📄 longgetutadlg.cpp
字号:
}
//////////////////////////////////////////////////////////////////////////
r=sqrt(XYZ[0]*XYZ[0]+XYZ[1]*XYZ[1]+XYZ[2]*XYZ[2]);
//////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////
if (CoorTime<=NowTime)
{
h=0.1;
//////////////////////////////////////////////////////////////////////////
for (j=0;j<=NowTime-CoorTime;j=j+h)
{
if (NowTime-CoorTime-j<h)
h=NowTime-CoorTime-j;
//////////////////////////////////////////////////////////////////////////
k1[0] =(-(GM)*XYZ[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ[0]*(1-(5*XYZ[2]*XYZ[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ[0]+2*we*XYZDot[1];
k1[1] =(-(GM)*XYZ[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ[1]*(1-(5*XYZ[2]*XYZ[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ[1]-2*we*XYZDot[0];
k1[2] =(-(GM)*XYZ[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ[2]*(3-(5*XYZ[2]*XYZ[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
XYZ2[0] =XYZ[0]+XYZDot[0]*h/2+k1[0]*h*h/8;
XYZDot2[0] =XYZDot[0]+k1[0]*h/2;
XYZ2[1] =XYZ[1]+XYZDot[1]*h/2+k1[1]*h*h/8;
XYZDot2[1] =XYZDot[1]+k1[1]*h/2;
XYZ2[2] =XYZ[2]+XYZDot[2]*h/2+k1[2]*h*h/8;
XYZDot2[2] =XYZDot[2]+k1[2]*h/2;
r =sqrt(XYZ2[0]*XYZ2[0]+XYZ2[1]*XYZ2[1]+XYZ2[2]*XYZ2[2]);
k2[0] =(-(GM)*XYZ2[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[0]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ2[0]+2*we*XYZDot2[1];
k2[1] =(-(GM)*XYZ2[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[1]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ2[1]-2*we*XYZDot2[0];
k2[2] =(-(GM)*XYZ2[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[2]*(3-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
XYZ2[0] =XYZ[0]+XYZDot[0]*h/2+k2[0]*h*h/8;
XYZDot2[0] =XYZDot[0]+k2[0]*h/2;
XYZ2[1] =XYZ[1]+XYZDot[1]*h/2+k2[1]*h*h/8;
XYZDot2[1] =XYZDot[1]+k2[1]*h/2;
XYZ2[2] =XYZ[2]+XYZDot[2]*h/2+k2[2]*h*h/8;
XYZDot2[2] =XYZDot[2]+k2[2]*h/2;
r =sqrt(XYZ2[0]*XYZ2[0]+XYZ2[1]*XYZ2[1]+XYZ2[2]*XYZ2[2]);
k3[0] =(-(GM)*XYZ2[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[0]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ2[0]+2*we*XYZDot2[1];
k3[1] =(-(GM)*XYZ2[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[1]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ2[1]-2*we*XYZDot2[0];
k3[2] =(-(GM)*XYZ2[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[2]*(3-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
XYZ2[0] =XYZ[0]+XYZDot[0]*h+k3[0]*h*h/2;
XYZDot2[0] =XYZDot[0]+k3[0]*h;
XYZ2[1] =XYZ[1]+XYZDot[1]*h+k3[1]*h*h/2;
XYZDot2[1] =XYZDot[1]+k3[1]*h;
XYZ2[2] =XYZ[2]+XYZDot[2]*h+k3[2]*h*h/2;
XYZDot2[2] =XYZDot[2]+k3[2]*h;
r =sqrt(XYZ2[0]*XYZ2[0]+XYZ2[1]*XYZ2[1]+XYZ2[2]*XYZ2[2]);
k4[0] =(-(GM)*XYZ2[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[0]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ2[0]+2*we*XYZDot2[1];
k4[1] =(-(GM)*XYZ2[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[1]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ2[1]-2*we*XYZDot2[0];
k4[2] =(-(GM)*XYZ2[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[2]*(3-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
for (i=0;i<3;i++){
XYZ[i] =XYZ[i]+XYZDot[i]*h+h*h*(k1[i]+2*k2[i]+2*k3[i]+k4[i])/12;
XYZDot[i] =XYZDot[i]+h*(k1[i]+2*k2[i]+2*k3[i]+k4[i])/6;
}
r =sqrt(XYZ[0]*XYZ[0]+XYZ[1]*XYZ[1]+XYZ[2]*XYZ[2]);
//////////////////////////////////////////////////////////////////////////
if (h==NowTime-CoorTime-j)
break;
}
}
else
{
h=-0.1;
//////////////////////////////////////////////////////////////////////////
for (j=0;j>=NowTime-CoorTime;j=j+h)
{
if (NowTime-CoorTime-j>h)
h=NowTime-CoorTime-j;
//////////////////////////////////////////////////////////////////////////
k1[0] =(-(GM)*XYZ[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ[0]*(1-(5*XYZ[2]*XYZ[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ[0]+2*we*XYZDot[1];
k1[1] =(-(GM)*XYZ[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ[1]*(1-(5*XYZ[2]*XYZ[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ[1]-2*we*XYZDot[0];
k1[2] =(-(GM)*XYZ[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ[2]*(3-(5*XYZ[2]*XYZ[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
XYZ2[0] =XYZ[0]+XYZDot[0]*h/2+k1[0]*h*h/8;
XYZDot2[0] =XYZDot[0]+k1[0]*h/2;
XYZ2[1] =XYZ[1]+XYZDot[1]*h/2+k1[1]*h*h/8;
XYZDot2[1] =XYZDot[1]+k1[1]*h/2;
XYZ2[2] =XYZ[2]+XYZDot[2]*h/2+k1[2]*h*h/8;
XYZDot2[2] =XYZDot[2]+k1[2]*h/2;
r =sqrt(XYZ2[0]*XYZ2[0]+XYZ2[1]*XYZ2[1]+XYZ2[2]*XYZ2[2]);
k2[0] =(-(GM)*XYZ2[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[0]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ2[0]+2*we*XYZDot2[1];
k2[1] =(-(GM)*XYZ2[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[1]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ2[1]-2*we*XYZDot2[0];
k2[2] =(-(GM)*XYZ2[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[2]*(3-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
XYZ2[0] =XYZ[0]+XYZDot[0]*h/2+k2[0]*h*h/8;
XYZDot2[0] =XYZDot[0]+k2[0]*h/2;
XYZ2[1] =XYZ[1]+XYZDot[1]*h/2+k2[1]*h*h/8;
XYZDot2[1] =XYZDot[1]+k2[1]*h/2;
XYZ2[2] =XYZ[2]+XYZDot[2]*h/2+k2[2]*h*h/8;
XYZDot2[2] =XYZDot[2]+k2[2]*h/2;
r =sqrt(XYZ2[0]*XYZ2[0]+XYZ2[1]*XYZ2[1]+XYZ2[2]*XYZ2[2]);
k3[0] =(-(GM)*XYZ2[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[0]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ2[0]+2*we*XYZDot2[1];
k3[1] =(-(GM)*XYZ2[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[1]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ2[1]-2*we*XYZDot2[0];
k3[2] =(-(GM)*XYZ2[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[2]*(3-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
XYZ2[0] =XYZ[0]+XYZDot[0]*h+k3[0]*h*h/2;
XYZDot2[0] =XYZDot[0]+k3[0]*h;
XYZ2[1] =XYZ[1]+XYZDot[1]*h+k3[1]*h*h/2;
XYZDot2[1] =XYZDot[1]+k3[1]*h;
XYZ2[2] =XYZ[2]+XYZDot[2]*h+k3[2]*h*h/2;
XYZDot2[2] =XYZDot[2]+k3[2]*h;
r =sqrt(XYZ2[0]*XYZ2[0]+XYZ2[1]*XYZ2[1]+XYZ2[2]*XYZ2[2]);
k4[0] =(-(GM)*XYZ2[0]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[0]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[0]+we*we*XYZ2[0]+2*we*XYZDot2[1];
k4[1] =(-(GM)*XYZ2[1]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[1]*(1-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[1]+we*we*XYZ2[1]-2*we*XYZDot2[0];
k4[2] =(-(GM)*XYZ2[2]/(r*r*r))+3*C20*GM*Ae*Ae*XYZ2[2]*(3-(5*XYZ2[2]*XYZ2[2]/(r*r)))/(2*(r*r*r*r*r))\
+XYZDDot[2];
//////////////////////////////////////////////////////////////////////////
for (i=0;i<3;i++){
XYZ[i] =XYZ[i]+XYZDot[i]*h+h*h*(k1[i]+2*k2[i]+2*k3[i]+k4[i])/12;
XYZDot[i] =XYZDot[i]+h*(k1[i]+2*k2[i]+2*k3[i]+k4[i])/6;
}
r =sqrt(XYZ[0]*XYZ[0]+XYZ[1]*XYZ[1]+XYZ[2]*XYZ[2]);
//////////////////////////////////////////////////////////////////////////
if (h==NowTime-CoorTime-j)
break;
}
}
AfxMessageBox("0k");
m_X=XYZ[0];m_Y=XYZ[1];m_Z=XYZ[2];
UpdateData(false);
}
void CLonggetutaDlg::OnButton1()
{
// TODO: Add your control notification handler code here
double VectorNow[StataN]={0};
double T;
double VOmiga=0;
double SOmiga=0;
long Seed;
for (int i=0;i<5000;i++)
{
VectorNow[0]=VectorNow[0]+VectorNow[2]*T;VectorNow[1]=VectorNow[1]+VectorNow[3]*T;
CarPoint[i][0]=VectorNow[0]; CarPoint[i][1]=VectorNow[1];
CarPoint[i][2]=VectorNow[2]; CarPoint[i][3]=VectorNow[3];
CarPoint[i][4]=VectorNow[4]; CarPoint[i][5]=VectorNow[5];
VOmiga=gauss(0,0.005*0.005,&Seed);
SOmiga=gauss(0,0.7*0.7,&Seed);
DRw[i-1]=(VectorNow[4]*VectorNow[3]-VectorNow[5]*VectorNow[2])/(VectorNow[2]*VectorNow[2]+VectorNow[3]*VectorNow[3])+VOmiga;
DRV[i-1]=T*sqrt(VectorNow[2]*VectorNow[2]+VectorNow[3]*VectorNow[3])+SOmiga;
}
}
//////////////////////////////////////////////////////////////////////////
//均匀分布随机数
//a 下限
//b 上限
//////////////////////////////////////////////////////////////////////////
double CLonggetutaDlg::uniform(double a,double b,long int *seed)
{
double t;
*seed=2045*(*seed)+1;
*seed=*seed-(*seed/1048576)*1048576;
t=(*seed)/1048576.0;
t=a+(b-a)*t;
return(t);
}
//////////////////////////////////////////////////////////////////////////
// 均匀分布生成高斯分布随机数
//mean 为数学期望
//sigma 方差因子
//s 随机种子 1-2147483647
//////////////////////////////////////////////////////////////////////////
double CLonggetutaDlg::gauss(double mean,double sigma,long int *s)
{
int i;
double x,y;
for(x=0,i=0;i<12;i++)
x+=uniform(0.0,1.0,s);
x=x-6.0;
y=mean+x*sigma;
return(y);
}⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -