📄 标准弹道.cpp
字号:
/******************************************************************************/
/* 三自由度包络计算模块 */
/* */
/* 作者：胡锐 */
/* */
/* 时间：2006年5月22日 */
/******************************************************************************/
#include <math.h>
#include <stdio.h>
#include <malloc.h>
#include <fstream.h> //输入输出流
double ACX[7][5]={{0.1659 , 0.1725 , 0.1963 , 0.3507 , 0.5381}, //0
{0.1809 , 0.1859 , 0.2092 , 0.3674 , 0.5486}, //2
{0.2066 , 0.2141 , 0.2375 , 0.4014 , 0.5841}, //4
{0.2525 , 0.2641 , 0.2902 , 0.4536 , 0.6408}, //6
{0.3140 , 0.3296 , 0.3636 , 0.5310 , 0.7293}, //8
{0.4661 , 0.4943 , 0.5395 , 0.7143 , 0.9184}, //11
{0.6846 , 0.7332 , 0.7911 , 0.9770 , 1.1594}}; //14 阻力系数CX
double ACY0[2][5]={{ 0.6, 0.8, 0.9, 1.0, 1.15},
{ 0.08028, 0.08316,0.08506 ,0.08767 , 0.09362 }}; //升力系数CY0
double ACY1[2][5]={{ 0.6, 0.8, 0.9, 1.0, 1.15},
{ 2.24162e-4, 2.39065e-4, 2.47800e-4, 2.36159e-4, 1.81134e-4}}; //升力系数CY1
double Mzdetla1[2][5]={{ 0.6, 0.8, 0.9, 1.0, 1.1}, //俯仰诱导系数
{0.1306,0.1347,0.1379,0.1560,0.1421}};
double Mzalpha[7][5]={{-0.0484,-0.0486,-0.0498,-0.0579,-0.0665}, //攻角0度
{-0.06357,-0.06277,-0.06462,-0.07577,-0.08510}, //2
{-0.08634,-0.09101,-0.08981,-0.10069,-0.1172}, //4
{-0.11237,-0.11807,-0.11624,-0.12293,-0.14747}, //6
{-0.13840,-0.13884,-0.13876,-0.13802,-0.16468}, //8
{-0.15742,-0.15094,-0.15642,-0.14492,-0.17674}, //10
{-0.16331,-0.16147,-0.16290,-0.15236,-0.18441}}; //12
double Mzdetla[7][5]={{-0.1306,-0.1347,-0.1379,-0.156,-0.1421}, //攻角0
{-0.12911,-0.13781,-0.14806,-0.16250,-0.14855}, //2
{-0.13314,-0.14157,-0.15050,-0.16712,-0.15432}, //4
{-0.13457,-0.14461,-0.15454,-0.17047,-0.15807}, //6
{-0.13870,-0.14912,-0.15816,-0.17657,-0.16393}, //8
{-0.14348,-0.15382,-0.16144,-0.18164,-0.17088}, //10
{-0.14996,-0.15897,-0.16849,-0.18717,-0.17770}}; //12
double Y[6],SONIC,RHO;
double GJ[7]={0.000,2.000,4.000,6.000,8.000,10.000,12.000}; //alpha攻角
double MHS[5]={0.6,0.8,0.9,1.0,1.15}; //Ma值
double MA,ABC_ALPHA,CX,CY0,CY1,ALPHA;
double dX=0, dY=0, dTheta=0,guozai=0;
double deta_X,deta_Y,deta_V; //投弹时的误差
FILE *fo;
void rk(int ,double,double ); //四阶龙格库塔法子程序
void dery(double* ,double*,double ); //计算微分方程组右端函数值
void interp(); //调用插值函数
double int32(double ,double ,int ,int ,double GJ[],double MHS[],double AC[][5]); //双变元不等距线性插值
double int31(double A[],double XS[][5],int N1,int N2,int ,double ); //单变元不等距抛物线插值
double int11(double yy[][5],int n, double x);
void solve(double V, double H, double ALPHA, double deta_x, double deta_y, double deta_v);
void result(double);
/*************************************************************************************/
/* 主函数设置参数计算弹道 */
/* main() */
/*************************************************************************************/
void main(void)
{
/*打开文件输出结果的目录*/
fo=fopen("标准弹道计算.txt","w");
/* fprintf(fo,"%s"," Time Ma ALPHA X dX Y dY V THETA dTHETA guozai");
fprintf(fo,"\n");
fprintf(fo,"%s"," s 。 m m/s m m/s m/s rad rad/s ");
fprintf(fo,"\n");
fprintf(fo,"\n");
fprintf(f1,"%s"," 攻角 时间 射程 垂直速度 过载 合速度");
fprintf(f1,"\n");
fprintf(f1,"%s"," 。 s m m/s ");
fprintf(f1,"\n");
*/
solve(800,4000, 0, 0, 0, 0);
fclose(fo);
}
/*************************************************************************************/
/* 调用进行弹道计算 */
/* solve() */
/* 无返回值 */
/*************************************************************************************/
void solve(double V,double H, double ALPHA, double deta_x, double deta_y, double deta_v)
{
double h=0.02,t2000,t75;
double alpha1[4]={0,-6.0,0,0};
Y[0]=0; //时间 *
Y[1]=sqrt(((V+deta_v)/3.6)*((V+deta_v)/3.6)+16); //合成速度 *
Y[2]=-atan(14.4/(V+deta_v)); //弹道倾角 ***置初始值
Y[3]=deta_x; //射程 *
Y[4]=H+deta_y; //高度 *
while(Y[4]>=2000)
{
// for(int i=0;i<=10000;i++) //可进行不同步长时的定时间值取出
// {
// if(fabs(Y[0]-0.02*i)<=0.00000001)
// result(ALPHA);
// }
rk(5,h,ALPHA); //调用龙格库塔法
result(ALPHA);
// ALPHA=ALPHA0;
}
t2000=Y[0];
while(Y[2]>=-1.1&& Y[4]>=0.0)
{
ALPHA=alpha1[1]+(alpha1[0]-alpha1[1])*exp(-2.0*(Y[0]-t2000));
rk(5,h,ALPHA);
result(ALPHA);
}
t75=Y[0];
alpha1[2]=ALPHA;
while(Y[4]>=0)
{
ALPHA=alpha1[3]+(alpha1[2]-alpha1[3])*exp(-2.0*(Y[0]-t75));
rk(5,h,ALPHA);
result(ALPHA);
}
}
/*************************************************************************************/
/* 四阶龙格库塔法子程序 */
/* rk() */
/* 无返回值 */
/*************************************************************************************/
void rk(int n,double h,double ALPHA)
{
double a[4],old_y[5],Y1[5],*dy;
int i,j;
dy=(double*)calloc(n,sizeof(double)); //定义一个n长的空间
a[0]=a[1]=h/2;
a[2]=a[3]=h;
dery(dy,Y,ALPHA);
for(i=0;i<n;i++) {old_y[i]=Y[i];}
for(j=0;j<3;j++)
{
for(i=0;i<n;i++)
{
Y1[i]=old_y[i]+a[j]*dy[i];
Y[i]=Y[i]+a[j+1]*dy[i]/3;
}
dery(dy,Y1,ALPHA);
}
for(i=0;i<n;i++)
Y[i]=Y[i]+a[0]*dy[i]/3;
free(dy);
return;
}
/*************************************************************************************/
/* 计算微分方程组右端函数值 */
/* dery() */
/* 无返回值，用于获取龙格库塔法中的K1-K4 */
/*************************************************************************************/
void dery(double dy[],double Y[],double ALPHA)
{
double q,a[2],XF,YF,g,S,L;
L=2.2; //弹长
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -