Soy un estudiante de doctorado que intenta reproducir resultados similares a los encontrados en Mizuno, Nakazawa y Hayashi, 1978 . La idea es integrar el siguiente conjunto de ecuaciones diferenciales de primer orden
$\frac{d \rho}{dr} = - \left(\frac{G\mu}{k_{b} T}\right)\frac{M_{r} \rho}{r^{2}}$
$\frac{d M_{r}}{dr} = - 4\pi r^{2}\rho$
Dadas las condiciones iniciales $ M_{r}\left(r=r_{p}\right) = M_{0}$ , $ \rho\left(r=_{p}\right) = \rho_{0}$ . Mi dominio es la esfera de Hill de un planeta inmerso en un disco protoplanetario, desde el radio del planeta ( $r_{p}$ ) hasta el límite de la Esfera de la Colina. $G$ es la constante de gravitación, $k_{b}$ es la constante de Boltzmann, $T$ es la temperatura y $\mu$ es la masa de la partícula media. En particular, $M_{0}$ equivale a la masa de un planeta tipo Júpiter. Exploro el otro parámetro $\rho_{0}$ la densidad en el límite superior del planeta.
He explorado varios esquemas numéricos que son bastante sencillos para mí y probablemente bien conocidos por todos ustedes; Euler explícito, Euler implícito y Runge Kutta de 4º orden .
Configuración $M_{0}$ a su valor adecuado parece desencadenar una inestabilidad haciendo que los coeficientes en los esquemas anteriores diverjan en los primeros pasos . Aumentar la discretización del dominio parecía una buena solución para rodear este problema, pero el nivel de discretización tiene que ser más alto de lo que mis capacidades computacionales son capaces de soportar. El Euler implícito, por otro lado, parece ser estable pero produce un resultado no físico (oscilante). Si $M_{0}= 0$ las soluciones están bien pero no representan lo que quiero modelar
Siento que me falta algo al tratar este sistema, y no sé qué otro enfoque puedo adoptar para abordar este problema. El documento anterior no menciona nada relacionado con los esquemas numéricos utilizados. ¿Tienen alguna sugerencia?
Muchas gracias.
