simulation - ¿Es posible hacer una simulación realista del sistema solar de n cuerpos en materia de tamaño y masa?

physics game-physics

Question

Nota importante: esta pregunta no tiene ninguna relación con "PhysX", que es un sistema de física de juegos de computadora (útil para la física en juegos de arcade como juegos de pelota, etc.); PhysX es un sistema integrado en Unity3D y otros motores de juegos; PhysX es totalmente irrelevante aquí.

//////////////////// ACTUALIZACIÓN (lea la parte inferior primero) /////////////////////

He estado registrando los valores y buscando dónde está el problema exacto, y creo que lo encontré. Tengo algo como esto en mi código

Velocity += Acceleration * Time.deltaTime;
position += Velocity * Time.deltaTime;

La aceleración es algo así como 0,0000000000000009 .. en este momento. A medida que el código fluye, la velocidad aumenta como debería ser, no hay problema con el flotador. Pero al principio, la posición inicial de la tierra es (0,0,23500f). Puedes ver esto en la tabla que di al final.

Bueno, ahora cuando agrego speed * timedelta (que es algo así como 0,00000000000000005 en este punto) a la posición que es 23500, básicamente no lo agrega. la posición sigue siendo (0, 0, 23500) no algo así como (0,0, 23500.00000000000005), por lo tanto, la tierra no se mueve, por lo tanto, la aceleración no cambia.

Si configuré la posición inicial de la tierra en 0,0,0 y aún así, configuré la aceleración en 0,0000000000000000009 para suponer que su posición es (0,0,23500). Luego "AGREGA" la velocidad * timedelta. Se convierte en algo así como (0,0,000000000000000000005) y sigue aumentando. Cuando el flotante es 0, no hay ningún problema al agregar un valor tan pequeño. Pero si el flotante es algo así como 23500, entonces no suma los valores pequeños.

No sé si es exactamente el problema de la unidad o el flotador de C #.

Y es por eso que no puedo hacer que funcione con valores pequeños. Si puedo superar esto, mi problema se resolverá.

//////////////////////////////////////////////////////////////////////////////////////////////////////////// /////////////////////////////

He estado desarrollando la física de los cuerpos n para simular nuestro sistema solar, por lo que he estado recopilando datos para que sea lo más realista posible. Pero hay un problema con el tamaño de los datos. Busqué en Internet y no pude encontrar una sola explicación de cómo la gente supera esto. (Si es así) Así que estoy intentando mi tiro aquí.

Entonces, para mantener fija la relación de distancia, radio y "masa" entre planetas, creé un archivo de Excel para calcular todos los datos. (¿Por qué diablos alguien pondría "cuál sería la masa de la Tierra si tuviera" ese "gráfico de radio" en Internet?) Le daré el ss como archivo adjunto. Básicamente "normaliza" o, en otras palabras, "escala" cada propiedad de un planeta a una referencia dada. En este caso, tomé la referencia como "radio de la tierra".

Trabajo en unidad, y sabes, no puedes trabajar con valores "demasiado grandes" o "demasiado pequeños" en unidad. Así que tuve que reducir el sistema solar, "¡mucho!"

Así que uso la ley de Newton de la gravitación universal, que es F = GMm / r ^ 2, para hacerlo simple, estoy calculando directamente a = GM / r ^ 2, para un cuerpo dado a partir de todos los demás cuerpos.

Entonces, el valor real de la aceleración gravitacional de la Tierra "hacia el sol" es de aproximadamente 0,000006 km / s ^ 2, que es un valor increíblemente pequeño para trabajar en unidad, pero podría funcionar. Sin embargo, para obtener este valor, 1 necesito establecer el radio (escala) de la Tierra en 6371 unidades, y el sol en una escala de 696,342 !, lo cual es DEMASIADO grande para hacerlo unitario.

Entonces dije, que el radio de la Tierra sea 1, en unidades de unidad. Entonces, cuando el radio cambia, todo cambia, la masa, la distancia ... Mantuve la densidad del planeta y calculo la masa del nuevo volumen con el nuevo radio. Todos los cálculos están en el archivo adjunto.

Entonces, cuando tomo el radio de la Tierra como 1, la aceleración gravitacional hacia el sol se convierte en algo así como 0,0000000000009, que es ridículamente pequeño. Y, por supuesto, Unity no funciona con ese valor.

Entonces, si aumento el radio de la Tierra, entonces la masa y el radio del Sol se vuelven ridículamente grandes y, de nuevo, no puedo trabajar con él.

