system - runge kutta 四阶求解微分方程组

Question

dT/dt=(1.344-1.025T)/h (1)

dh/dt=0.025-(3.5*10^-4)*sqrt(h) (2)

h(0)=1

T(0)=1

我必须在 fortran 中求解这个方程组。我在matlab中解决了这个问题，但我不知道fortran编程，所以如果有人可以帮助我或者有人有fortran代码来帮助我，请非常感谢

Answer 1

尝试使用欧拉积分。先做点简单的。你有一个优势：你已经解决了一次，所以当你得到答案时你知道答案是什么样子的。

由于版主坚持认为这是一个低质量的答案，因为篇幅很短，我将提供一个 Java 中的有效答案，它应该会激发您的一些想法。我使用了 Apache Commons 数学库；它有几种不同的 ODE 集成方案，包括 Euler 和 Runge Kutta。

我在使用 JDK 8 的 Windows 7 机器上运行它。您可以使用命令行在 Euler 和 Runge-Kutta 之间切换：

package math.ode;

import org.apache.commons.math3.exception.DimensionMismatchException;
import org.apache.commons.math3.exception.MaxCountExceededException;
import org.apache.commons.math3.ode.FirstOrderDifferentialEquations;
import org.apache.commons.math3.ode.FirstOrderIntegrator;
import org.apache.commons.math3.ode.nonstiff.ClassicalRungeKuttaIntegrator;
import org.apache.commons.math3.ode.nonstiff.EulerIntegrator;

/**
 * IntegrationExample solves coupled ODEs using Euler and Runge Kutta
 * Created by Michael
 * Creation date 12/20/2015.
 * @link https://stackoverflow.com/questions/20065521/dependencies-for-jama-in-maven
 */
public class IntegrationExample {

    public static final double DEFAULT_STEP_SIZE = 0.001;
    private static final double DEFAULT_MAX_TIME = 2.0;

    public static void main(String[] args) {
        // Problem set up
        double step = (args.length > 0) ? Double.valueOf(args[0]) : DEFAULT_STEP_SIZE;
        double maxTime = (args.length > 1) ? Double.valueOf(args[1]) : DEFAULT_MAX_TIME;
        String integratorName = (args.length > 2) ? args[2] : "euler";
        // Choose different integration schemes here.
        FirstOrderIntegrator firstOrderIntegrator = getFirstOrderIntegrator(step, integratorName);
        // Equations to solve here; see class below
        FirstOrderDifferentialEquations odes = new CoupledOdes();
        double [] y = ((CoupledOdes) odes).getInitialConditions();
        double t = 0.0;
        int i = 0;
        while (t <= maxTime) {
            System.out.println(String.format("%5d %10.6f %10.6f %10.6f", i, t, y[0], y[1]));
            firstOrderIntegrator.integrate(odes, t, y, t+step, y);
            t += step;
            ++i;
        }
    }

    private static FirstOrderIntegrator getFirstOrderIntegrator(double step, String integratorName) {
        FirstOrderIntegrator firstOrderIntegrator;
        if ("runge-kutta".equalsIgnoreCase(integratorName)) {
            firstOrderIntegrator = new ClassicalRungeKuttaIntegrator(step);
        } else {
            firstOrderIntegrator = new EulerIntegrator(step);
        }
        return firstOrderIntegrator;
    }
}

class CoupledOdes implements FirstOrderDifferentialEquations {

    public double [] getInitialConditions() {
        return new double [] { 1.0, 1.0 };
    }

    @Override
    public int getDimension() {
        return 2;
    }

    @Override
    public void computeDerivatives(double t, double[] y, double[] yDot) throws MaxCountExceededException, DimensionMismatchException {
        yDot[0] = (1.344-1.025*y[0])/y[1];
        yDot[1] = 0.025-3.5e-4*Math.sqrt(y[1]);
    }
}

你没有说你需要多长时间才能及时整合，所以我假设 2.0 作为最长时间。您也可以在命令行上更改它。

这是 Excel 中结果与时间的关系图。如您所见，响应流畅且表现良好。欧拉对这样的方程组没有任何问题。