Какие алгоритмы позволят мне симулировать физику планет?

Question

Я заинтересован в создании симулятора "Солнечная система", который позволит мне моделировать вращательные и гравитационные силы планет и звезд.

Я бы хотел сказать, моделировать нашу солнечную систему и моделировать ее на разных скоростях (то есть наблюдать, как Земля и другие планеты вращаются вокруг Солнца в течение дней, лет и т. Д.). Я хотел бы иметь возможность добавлять планеты и изменять массу планет и т. Д., Чтобы посмотреть, как это повлияет на систему.

Есть ли у кого-нибудь ресурсы, которые укажут мне правильное направление для написания симулятора такого типа?

Существуют ли какие-либо физические движки, разработанные для этой цели?

Answer 1

Здесь все здесь и вообще все, что написала Джин Миус.

Answer 2

Вам нужно знать и понимать Закон всемирного тяготения Ньютона и Законы движения Кеплера . Эти два просты, и я уверен, что вы слышали о них, если не изучали их в средней школе. Наконец, если вы хотите, чтобы ваш симулятор был максимально точным, вам следует ознакомиться с проблемой n-Body .

Вы должны начать с простого. Попробуйте создать объект Sun и объект Earth, которые вращаются вокруг него. Это должно дать вам очень твердое начало, и оттуда довольно легко расширяться. Объект планеты будет выглядеть примерно так:

Class Planet {
  float x;
  float y;
  float z; // If you want to work in 3D
  double velocity;
  int mass;
}

Просто помните, что F = MA, а остальное просто скучная математика: P

Answer 3

Это отличный учебник по проблемам N-тела в целом.

http://www.artcompsci.org/#msa

Он написан с использованием Ruby, но довольно легко отображать на другие языки и т. Д. Он охватывает некоторые из общих подходов интеграции; Форвард-Эйлер, чехарда и эрмит.

Answer 4

Возможно, вы захотите взглянуть на Celestia , симулятор свободного пространства. Я считаю, что вы можете использовать его для создания фиктивных солнечных систем, и это с открытым исходным кодом .

Answer 5

В физике это известно как проблема N-тела . Он известен тем, что вы не можете решить это вручную для системы с более чем тремя планетами. К счастью, вы можете получить приблизительные решения с помощью компьютера очень легко.

Хорошую статью о написании этого кода с нуля можно найти здесь .

Тем не менее, я чувствую слово предупреждения здесь важно. Вы можете не получить ожидаемых результатов. Если вы хотите увидеть, как:

1. масса планеты влияет на ее орбитальную скорость вокруг Солнца, круто. Вы увидите это.
2. разные планеты взаимодействуют друг с другом, вы будете огорчены.

Проблема в следующем.

Да, современные астрономы обеспокоены тем, как масса Сатурна меняет орбиту Земли вокруг Солнца. Но это ОЧЕНЬ незначительный эффект. Если вы собираетесь проложить путь планеты вокруг Солнца, вряд ли будет иметь значение, что в Солнечной системе есть другие планеты. Солнце настолько велико, что заглушит всю остальную гравитацию. Единственные исключения из этого:

1. Если ваши планеты имеют очень эллиптические орбиты. Это приведет к тому, что планеты потенциально станут ближе друг к другу, поэтому они будут больше взаимодействовать.
2. Если ваши планеты находятся практически на одинаковом расстоянии от Солнца. Они будут больше взаимодействовать.
3. Если вы сделаете свои планеты настолько комичными, что они будут конкурировать с Солнцем за гравитацию во внешней Солнечной системе.

Чтобы было ясно, да, вы сможете рассчитать некоторые взаимодействия между планетами. Но нет, эти взаимодействия не будут значимыми для невооруженного глаза, если вы создадите реалистичную Солнечную систему.

Попробуйте и узнайте!

Answer 6

Все, что вам нужно реализовать, - это правильное дифференциальное уравнение (закон Кеплера) и использование Рунге-Кутты. (пока это не сработало для меня, но, возможно, есть лучшие методы)

В сети есть множество таких симуляторов.

Вот один простой пример, реализованный в 500 строках кода c. (Монитор Алгоритм намного меньше) http://astro.berkeley.edu/~dperley/programs/ssms.html.

Также проверьте это:
http://en.wikipedia.org/wiki/Kepler_problem
http://en.wikipedia.org/wiki/Two-body_problem
http://en.wikipedia.org/wiki/N-body_problem

Answer 7