No sé cómo otras personas arreglaron esto, qué hicieron para superar este problema, pero como veo desde aquí, parece imposible hacer una simulación realista del sistema solar con n cuerpos. (en unidad al menos)

Así que necesito tener 10 representantes para publicar imágenes -_-, en su lugar, daré un enlace. http://berkaydursun.com/solar_system_simulator/data.PNG También un directorio arriba es la simulación experimental del sistema solar en funcionamiento con cálculos de n cuerpos pero con valores NO REALISTAS. Funciona bastante bien, e incluso se ve de alguna manera casi real, pero no, no tiene las proporciones correctas ^^ Puede probarlo aquí si lo desea http://berkaydursun.com/solar_system_simulator/

Editar: WoW casi comencé cada párrafo con "So" ^^

Answer 1

Reducir las cosas no necesariamente ayudará, como has descubierto. Aquí hay una buena lectura sobre cosas a considerar cuando se usan números de coma flotante: http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html

Básicamente, hacer una simulación a partir de los primeros principios (las leyes de Newton) es malo para la precisión numérica, porque no imbuyes a los métodos numéricos con una idea de la escala de los efectos importantes, por lo que terminas arrojando un montón de efectos diferentes a diferentes escalas juntos y el resultado es baja precisión.

Por lo general, cosas como efemérides para planetas, satélites, etc. no comienzan con la ley de Newton, comienzan suponiendo que las órbitas son keplerianas, y luego pequeñas correcciones perturbativas.

Aquí hay un algoritmo que calcula la posición de los planetas (semiimperialmente). http://ssd.jpl.nasa.gov/txt/aprx_pos_planets.pdf

Si quieres hacer una simulación de N-body, parece que necesitarás más precisión. Si la unidad le impide usar doble precisión, entonces sugiero hacer los cálculos en C # simple, luego convertir a precisión simple cuando el trabajo esté terminado.

Answer 2

También programé la simulación del sistema sol, así que aquí están las ideas sobre las minas:

1. representación

   Yo uso OpenGL con escala 1: 1 . Todas las unidades están en SI así que [m, s, kg, ...] . El problema comienza con Z-buffer . El ancho habitual del búfer Z es de 16/24/32 bit que no se acerca a lo que necesita. Estoy renderizando desde 0.1m hasta 1000 AU entonces, ¿cómo superar esto?

   Lo logré renderizando con 3 frustrums a la vez combinando Z-sorting y Z-buffering (Z-sort es necesario debido a anillos transparentes ... y otros efectos). Entonces, primero renderizo las partes más distantes hasta zfar=1000AU . El domo del cielo se proyecta a z=750AU distancia z=750AU luego despeja el búfer Z y renderiza objetos hasta zfar=0.1AU . Luego, borre el Z-buffer nuevamente y renderice objetos cercanos hasta zfar=100000 m .

   Para realizar este trabajo, debe tener una matriz de proyección lo más precisa posible. El gluPerspective tiene cotangens poco precisos, por lo que necesita reparar los elementos afectados (me lleva mucho tiempo detectarlo). Z near valor Z near depende del ancho de bits del búfer Z. Cuando se codifica correctamente, esto funciona bien incluso con zoom 10000x . Utilizo este programa como navegador / buscador de objetos para el telescopio de minas :) en tiempo real desde la vista de inicio de la mina. Combino estrellas 3D, cuerpos astronómicos, barcos, tierra real (a través de DTM y textura satelital). Es capaz incluso de salida anaglifo rojo-cian :). Puede renderizar desde la superficie, la atmósfera, el espacio ... (no solo bloqueado a la Tierra). No se utiliza ninguna otra biblioteca de terceros, entonces OpenGL. Así es como se ve:

   Como puede ver, funciona bien en cualquier altitud o zoom, la atmósfera se hace como este sombreador de dispersión de atmósfera

2. simulación

   No estoy usando la simulación de gravedad de cuerpo n porque para eso necesitas una gran cantidad de datos que es muy difícil de obtener (y casi imposible en la precisión deseada). Los cálculos deben hacerse con mucha precisión .

   En su lugar, uso la ecuación de Kepler, así que mira esto:

   • Resolviendo la ecuación de Kepler
   • Implementación de C ++ .

   Si todavía quieres usar el modelo de gravedad, usa los horizontes JPL de la NASA . Creo que también tienen códigos fuente en C / C ++, pero están usando un marco de referencia incompatible con los mapas de minas, por lo que es inutilizable para mí.

   En general, la ecuación de Kepler tiene un error mayor, pero no aumenta demasiado con el tiempo. El modelo de gravedad es más preciso, pero su error aumenta con el tiempo y debe actualizar los datos del cuerpo del astro continuamente para que funcione ...

