Ранг Кутта 4-го порядка для решения системы дифференциальных уравнений - PullRequest
Новая
       24

Ранг Кутта 4-го порядка для решения системы дифференциальных уравнений

0 голосов
/ 30 марта 2011 
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

Я должен решить эту систему уравнений в Фортране.Я решил проблему в Matlab, но я не знаю фортрановское программирование, так что, ребята, если кто-то может мне помочь или у кого-то есть код на фортране, помогите мне, пожалуйста, пожалуйста, спасибо большое

1 Ответ

3 голосов
/ 30 марта 2011

Попробуйте это с интеграцией Эйлера. Сначала сделай что-нибудь простое. У вас есть одно преимущество: вы решили это один раз, так что вы знаете, как выглядит ответ, когда получите его.

Поскольку модераторы настаивают на том, что это ответ низкого качества из-за короткой длины, я приведу рабочий ответ на Java, который должен вызвать у вас некоторые мысли. Я использовал математическую библиотеку Apache Commons; у него есть несколько различных схем интеграции ODE, включая Euler и Runge Kutta.

Я запустил это на компьютере с Windows 7, используя JDK 8. Вы можете переключаться между 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. Как видите, ответы гладкие и хорошо себя ведут. У Эйлера нет проблем с такими системами уравнений.

enter image description here

...