## Sistema de ecuaciones diferenciales de segundo orden

$\begin{array}{l}\frac{{d}^{2}x}{d{t}^{2}}=-4{\pi }^{2}\left(\frac{x}{{\left({x}^{2}+{y}^{2}\right)}^{3/2}}+\frac{\alpha \left(x-d\right)}{{\left({\left(x-d\right)}^{2}+{y}^{2}\right)}^{3/2}}\right)\\ \frac{{d}^{2}y}{d{t}^{2}}=-4{\pi }^{2}\left(\frac{y}{{\left({x}^{2}+{y}^{2}\right)}^{3/2}}+\frac{\alpha y}{{\left({\left(x-d\right)}^{2}+{y}^{2}\right)}^{3/2}}\right)\end{array}$

Definimos una nueva energía e por unidad de masa en este sistema de unidades

$e=\frac{1}{2}\left({v}_{x}^{2}+{v}_{y}^{2}\right)-4{\pi }^{2}\left(\frac{1}{\sqrt{{x}^{2}+{y}^{2}}}+\frac{\alpha }{\sqrt{{\left(x-d\right)}^{2}+{y}^{2}}}\right)$

Con las condiciones iniciales siguientes: en el instante t=0, x=x0, y=y0, vx=v0x, vy=v0y.

 public class Sistema extends RungeKutta{ double alfa; double d; public Sistema(double alfa, double d, double h){ super(h); this.alfa=alfa; this.d=d; } public double f(double x, double y, double vx, double vy, double t){ double r1=Math.sqrt(x*x+y*y); double r2=Math.sqrt((x-d)*(x-d)+y*y); double z=-4*Math.PI*Math.PI*(x/(r1*r1*r1)+alfa*(x-d)/(r2*r2*r2)); return z; } public double g(double x, double y, double vx, double vy, double t){ double r1=Math.sqrt(x*x+y*y); double r2=Math.sqrt((x-d)*(x-d)+y*y); double z=-4*Math.PI*Math.PI*y*(1.0/(r1*r1*r1)+alfa/(r2*r2*r2)); return z; } public double energia(double x, double y, double vx, double vy){ double r1=Math.sqrt(x*x+y*y); double r2=Math.sqrt((x-d)*(x-d)+y*y); double z=(vx*vx+vy*vy)/2-4*Math.PI*Math.PI*(1.0/r1+alfa/r2); return z; } } //Objeto de la clase Sistema void setNuevo(double alfa, double d, double x0, double y0, double V0x, double V0y){ this.d=d; estado=new Estado(0.0, x0, y0, V0x, V0y); sis=new Sistema(alfa, d, calculaPaso()); eInicial=oscilador.energia(x0, y0, V0x, V0y); error=0.0; } double calculaPaso(){ double r1=Math.sqrt(estado.x*estado.x+estado.y*estado.y); double r2=Math.sqrt((estado.x-d)*(estado.x-d)+estado.y*estado.y); double distMed=Math.sqrt(r1*r2); double dt=0.0001*Math.pow((distMed/(d*0.1)), 1.5); return dt; } void mover(){ sis.resolver(estado); sis.setPaso(calculaPaso()); double energia=oscilador.energia(estado.x, estado.y, estado.vx, estado.vy); error=Math.abs((energia-eInicial)*100/eInicial); if(error>1.0){ //se para; } }