[edit1] precisión de integración

su implementación actual es así:

// object variables
double acc[3],vel[3],pos[3];
// timer iteration
double dt=timer.interval;
for (int i=0;i<3;i++)
{
vel[i]+=acc[i]*dt;
pos[i]+=vel[i]*dt;
}

El problema es que cuando agregas un valor muy pequeño y muy grande, se cambian al mismo exponente antes de la suma, lo que redondeará los datos significativos. Para evitar esto, simplemente cámbialo a esto:

// object variables
double vel0[3],pos0[3]; // low
double vel1[3],pos1[3]; // high
double acc [3],vel [3],pos [3]; // full
// timer iteration
double dt =timer.interval;
double max=10.0; // precision range constant
for (int i=0;i<3;i++)
{
vel0[i]+=acc[i]*dt; if (fabs(vel0[i]>=max)) { vel1[i]+=vel0[i]; vel0[i]=0.0; } vel[i]=vel0[i]+vel1[i];
pos0[i]+=vel[i]*dt; if (fabs(pos0[i]>=max)) { pos1[i]+=pos0[i]; pos0[i]=0.0; } pos[i]=pos0[i]+pos1[i];
}

Ahora xxx0 está integrado hasta el max y todo se agrega a xxx1

El redondeo sigue ahí, pero ya no es acumulativo. xxx0+xxx1 seleccionar el valor max que la integración en sí es segura y también la adición xxx0+xxx1 debe ser segura. Entonces, si los números son demasiado diferentes para una división, divídalos dos veces o más ...

• como: xxx0+=yyy*dt; if (fabs(xxx0>max0))... if (fabs(xxx1>max1))... xxx0+=yyy*dt; if (fabs(xxx0>max0))... if (fabs(xxx1>max1))...

[Editar2] Estrellas

• Los colores de las estrellas La mejor visualización de estrellas que he visto
• Índice de color Star BV a color RGB aparente, todos los catálogos de estrellas utilizan el índice BV
• Usando los catálogos de Stelar también hay enlaces de referencias cruzadas de nombres de estrellas
• Skybox: combina diferentes datos de estrellas

[Edit3] Mejorando aún más la precisión de integración de Newton D''ALembert

El problema básico con la integración iterativa es que aplicar una aceleración basada en la gravedad basada en la posición actual del cuerpo dará como resultado órbitas más grandes porque durante el paso de integración dt la posición cambia un poco, lo que no se tiene en cuenta en la integración ingenua. Para remediar esto, eche un vistazo a esta imagen:

Supongamos que nuestro cuerpo está en órbita circular y en la posición de 0 grados. En lugar de usar la dirección de aceleración basada en la posición actual, usé la posición después de 0.5*dt . Esto aumenta el bit pequeño de aceleración que resulta en una precisión mucho más alta (correspondencia con las órbitas de Kepler). Con este ajuste pude convertir con éxito de la órbita de Kepler a Newton D''Alembert para un sistema de 2 cuerpos. (hacer esto para n-body es el siguiente paso). La correlación aproximada con datos reales de nuestro sistema solar solo es posible para un sistema de 2 cuerpos no afectado por los efectos de las mareas y / o lunas. Para construir datos ficticios propios, puede usar la órbita circular de Kepler y la fuerza contripedal que iguala la gravedad:

G = 6.67384e-11;
v = sqrt(G*M/a); // orbital speed
T = sqrt((4.0*M_PI*M_PI*a*a*a)/(G*(m+M))); // orbital period

donde a es el radio circular de la órbita m es la masa corporal, M es la masa focal del cuerpo (sol). Para mantener la precisión en una tolerancia aceptable (para mí), el paso de integración dt debería ser:

dt = 0.000001*T

Entonces, para poner un nuevo cuerpo a prueba, simplemente póngalo en:

pos = (a,0,0)
vel = (0,sqrt(G*M/a),0)

Mientras el cuerpo focal principal (Sol) está en:

pos = (0,0,0)
vel = (0,0,0)

Esto colocará su cuerpo en órbita circular para que pueda comparar Kepler versus Newton D''Alembert para evaluar la precisión de su simulación.