Похоже, что это очень сложно и требует глубоких знаний физики, но на самом деле это очень просто, вам нужно знать только 2 формулы и базовое понимание векторов:

Сила притяжения (или гравитационная сила) между планетой 1 и планетой 2 с массой m1 и m2 и расстоянием между ними d: Fg = G * m1 * m2 / d ^ 2; Fg = m * a. G - это константа, найдите ее, подставив случайные значения так, чтобы ускорение «а» не было слишком маленьким и не слишком большим, приблизительно «0,01» или «0,1».

Если у вас есть общая векторная сила, действующая на текущую планету в этот момент времени, вы можете найти мгновенное ускорение a = (общая сила) / (масса текущей планеты). И если у вас есть текущее ускорение, текущая скорость и текущая позиция, вы можете найти новую скорость и новую позицию

Если вы хотите выглядеть реально, вы можете использовать следующий сверхлегкий алгоритм (псевдокод):

int n; // # of planets
Vector2D planetPosition[n]; 
Vector2D planetVelocity[n]; // initially set by (0, 0)
double planetMass[n];

while (true){
    for (int i = 0; i < n; i++){
        Vector2D totalForce = (0, 0); // acting on planet i
        for (int j = 0; j < n; j++){
            if (j == i)
                continue; // force between some planet and itself is 0
            Fg = G * planetMass[i] * planetMass[j] / distance(i, j) ^ 2;
        // Fg is a scalar value representing magnitude of force acting
        // between planet[i] and planet[j]
        // vectorFg is a vector form of force Fg
        // (planetPosition[j] - planetPosition[i]) is a vector value
        // (planetPosition[j]-planetPosition[i])/(planetPosition[j]-plantetPosition[i]).magnitude() is a
        // unit vector with direction from planet[i] to planet[j]
            vectorFg = Fg * (planetPosition[j] - planetPosition[i]) / 
                  (planetPosition[j] - planetPosition[i]).magnitude();
            totalForce += vectorFg;
        }
        Vector2D acceleration = totalForce / planetMass[i];
        planetVelocity[i] += acceleration;
    }

    // it is important to separate two for's, if you want to know why ask in the comments
    for (int i = 0; i < n; i++)
        planetPosition[i] += planetVelocity[i];

    sleep 17 ms;
    draw planets;
}

Answer 8

Алгоритмы для моделирования физики планет.

Вот реализация частей Keppler в моем приложении для Android. Основные части находятся на моем веб-сайте, вы можете скачать весь источник: http://www.barrythomas.co.uk/keppler.html

Это мой метод рисования планеты в «следующей» позиции на орбите. Подумайте о таких шагах, как обход круга, по одному градусу за раз, по кругу, который имеет тот же период, что и планета, которую вы пытаетесь отследить. Вне этого метода я использую глобальный двойной как счетчик шагов - dTime, который содержит несколько градусов поворота.

Ключевые параметры, передаваемые методу: dEccentricty, dScalar (масштабный коэффициент, позволяющий отображать всю орбиту на дисплее), dYear (продолжительность орбиты в земных годах) и ориентация орбиты так, чтобы перигелий находился на правильное место на циферблате, так сказать, dLongPeri - долгота перигелия.

drawPlanet:

