行星轨迹制作_行星运动轨迹的程序实现

àì夳堔傛蜴生んèń 2023-01-10 03:49 109阅读 0赞

//三体问题的导数曲线

class ThreeBody : publicDifferentiationFunction

\{public:

ThreeBody()

m\_star1\_x= -10.0f;

m\_star1\_y= 0.0f;

m\_star1\_z= 0.0f;

m\_star1\_weight= 1.0f;

m\_star2\_x= 10.0f;

m\_star2\_y= 0.0f;

m\_star2\_z= 0.0f;

m\_star2\_weight= 1.0f;

m\_planet\_x= 5.0f;

m\_planet\_y= 5.0f;

m\_planet\_z= 0.1f;

m\_planet\_weight= 0.1f;

m\_planet\_speed\_x= 0.0f;

m\_planet\_speed\_y= 2.0f;

m\_planet\_speed\_z= 0.0f;

m\_g= 50.0f;

m\_ek= 0.5f\*m\_planet\_weight\*(m\_planet\_speed\_x\*m\_planet\_speed\_x + m\_planet\_speed\_y\*m\_planet\_speed\_y + m\_planet\_speed\_z\*m\_planet\_speed\_z); //1/2\*m\*v\*v

float d1x = m\_star1\_x -m\_planet\_x;float d1y = m\_star1\_y -m\_planet\_y;float d1z = m\_star1\_z -m\_planet\_z;float sqd1 = d1x\*d1x + d1y\*d1y + d1z\*d1z;float d1 =sqrt(sqd1);

m\_ep1= -m\_g\*m\_star1\_weight\*m\_planet\_weight/d1;float d2x = m\_star2\_x -m\_planet\_x;float d2y = m\_star2\_y -m\_planet\_y;float d2z = m\_star2\_z -m\_planet\_z;float sqd2 = d2x\*d2x + d2y\*d2y + d2z\*d2z;float d2 =sqrt(sqd2);

m\_ep2= -m\_g\*m\_star2\_weight\*m\_planet\_weight/d2;

m\_e= m\_ek + m\_ep1 +m\_ep2;

\}void Differentiate(float x, float y, float z, floatt,float& outX, float& outY, float&outZ)

t= t\*20.0f;float d1x = m\_star1\_x -x;float d1y = m\_star1\_y -y;float d1z = m\_star1\_z -z;float sqd1 = d1x\*d1x + d1y\*d1y + d1z\*d1z;float d1 =sqrt(sqd1);float d2x = m\_star2\_x -x;float d2y = m\_star2\_y -y;float d2z = m\_star2\_z -z;float sqd2 = d2x\*d2x + d2y\*d2y + d2z\*d2z;float d2 =sqrt(sqd2);float a1 = m\_g\*m\_star1\_weight/sqd1;float a1x = a1\*d1x/d1;float a1y = a1\*d1y/d1;float a1z = a1\*d1z/d1;float a2 = m\_g\*m\_star2\_weight/sqd2;float a2x = a2\*d2x/d2;float a2y = a2\*d2y/d2;float a2z = a2\*d2z/d2;

outX= x + m\_planet\_speed\_x\*t + 0.5f\*(a1x + a2x)\*t\*t;

outY= y + m\_planet\_speed\_y\*t + 0.5f\*(a1y + a2y)\*t\*t;

outZ= z + m\_planet\_speed\_z\*t + 0.5f\*(a1z + a2z)\*t\*t;

m\_planet\_speed\_x+= (a1x + a2x)\*t;

m\_planet\_speed\_y+= (a1y + a2y)\*t;

m\_planet\_speed\_z+= (a1z + a2z)\*t;

\{float d1x = m\_star1\_x -outX;float d1y = m\_star1\_y -outY;float d1z = m\_star1\_z -outZ;float sqd1 = d1x\*d1x + d1y\*d1y + d1z\*d1z;float d1 =sqrt(sqd1);

m\_ep1= -m\_g\*m\_star1\_weight\*m\_planet\_weight/d1;float d2x = m\_star2\_x -outX;float d2y = m\_star2\_y -outY;float d2z = m\_star2\_z -outZ;float sqd2 = d2x\*d2x + d2y\*d2y + d2z\*d2z;float d2 =sqrt(sqd2);

m\_ep2= -m\_g\*m\_star2\_weight\*m\_planet\_weight/d2;

m\_ek= m\_e - m\_ep1 -m\_ep2;if (m\_ek < 0.0f)

m\_ek= 0.0f;

\}float v = sqrt(2\*m\_ek/m\_planet\_weight);float w = sqrt(m\_planet\_speed\_x\*m\_planet\_speed\_x + m\_planet\_speed\_y\*m\_planet\_speed\_y + m\_planet\_speed\_z\*m\_planet\_speed\_z);

m\_planet\_speed\_x\*= v/w;

m\_planet\_speed\_y\*= v/w;

m\_planet\_speed\_z\*= v/w;

\}float GetExtendT() const\{return 20.0f;

\}float GetStartX() const\{returnm\_planet\_x;

\}float GetStartY() const\{returnm\_planet\_y;

\}float GetStartZ() const\{returnm\_planet\_z;

\}public:floatm\_star1\_x;floatm\_star1\_y;floatm\_star1\_z;floatm\_star1\_weight;floatm\_star2\_x;floatm\_star2\_y;floatm\_star2\_z;floatm\_star2\_weight;floatm\_planet\_x;floatm\_planet\_y;floatm\_planet\_z;floatm\_planet\_weight;floatm\_planet\_speed\_x;floatm\_planet\_speed\_y;floatm\_planet\_speed\_z;float m\_g; //万有引力系数

floatm\_e;float m\_ek; //动能

float m\_ep1; //引力势能

float m\_ep2; //引力势能

\};