Siguiente Anterior

Ecuaciones diferenciales de primer orden

Descripción

Sea una ecuación diferencial de primer orden, con la condición inicial

x(t0)=x0

Se elige una anchura de paso h, y se calculan cuatro números k1, k2, k3, k4 de acuerdo con el procedimiento esquematizado en la tabla adjunta. Según el procedimiento ordinario de Runge-Kutta, a partir del valor de x en el instante t se determina el valor de x en el instante t+h mediante la fórmula que figura en la última fila de dicha tabla.

dx dt =f(x,t)
k 1 =h·f(x,t) k 2 =h·f( x+ 1 2 k 1 ,t+ 1 2 h ) k 3 =h·f( x+ 1 2 k 2 ,t+ 1 2 h ) k 4 =h·f( x+ k 3 ,t+h )
x(t+h)=x(t)+ 1 6 ( k 1 +2 k 2 +2 k 3 + k 4 )

La jerarquía de clases

Definimos una clase base abstracta denomina RungeKutta, cuya función miembro resolver nos calcula la solución aproximada de la ecuación diferencial x en el instante t+h, cuando le pasamos el estado x del sistema en el instante t. La función devuelve el estado del sistema que puede ser la posición de un móvil, la intensidad de una corriente eléctrica en un circuito, etc.

public abstract class RungeKutta {
      double h;
  public RungeKutta(double h) {
     this.h=h;
  }
  public void setPaso(double h){
        this.h=h;
  }
  public void resolver(Estado e){
        double k1, k2, k3, k4;
     //estado inicial en el instante t
         double x=e.x;
	       double t=e.t;

            k1=h*f(x, t);
            k2=h*f(x+k1/2, t+h/2);
            k3=h*f(x+k2/2, t+h/2);
            k4=h*f(x+k3, t+h);
            x+=(k1+2*k2+2*k3+k4)/6;
            t+=h;
     // estado en el instante t+h
	    e.x=x;
	    e.t=t;   
    }
    abstract public double f(double x, double t);
}

La función resolver recibe el estado inicial e del sistema. Ya que los objetos que se pasan a una función se pueden modificar en el curso de su llamada. La función resolver modifica dicho objeto, devolviendo el estado final cuando se le pasa el estado inicial.

La definición de la clase Estado cuyos objetos guardan el estado x del sistema en el instante t

public class Estado {
   public double t;
   public double x;
   public Estado(double t, double x) {
      this.t=t;
      this.x=x;
  }
}

Para obtener el estado e (objeto de la clase Estado) accedemos a sus miembros dato públicos t, x

	System.out.println("t " +e.t);
	System.out.println("x " +e.x);

Para crear el estado e a partir de los valores de t, x llamamos al constructor de la clase Estado.

        double t=0;
        double x=5000;

	    Estado e=new Estado(t, x);

En la clase Sistema derivada de RungeKutta describimos el modelo físico particular, redefiniendo la función f(x, t). Consideremos la ecuación diferencial que describe la desintegración de una sustancia radioactiva en otra estable, o su fenómeno análogo la descarga de un tubo-capilar.

dx dt =− a x  x= x 0 exp⁡(− a t)

Donde a es la constante de desintegración radioactiva. A la izquierda tenemos la ecuación diferencial y a la derecha su solución analítica.

La clase Sistema define la función f(x,t), tomando el valor de la constante de desintegración a igual a 0.1

public class Sistema extends RungeKutta{
  public Sistema(double h) {
     super(h);
  }
   public double f(double x, double t){
        return (-0.1*x);
    }
}

La llamada a la función resolver se efectuará desde un objeto de la clase derivada Sistema. Se nos pedirá los siguientes datos: el estado inicial, es decir, el número de núcleos x0 en el instante inicial t=0, el instante final t en el que queremos calcular el nuevo estado del sistema x, y el paso h para resolver numéricamente la ecuación diferencial. 

public class Aplicacion {
  public static void main(String[] args) {
        double h=0.5;       //paso
        Estado e=new Estado(0.0, 5000.0);
        Sistema sis=new Sistema(h);
        double t=20.0; //instante final
//solución aproximada
        for(int i=0; i<t/h; i++){
            sis.resolver(e);
        }
        System.out.println("valor aproximado de x "+(int)e.x);
// valor exacto
        double x=(int)(5000.0*Math.exp(-0.1*t));
        System.out.println("valor exacto de x "+(int)x);
  }
}

Comparamos el valor exacto y el valor aproximado, tomando como paso h el valor 0.5, el número inicial de núcleos, 5000, y el instante t, 20. Se obtiene para el resultado aproximado 676 y para el resultado exacto 676.

SiguienteAnterior