public void drawPlanet (double dEccentricity, double dScalar, double dYear, Canvas canvas, Paint paint, 
            String sName, Bitmap bmp, double dLongPeri)
{
        double dE, dr, dv, dSatX, dSatY, dSatXCorrected, dSatYCorrected;
        float fX, fY;
        int iSunXOffset = getWidth() / 2;
        int iSunYOffset = getHeight() / 2;

        // get the value of E from the angle travelled in this 'tick'

        dE = getE (dTime * (1 / dYear), dEccentricity);

        // get r: the length of 'radius' vector

        dr = getRfromE (dE, dEccentricity, dScalar);

        // calculate v - the true anomaly

        dv = 2 * Math.atan (
                Math.sqrt((1 + dEccentricity) / (1 - dEccentricity))
                *
                Math.tan(dE / 2)
                ); 

        // get X and Y coords based on the origin

        dSatX = dr / Math.sin(Math.PI / 2) * Math.sin(dv);
        dSatY = Math.sin((Math.PI / 2) - dv) * (dSatX / Math.sin(dv));

        // now correct for Longitude of Perihelion for this planet

        dSatXCorrected = dSatX * (float)Math.cos (Math.toRadians(dLongPeri)) - 
            dSatY * (float)Math.sin(Math.toRadians(dLongPeri));
        dSatYCorrected = dSatX * (float)Math.sin (Math.toRadians(dLongPeri)) + 
            dSatY * (float)Math.cos(Math.toRadians(dLongPeri));

        // offset the origin to nearer the centre of the display

        fX = (float)dSatXCorrected + (float)iSunXOffset;
        fY = (float)dSatYCorrected + (float)iSunYOffset;

        if (bDrawOrbits)
            {
            // draw the path of the orbit travelled
            paint.setColor(Color.WHITE);
            paint.setStyle(Paint.Style.STROKE);
            paint.setAntiAlias(true);

            // get the size of the rect which encloses the elliptical orbit

            dE = getE (0.0, dEccentricity);
            dr = getRfromE (dE, dEccentricity, dScalar);
            rectOval.bottom = (float)dr;
            dE = getE (180.0, dEccentricity);
            dr = getRfromE (dE, dEccentricity, dScalar);
            rectOval.top = (float)(0 - dr);

            // calculate minor axis from major axis and eccentricity
            // http://www.1728.org/ellipse.htm

            double dMajor = rectOval.bottom - rectOval.top;
            double dMinor = Math.sqrt(1 - (dEccentricity * dEccentricity)) * dMajor;

            rectOval.left = 0 - (float)(dMinor / 2);
            rectOval.right = (float)(dMinor / 2);

            rectOval.left += (float)iSunXOffset;
            rectOval.right += (float)iSunXOffset;
            rectOval.top += (float)iSunYOffset;
            rectOval.bottom += (float)iSunYOffset;

            // now correct for Longitude of Perihelion for this orbit's path

            canvas.save();
                canvas.rotate((float)dLongPeri, (float)iSunXOffset, (float)iSunYOffset);
                canvas.drawOval(rectOval, paint);
            canvas.restore();
            }

        int iBitmapHeight = bmp.getHeight();

        canvas.drawBitmap(bmp, fX - (iBitmapHeight / 2), fY - (iBitmapHeight / 2), null);

        // draw planet label

        myPaint.setColor(Color.WHITE);
        paint.setTextSize(30);
        canvas.drawText(sName, fX+20, fY-20, paint);
}

Приведенный выше метод вызывает два дополнительных метода, которые обеспечивают значения E (средняя аномалия) и r, длины вектора, в конце которого находится планета.

Gete:

public double getE (double dTime, double dEccentricity)
    {
    // we are passed the degree count in degrees (duh) 
    // and the eccentricity value
    // the method returns E

    double dM1, dD, dE0, dE = 0; // return value E = the mean anomaly
    double dM; // local value of M in radians

    dM = Math.toRadians (dTime);

    int iSign = 1;

    if (dM > 0) iSign = 1; else iSign = -1;

    dM = Math.abs(dM) / (2 * Math.PI); // Meeus, p 206, line 110
    dM = (dM - (long)dM) * (2 * Math.PI) * iSign; // line 120
    if (dM < 0)
        dM = dM + (2 * Math.PI); // line 130
    iSign = 1;
    if (dM > Math.PI) iSign = -1; // line 150
    if (dM > Math.PI) dM = 2 * Math.PI - dM; // line 160

    dE0 = Math.PI / 2; // line 170
    dD = Math.PI / 4; // line 170

    for (int i = 0; i < 33; i++) // line 180 
        {
        dM1 = dE0 - dEccentricity * Math.sin(dE0); // line 190
        dE0 = dE0 + dD * Math.signum((float)(dM - dM1));
        dD = dD / 2; 
        }

    dE = dE0 * iSign;

    return dE;
    }

getRfromE:

public double getRfromE (double dE, double dEccentricty, double dScalar)
{
    return Math.min(getWidth(), getHeight()) / 2 * dScalar * (1 - (dEccentricty * Math.cos(dE)));
}

Answer 9

Проверьте nMod , набор инструментов для моделирования n-body, написанный на C ++ и использующий OpenGL. Он имеет довольно густонаселенную модель солнечной системы, которая поставляется с ним, и ее должно быть легко изменить. Кроме того, у него есть довольно хорошая вики о симуляции n-тела в целом. Тот же самый парень, который создал это, также создает новую программу под названием Moody , но, похоже, она не так уж и далеко.

