内弹道解算程序C++

内弹道解算程序C++
附录程序
#define _CRT_SECURE_NO_DEPRECA TE
#include
#include
#include
#include
#include
#include //C head file
#include
#include
#include
#include
using namespace std;
/*初始参数*/
//////////燃烧室内参数
//火药燃烧的参数
#define U1 0.2 //燃速系
#define an 0.82 //燃速指数
#define omga 1.76 //装药质量
#define am 21.76 //弹丸质量
#define f 980000 //火药力
#define F1 1.1 //阻力系数
#define STA 0.2 //（k-1）
#define DP 1600 //火药装填密度
#define ALF 1.0e-3 //火药燃气余容double PG=30e6 ; //挤进压力
//火药参数（单基9/7）
#define HE1 0.5e-3 //第二种火药半弧厚
#define HD0 0.5e-3 //第二种火药内孔直径
#define HD2C 12e-2 //第二种火药火药长
//火药形状特征量Π1、Q1、β
double XHT,XHQ,XHB;
//double XHB;
//碎粒燃烧结束时的Zb，火药形状参数χ、λ和μdouble Zb0,CI,ALM,AMU;
//碎粒燃烧阶段的形状参数χs、λs
double CIS,ALMS;
//////////火炮参数
#define V0 4.2e-3
#define omga2 0.34 //第一种药质量
#define f2 980000 //第一种药火药力
#define DP2 1600 //第一种药密度
#define ALF2 1e-3 //第一种药余容
//火药参数（选用杆状药）
#define HE12 0.24e-3 //火药半弧厚
#define HD02 0.3e-3 //内孔直径
#define HD2C2 6.5e-3 //火药长
//#define f2 980000
double Z2,PSI2,dZ2,Zb2,Zm2;
double DLT2;
//火药形状参数χ、λ和μ
double CI2,ALM2,AMU2;
double XHB2;
#define pi 3.1415926 //圆周率
#define atm 1.013e5 //大气压
#define h 1e-7 //步长
int i,i0,j,k,ii;
//次要功系数、身管截面积，燃烧室当量长度double FI,s,L0;
//火药的装填密度
double DLT;
double aaa;
double NST,HIB;
//火药燃烧时间
double TIME,TIMm,TIMb;
//PTDD：膛底压力
double PTDD;
//弹后空间平均压力、弹底压力double P,PD;
double PSI,ZJSN;
//double PSI2;
double PV,LPSI;
long NN;
//弹丸运动距离（不大于身管长度）double L,Lb,Lm,dL;
//弹丸速度
double v,vm,vb,dv;
int nn;
//////////////微分过程参数double a1[6];
double Z,dZ,Zb,Zm;
void ranshaoshi();
void chushitiaojian();
void ranshao_0();
void ranshao_1();
void thy();
void main()
{
FILE *outhy;
ranshaoshi();
NST=5;
HIB=NST*h;
outhy=fopen("hy0.txt","w");
fprintf(outhy,"omga(kg)\tFI\tDLT\n");
fprintf(outhy,"%f\t%f\t%f\n",
omga,FI,DLT*1e-3);
fprintf(outhy,"t（ms）\tZ\tψ\tZ2\tψ2\tPTDD（MPa）\tP （MPa）\tPD（MPa）\tv（m/s）\tLψ（m）\tL（m）\n");
//初始燃烧时条件
P=5e6;
TIME=0;
chushitiaojian();
ranshao_0();
while(L<=2.384)
{
thy();
fprintf(outhy,"%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\t%f\ n",
TIME*1e3,Z,PSI,Z2,PSI2,PTDD*1e-6,P*1e-6,PD*1e-6,v,LPSI,L);
}
fclose(outhy);
printf("FI is %f\n",FI);
printf("v is %f\n",v);
printf("Finished!");
getchar();
}
void ranshaoshi()
{
FI=1.132;
s=119.6e-4;
L0=V0/s;
DLT=omga/V0;
DLT2=omga2/V0;
//////9/7多孔药
Zb0=1+0.2956*(1+HD0/2/HE1);
XHT=((3*HD0+8*HE1)+7*HD0)/HD2C;
XHQ=(pow((3*HD0+8*HE1),2)-7*HD0*HD0)/(HD2C*HD2C); XHB=2*HE1/HD2C;
CI=(2*XHT+XHQ)/XHQ*XHB;
ALM=(6-2*XHT)/(2*XHT+XHQ)*XHB;
AMU=-1*6/(2*XHT+XHQ)*XHB*XHB;
CIS=(1-CI*(1+ALM+AMU)*Zb0*Zb0)/Zb0/(1-Zb0);
ALMS=CI*(1+ALM+AMU)/CIS-1;
////////杆状药
XHB2=2*HE12/HD2C2;
CI2=1+XHB2;
ALM2=-1*XHB2/(1+XHB2);
AMU2=0;
}
void chushitiaojian()
{
double segma,segma2;
PSI=(1/DLT-1/DP)/(f/P+ALF-1/DP);
segma=sqrt(1+4*ALM/CI*PSI);
Z=2*PSI/CI/(1+segma);
PSI2=(1/DLT2-1/DP)/(f2/P+ALF-1/DP);
segma2=sqrt(1+4*ALM2/CI2*PSI2);
Z2=2*PSI2/CI2/(1+segma2);
}
void thy()
{
a1[1]=0.5*HIB;
a1[2]=a1[1];
a1[3]=HIB;
a1[4]=HIB;
a1[5]=a1[1];
Zm=Z;
Zb=Z;
Zm2=Z2;
Zb2=Z2;
Lm=L;
Lb=L;
vm=v;
vb=v;
TIMm=TIME;
TIMb=TIME;
for (j=1; j<5; j++)
{
ranshao_1();
Z=Zm+a1[j]*dZ;
Zb=Zb+a1[j+1]*dZ/3;
Z2=Zm2+a1[j]*dZ2;
Zb2=Zb2+a1[j+1]*dZ2/3; L=Lm+a1[j]*dL;
Lb=Lb+a1[j+1]*dL/3;
v=vm+a1[j]*dv;
vb=vb+a1[j+1]*dv/3; TIME=TIMm+a1[j]; TIMb=TIMb+a1[j+1]/3;
if (j==4)
{
Z=Zb;
Z2=Zb2;
L=Lb;
v=vb;
TIME=TIMb;
}
ranshao_0();
}
}
void ranshao_0()
{
PSI=1;
//主体燃烧阶段
if (Z<1)
{
PSI=CI*Z*(1+ALM*Z+AMU*Z*Z);
}
//碎粒燃烧阶段
if (Z>=1 && Z<=Zb0)
{
PSI=CIS*Z*(1+ALMS*Z);
}
if (Z2<1)
{
PSI2=CI2*Z2*(1+ALM2*Z2+AMU2*Z2*Z2);
}
LPSI=L0*(1-DLT/DP*(1-PSI)-DLT2/DP2*(1-PSI2)-ALF*DLT*PSI-ALF2*DLT2*PSI2);
PV=f*omga*PSI+f2*omga2*PSI2-STA*FI*am*v*v/2; P=PV/(s*(L+LPSI));
PTDD=P*(1+(omga+omga2+0.03)/(2*F1*am))/FI; }
void ranshao_1()
{
dZ=0;
if (Z<zb0)< p="">
{
dZ=U1*pow(P,an)/pow(1e6,an)/HE1/100;
}
dZ2=0;
if (Z2<1)
{
dZ2=U1*pow(P,an)/pow(1e6,an)/HE12/100;
}
dL=v;
dv=0;
PD=P*F1/FI;
if (PD>PG)
{
dv=s*P/(FI*am);
nn=1;
}
if(nn==1)
{
dv=s*P/(FI*am);
}
}
</zb0)<>。