#include #include #include #include #include #include "GPM.h" // 用户新增头文件 #include "OrionShip.h" using namespace std; using namespace arma; int main_Orion() { // 计时器开始 double dur; clock_t start, end; start = clock(); ////////////////////////////////////参数初始化////////////////////////////////////////////// input0 input; output0* output; auxdata1 auxdata; // 伪谱法条件初始化 Intput_Initial(&input); // 人为编辑定义部分 auxdata.Sref = 23.8; auxdata.M0 = 9500; auxdata.C1_v = 2e-8; // 与热流密度相关参数 auxdata.Rd_v = 0.006; // 与热流密度相关参数 auxdata.n_v = 0.5; // 与热流密度相关参数 auxdata.m_v = 3.15; auxdata.g = 9.80665; auxdata.Re = 6371110; input.setup.auxdata = auxdata; // 归一化单位向量 double pi = datum::pi; double g = auxdata.g; double Re = auxdata.Re; double tc = sqrt(Re / g); double rc = Re; double vc = sqrt(Re * g); vec state_c; state_c << vc << 1.0 << 1.0 << rc << 1.0 << 1.0 << endr; // 范围赋值 int length_state = 6; int length_control = 1; double initialtime_lower = 0; double initialtime_upper = 0; double finaltime_lower = 0; double finaltime_upper = 4500.0 / tc; vec initialstate_lower; initialstate_lower << 10600.0 << -5.9 / 180.0 * pi << 90.0 / 180.0 * pi << 120.0 * 1e3 + Re << 0.0 / 180.0 * pi << 0.0 / 180.0 * pi << endr; vec initialstate_upper; initialstate_upper << 10600.0 << -5.9 / 180.0 * pi << 90.0 / 180.0 * pi << 120.0 * 1e3 + Re << 0.0 / 180.0 * pi << 0.0 / 180.0 * pi << endr; vec state_lower; state_lower << 100.0 << -80.0 / 180.0 * pi << -90.0 / 180.0 * pi << 10.0 * 1e3 + Re << -1.0 / 180.0 * pi << -90.0 / 180.0 * pi << endr; vec state_upper; state_upper << 11700.0 << +20.0 / 180.0 * pi << 270.0 / 180.0 * pi << 200.0 * 1e3 + Re << 180.0 / 180.0 * pi << +90.0 / 180.0 * pi << endr; vec finalstate_lower; finalstate_lower << 100.0 << -80.0 / 180.0 * pi << -90.0 / 180.0 * pi << 10.0 * 1e3 + Re << 90.0 / 180.0 * pi << -90.0 / 180.0 * pi << endr; vec finalstate_upper; finalstate_upper << 300.0 << + 0.0 / 180.0 * pi << 270.0 / 180.0 * pi << 10.0 * 1e3 + Re << 90.0 / 180.0 * pi << +90.0 / 180.0 * pi << endr; vec control_lower; control_lower << -180.0 / 180.0 * pi << endr; vec control_upper; control_upper << +180.0 / 180.0 * pi << endr; vec path_lower; path_lower << -6.5 << endr; vec path_upper; path_upper << 6.5 << endr; rowvec integral_lower; integral_lower.reset(); rowvec integral_upper; integral_upper.reset(); input.setup.state_length = length_state; input.setup.control_length = length_control; input.setup.bounds.phase.initialtime.lower = initialtime_lower; input.setup.bounds.phase.initialtime.upper = initialtime_upper; input.setup.bounds.phase.finaltime.lower = finaltime_lower; input.setup.bounds.phase.finaltime.upper = finaltime_upper; input.setup.bounds.phase.initialstate.lower = initialstate_lower / state_c; input.setup.bounds.phase.initialstate.upper = initialstate_upper / state_c; input.setup.bounds.phase.state.lower = state_lower / state_c; input.setup.bounds.phase.state.upper = state_upper / state_c; // 默认状态上限为正无穷 input.setup.bounds.phase.finalstate.lower = finalstate_lower / state_c; input.setup.bounds.phase.finalstate.upper = finalstate_upper / state_c; input.setup.bounds.phase.control.lower = control_lower; input.setup.bounds.phase.control.upper = control_upper; input.setup.bounds.phase.path.lower = path_lower; input.setup.bounds.phase.path.upper = path_upper; input.setup.bounds.phase.integral.lower = integral_lower; input.setup.bounds.phase.integral.upper = integral_upper; input.setup.bounds.phase.event.lower.reset(); input.setup.bounds.phase.event.upper.reset(); // 猜测值 guess0 guess; int length_path = 1; vec guess_time; guess_time << 0.0 << 2500.0 / tc << endr; mat guess_state(2, length_state, fill::zeros); mat guess_control(2, length_control, fill::zeros); mat guess_path(2, length_path, fill::zeros); guess_state.row(0) = input.setup.bounds.phase.initialstate.lower.t(); guess_state.row(1) = input.setup.bounds.phase.finalstate.lower.t(); guess_control.row(0) = input.setup.bounds.phase.control.lower.t(); guess_control.row(1) = input.setup.bounds.phase.control.upper.t(); guess_path.row(0) = input.setup.bounds.phase.path.lower.t(); guess_path.row(1) = input.setup.bounds.phase.path.upper.t(); guess.phase.time = guess_time; guess.phase.state = guess_state; guess.phase.control = guess_control; guess.phase.path = guess_path; input.setup.guess = guess; // 配点设置 int inter_num = 5; input.setup.name = "OrionShip"; // 问题 mesh1 mesh; mesh.method = "p-h"; // 网格细化方法的名字,字符串表示 mesh.tolerance = 1e-5; // 网格细化精度的允许值 mesh.maxiterration = 6; // 网格细化最大迭代次数 mesh.colpointsmin = 2; // 网格细化产生新interval下的最小配点数目 mesh.colpointsmax = 20; // 网格细化产生新interval下的最大配点数目 mesh.splitmult = 1.2; // 比例因子,确定创建新网格的间隔数目 mesh.curveratio = 2; // 网格间隔中解的最大平均曲率比的阈值,以确定是否增加新间隔 mesh.phase.inter_num = inter_num; // 段p内的间隔数目, 整数 mesh.phase.colpoints = 2 * ones(inter_num, 1); // 段p内的配置点数目,整数数组,长度与间隔数目相等 mesh.phase.i_proportion = 2.0 / inter_num * ones(inter_num, 1); // 初始网格划分的比例(当前为均等划分,和为2) input.setup.mesh = mesh; // 组合setup.mesh结构体 // input.setup.Scal = 0; input.setup.warmstart = 0; input.setup.interp.method = "linear"; // 函数设置 functions1 functions; functions.continous1 = Reentry_Orion_Dynamic; functions.endpoint = OrionObjective; // 求解器设置 NLP1 NLP; NLP.ipoptoptions.print_level = 1; NLP.ipoptoptions.hessian_approximation = "limited-memory"; NLP.ipoptoptions.mu_strategy = "adaptive"; NLP.ipoptoptions.linear_solver = "ma27"; //mumps、pardiso NLP.ipoptoptions.tolerance = 1e-1; NLP.ipoptoptions.maxiterration = 3000; input.setup.NLP = NLP; input.setup.derivatives.supplier = "CD"; //input.setup.feature.FSE = 1; /////// // 伪谱法过程实现 int iteration_final; output = PM_Process(&input, functions, &iteration_final); // 计时结束 end = clock(); dur = (double)(end - start); #ifndef LUOYC20220721 printf("Use Time: %f \n", (dur / CLOCKS_PER_SEC)); #else printf_s("Use Time: %f \n", (dur / CLOCKS_PER_SEC)); #endif getchar(); return 1; }