Кроме того, если вы собираетесь проводить симуляции n-тел с несколькими объектами, вам действительно стоит взглянуть на метод быстрого мультиполя (также называемый алгоритмом быстрого мультиполя) Это может уменьшить количество вычислений от O (N ^ 2) до O (N), чтобы действительно ускорить ваше моделирование. Это также один из десяти самых успешных алгоритмов 20-го века , по мнению автора этой статьи.

Answer 10

Dear Friend here is the graphics code that simulate solar system

Kindly refer through it

/*Arpana*/

#include<stdio.h>
#include<graphics.h>
#include<conio.h>
#include<math.h>
#include<dos.h>
void main()
{
int i=0,j=260,k=30,l=150,m=90;
int n=230,o=10,p=280,q=220;
float pi=3.1424,a,b,c,d,e,f,g,h,z;
int gd=DETECT,gm;
initgraph(&gd,&gm,"c:\tc\bgi");
outtextxy(0,10,"SOLAR SYSTEM-Appu");
outtextxy(500,10,"press any key...");
circle(320,240,20);               /* sun */
setfillstyle(1,4);
floodfill(320,240,15);
outtextxy(310,237,"sun");
circle(260,240,8);
setfillstyle(1,2);
floodfill(258,240,15);
floodfill(262,240,15);
outtextxy(240,220,"mercury");
circle(320,300,12);
setfillstyle(1,1);
floodfill(320,298,15);
floodfill(320,302,15);
outtextxy(335,300,"venus");
circle(320,160,10);
setfillstyle(1,5);
floodfill(320,161,15);
floodfill(320,159,15);
outtextxy(332,150, "earth");
circle(453,300,11);
setfillstyle(1,6);
floodfill(445,300,15);
floodfill(448,309,15);
outtextxy(458,280,"mars");
circle(520,240,14);
setfillstyle(1,7);
floodfill(519,240,15);
floodfill(521,240,15);
outtextxy(500,257,"jupiter");
circle(169,122,12);
setfillstyle(1,12);
floodfill(159,125,15);
floodfill(175,125,15);
outtextxy(130,137,"saturn");
circle(320,420,9);
setfillstyle(1,13);
floodfill(320,417,15);
floodfill(320,423,15);
outtextxy(310,400,"urenus");
circle(40,240,9);
setfillstyle(1,10);
floodfill(38,240,15);
floodfill(42,240,15);
outtextxy(25,220,"neptune");
circle(150,420,7);
setfillstyle(1,14);
floodfill(150,419,15);
floodfill(149,422,15);
outtextxy(120,430,"pluto");
getch();
while(!kbhit())             /*animation*/
{
a=(pi/180)*i;
b=(pi/180)*j;
c=(pi/180)*k;
d=(pi/180)*l;
e=(pi/180)*m;
f=(pi/180)*n;
g=(pi/180)*o;
h=(pi/180)*p;
z=(pi/180)*q;
cleardevice();
circle(320,240,20);
setfillstyle(1,4);
floodfill(320,240,15);
outtextxy(310,237,"sun");

circle(320+60*sin(a),240-35*cos(a),8);
setfillstyle(1,2);
pieslice(320+60*sin(a),240-35*cos(a),0,360,8);
circle(320+100*sin(b),240-60*cos(b),12);
setfillstyle(1,1);
pieslice(320+100*sin(b),240-60*cos(b),0,360,12);
circle(320+130*sin(c),240-80*cos(c),10);
setfillstyle(1,5);
pieslice(320+130*sin(c),240-80*cos(c),0,360,10);
circle(320+170*sin(d),240-100*cos(d),11);
setfillstyle(1,6);
pieslice(320+170*sin(d),240-100*cos(d),0,360,11);
circle(320+200*sin(e),240-130*cos(e),14);
setfillstyle(1,7);
pieslice(320+200*sin(e),240-130*cos(e),0,360,14);
circle(320+230*sin(f),240-155*cos(f),12);
setfillstyle(1,12);
pieslice(320+230*sin(f),240-155*cos(f),0,360,12);
circle(320+260*sin(g),240-180*cos(g),9);
setfillstyle(1,13);
pieslice(320+260*sin(g),240-180*cos(g),0,360,9);
circle(320+280*sin(h),240-200*cos(h),9);
setfillstyle(1,10);
pieslice(320+280*sin(h),240-200*cos(h),0,360,9);
circle(320+300*sin(z),240-220*cos(z),7);
setfillstyle(1,14);
pieslice(320+300*sin(z),240-220*cos(z),0,360,7);
delay(20);
i++;
j++;
k++;
l++;
m++;
n++;
o++;
p++;
q+=2;
}
getch();
}