Un globo que asciende

Procedimiento numérico

Resolver la ecuación diferencial de primer orden

$\begin{array}{l} \frac{d y}{d t} = \sqrt{\frac{20}{3} g r_{0} {λ - \frac{M}{M_{A} n_{0}} \frac{1}{λ^{2}} {(1 - \frac{y}{y_{0}})}^{1 - η}}} \\ {(1 - \frac{y}{y_{0}})}^{η - 1} λ^{7} + b λ^{6} - a λ^{4} - b = 0 a = \frac{n}{n_{0}} b = \frac{16 π k r_{0}^{2}}{3 n_{0}} η = \frac{M_{A} g y_{0}}{R T_{0}} \end{array}$

por procedimientos numéricos con las siguientes condiciones iniciales t=0, y==0, y a la vez, la raíz de una ecuación por el procedimiento del punto medio.

Cuando el globo llega a la altura máxima, su velocidad es nula. En ese momento el argumento de la raiz cuadrada es cero. Un paso después el argumento es negativo y la función f de la clase Sistema lanza (thows) un excepción (un objeto de la clase NegativoExcepcion). La excepción es capturada en el bloque try..catch de la función resolver de la clase base RungeKutta. Esta reenvía (throws) el objeto ex de la clase NegativoExcepcion. La excepción es capturada en el bloque try..catch de la función mover de la clase MiCanvas. En dicho bloque catch se detiene el movimiento del globo.

public class Estado {
    double t;
    double x;
    double v;
    public Estado(double t, double x, double v) {
        this.t=t;
        this.x=x;
        this.v=v;
    }
}
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)throws NegativoExcepcion{
        double k1, k2, k3, k4;
//condiciones iniciales
        double x=e.x;
	    double t=e.t;
       try{
        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);
        }catch(NegativoExcepcion ex){
            throw ex;
        }

        x+=(k1+2*k2+2*k3+k4)/6;
        t+=h;
//cambia el estado de la partícula
	    e.x=x;
	    e.t=t;
    }
    abstract public double f(double x, double t)throws NegativoExcepcion;
}
public class Sistema extends RungeKutta{
    final double y0=49000;
    final double n0=12.5;
    final double r0=0.42;
    Funcion f;
    double masa;
    double eta;
    public Sistema(double eta, double masa, double moles, double paso){
        super(paso);
        this.masa=masa;
        this.eta=eta;
        f=new Funcion(eta, moles);
    }
    public double f(double y, double t) throws NegativoExcepcion{
      f.setAltura(y);
      double lambda=f.puntoMedio(1.0, 7.0);
      double temp=(20*9.8*r0/3)*(lambda-masa*Math.pow((1-y/y0), (1-eta))
/(0.0289*n0*lambda*lambda));
      if(temp<0) throw new NegativoExcepcion("Raíz cuadrada de un número negativo");
      return (Math.sqrt(temp));
    }

}

class NegativoExcepcion extends Exception {
  public NegativoExcepcion(String s) {
         super(s);
  }
}
public class Funcion extends Ecuacion{
  double a;
  final double n0=12.5;
  final double r0=0.42;
  double b;
  final double k=0.46235;
  final double y0=49000;
  double eta;
  double y=0.0;
  public Funcion(double eta, double n) {
    this.b=16*Math.PI*k*r0*r0/(3*n0);
    this.eta=eta;
    this.a=n/n0;
  }
  void setAltura(double y){
    this.y=y;
  }
  public double f(double x){
        double z=Math.pow((1-y/y0), (eta-1))
*x*x*x*x*x*x*x+b*x*x*x*x*x*x-a*x*x*x*x-b;
        return z;
    }
}
public class MiCanvas extends Canvas {
    final double eta;
     final double y0=49000;
     final double r0=0.42;
     final double n0=12.5;
     double radio;
     double masa;
//objeto sistema
   Sistema sistema=null;
   Estado estado=new Estado(0.0, 0.0, 0.0);
   Funcion f;
 
  void mover() {
     if(sistema!=null){
        try{
            sistema.resolver(estado);
            estado.v=sistema.f(estado.x, 0.0);
            f.setAltura(estado.x);
            radio=f.puntoMedio(1.0, 7.0)*r0*100;
        }catch(NegativoExcepcion ex){
            //el globo se detiene
            estado.x=yMax;
            estado.v=0.0;
        }
     }
  }