首页 > 分享 > c++模拟质点在万有引力下的运动轨迹

c++模拟质点在万有引力下的运动轨迹

萌宠菠菠乐园
2024-12-16 07:28

P h y s i c s I S t h e M o s t B e a u t i f u l mathcal{Physics IStheMostBeautiful} PhysicsIStheMostBeautiful

通过c++模拟质点在万有引力下的运动轨迹

0. F o r w a r d I t r o d u c t i o n mathcal{ForwardItroduction} ForwardItroduction

应某集团首领要求：

通过代码实现天体间万有引力作用的模拟运动轨迹
时限一晚上

通过模拟键盘和鼠标的点击在画图中绘制图像

1. P r e K n o w l e d g e 1 → T o C o n t r o l K e y b o a r d mathcal{PreKnowledge1to ToControlKeyboard} PreKnowledge1→ToControlKeyboard

直接上代码

keybd_event(num/*Key_Code*/, 0, 0, 0);//按下按键 keybd_event(num/*Key_Code*/, 0, KEYEVENTF_KEYUP, 0);//放开按键 12

num是你所要按下的键的键码，字母的键码就是本身的Ascll码，键码表如下（来自网络）：
cnblogs-万一

2. P r e K n o w l e d g e 2 → T o C o n t r o l M o u s e mathcal{PreKnowledge2to ToControlMouse} PreKnowledge2→ToControlMouse

仍然上代码

POINT p; GetCursorPos(&p);//获取鼠标坐标 printf("%d %dn", p.x, p.y); SetCursorPos(x, y);//移动到(x,y) mouse_event(MOUSEEVENTF_LEFTDOWN, x, y, 0, 0);//按下左键 mouse_event(MOUSEEVENTF_LEFTUP, x, y, 0, 0);//放开左键 mouse_event(MOUSEEVENTF_RIGHTDOWN, x, y, 0, 0);//按下右键 mouse_event(MOUSEEVENTF_RIGHTUP, x, y, 0, 0);//放开右键 123456789

鼠标的坐标 ( p . x , p . y ) (p.x,p.y) (p.x,p.y)以像素为单位，且一定为整数

3. P r e K n o w l e d g e 3 → U n i v e r s a l G r a v i t a t i o n ′ s L a w mathcal{PreKnowledge3to UniversalGravitation'sLaw} PreKnowledge3→UniversalGravitation′sLaw

根据中学知识知道
F g = G m 1 m 2 R 2 F_g=frac{Gm_1m_2}{R^2} Fg=R2Gm1m2
G ≈ 6.67 × 1 0 − 11 N ⋅ m 2 ⋅ k g − 2 Gapprox6.67times 10^{-11}Ncdot m^2cdot kg^{-2} G≈6.67×10−11N⋅m2⋅kg−2
因为 G G G的精度要求较大，所以可以适当缩放单位大小，如：

m → k m mto km m→km, k g → 10000 t kgto 10000t kg→10000t

4. P r e K n o w l e d g e 4 → N e w t o n ′ s S e c o n d L a w mathcal{PreKnowledge4to Newton'sSecondLaw} PreKnowledge4→Newton′sSecondLaw

F = m a F=ma F=ma

5. S o D o I t → P h y s i c s mathcal{SoDoItto Physics} SoDoIt→Physics

我们可以把时间划分成很小很小的单位，如 10 m s 10ms 10ms，每次以上一次的坐标近似表示这 10 m s 10ms 10ms内位置
F = G m 1 m 2 R 2 F=frac{Gm_1m_2}{R^2} F=R2Gm1m2用 O ( n 2 ) mathcal{O(n^2)} O(n2)的复杂度两两之间的引力
再算出每个质点的向量和 a = F m , Δ v = a t , Δ x = v t a=frac{F}{m},Delta v=at,Delta x=vt a=mF,Δv=at,Δx=vt( t t t取 0.01 s 0.01s 0.01s)

6. T u r n I n t o C o d e mathcal{TurnIntoCode} TurnIntoCode

此代码仅兼容部分设备qwq

#include <cstdio> #include <iostream> #include <cmath> #include <windows.h> #include <winuser.h> #define N 100 #define KEY_DOWN(VK_NONAME) ((GetAsyncKeyState(VK_NONAME) & 0x8000) ? 1:0) char st[100]; const int ox = 795, oy = 492, unt = 10; const long double G = 0.667408; struct Node{long double x, y, vx, vy, fx, fy, m; }a[N]; inline void cli_key(int id) {keybd_event(id, 0, 0, 0);keybd_event(id, 0, KEYEVENTF_KEYUP, 0); } inline void cli_mouse(int x, int y) {SetCursorPos(x, y);mouse_event(MOUSEEVENTF_LEFTDOWN, x, y, 0, 0);mouse_event(MOUSEEVENTF_LEFTUP, x, y, 0, 0); } inline void prt_st(char *st) {for (int i = 0; i < strlen(st); ++i) {int id;if (st[i] >= 'a' && st[i] <= 'z') id = st[i] - 'a' + 'A';else id = st[i];cli_key(id);} } inline void prt_line(int x1, int y1, int x2, int y2) {SetCursorPos(x1, y1);mouse_event(MOUSEEVENTF_LEFTDOWN, x1, y1, 0, 0);SetCursorPos(x2, y2);mouse_event(MOUSEEVENTF_LEFTUP, x2, y2, 0, 0); } inline void cli_pencil() {cli_mouse(269, 78); } inline void getxy() {Sleep(5000);POINT p;GetCursorPos(&p);printf("%d %dn", p.x, p.y); } int main() {int n, time;puts("time(s)?");scanf("%d", &time);puts("The scale is 100metres");printf("The number of particle(n <= 20):");scanf("%d", &n);for (int i = 1; i <= n; ++i) {printf("The X & Y(*100m) of %d-th particle:", i);std :: cin >> a[i].x >> a[i].y;printf("The Vx & Vy(m/s) of %d-th particle:", i);std :: cin >> a[i].vx >> a[i].vy;printf("The M(*10000000kg) of %d-th particle:", i);std :: cin >> a[i].m;}keybd_event(91, 0, 0, 0);keybd_event('R', 0, 0, 0);keybd_event('R', 0, KEYEVENTF_KEYUP, 0);keybd_event(91, 0, KEYEVENTF_KEYUP, 0);Sleep(1000);prt_st("mspaint");cli_key(13);Sleep(1000);keybd_event(18, 0, 0, 0);keybd_event(32, 0, 0, 0);keybd_event('X', 0, 0, 0);keybd_event('X', 0, KEYEVENTF_KEYUP, 0);keybd_event(32, 0, KEYEVENTF_KEYUP, 0);keybd_event(18, 0, KEYEVENTF_KEYUP, 0);Sleep(1000);keybd_event(17, 0, 0, 0);keybd_event('E', 0, 0, 0);keybd_event('E', 0, KEYEVENTF_KEYUP, 0);keybd_event(17, 0, KEYEVENTF_KEYUP, 0);Sleep(1000);prt_st("1580");cli_key(9);prt_st("660");cli_key(13);cli_pencil(); //getxy();prt_line(5, 492, 1585, 492);prt_line(795, 162, 795, 822);prt_line(ox + unt, oy, ox + unt, oy - 7);for (int k = 1; k <= time / 0.01; ++k) {for (int i = 1; i <= n; ++i) {a[i].fx = a[i].fy = 0;for (int j = 1; j <= n; ++j)if (i != j) {long double sqrR = (a[i].x - a[j].x) * (a[i].x - a[j].x) + (a[i].y - a[j].y) * (a[i].y - a[j].y);long double F = G * (a[i].m * a[j].m) / sqrR;if (F > 10000000) F = 10000000;a[i].fx += F * (a[j].x - a[i].x) / sqrt(sqrR);a[i].fy += F * (a[j].y - a[i].y) / sqrt(sqrR);}}for (int i = 1; i <= n; ++i) {printf("%.4Lf %.4Lf %.4Lf %.4Lf %.4Lf %.4Lf %.4Lfn", a[i].x, a[i].y, a[i].vx, a[i].vy, a[i].m, a[i].fx, a[i].fy);a[i].vx += a[i].fx / a[i].m * 0.01;a[i].vy += a[i].fy / a[i].m * 0.01;long double xx = a[i].x + a[i].vx * 0.01;long double yy = a[i].y + a[i].vy * 0.01;int xx1 = (int)(a[i].x * unt) + ox, yy1 = (int)(a[i].y * unt) + oy, xx2 = (int)(xx * unt) + ox, yy2 = (int)(yy * unt) + oy;int flg = 0;if (xx2 <= 10) {a[i].x = 78;flg = 1;}if (xx2 >= 1580) {a[i].x = -78;flg = 1;}if (yy2 <= 167) {a[i].y = 32;flg = 1;}if (yy2 >= 817) {a[i].y = -32;flg = 1;}if (flg) continue;prt_line(xx1, yy1, xx2, yy2);a[i].x = xx;a[i].y = yy;}for (int i = 1; i <= 1 * 100000; ++i) {if (KEY_DOWN(27)) {Sleep(300);while (!KEY_DOWN(27));Sleep(300);break;}} //Sleep(10);}return 0; }

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152 7. P i c t u r e s mathcal{Pictures} Pictures

data1

30
3
0 0
0 0
1000000
20 0
0 30
20000
-20 0
0 -30
20000

data2