C++ 4th Order Runge-Kutta Method

2010/11/12
(この記事の文字数: 948)

微分方程式を解くのに実用性の高い４次のルンゲクッタ(Runge-Kutta)法をC++で書いてみました。

Cで使える高階関数を利用して微分関数を指定できるようにしたので、それなりに汎用性はあると思います。

I created a program for 4th order Runge-Kutta method in C++. I used higher-order function in this program, so it is versatile.


void RungeKutta4th(double* dest, double st, double end,

     double step, double initValue, 

     double(*dfunc)(double t, double y))

{

        double y, t, k1, k2, k3, k4;

        uint itr_i, length;

        length = (uint)((end - st) / step);

        t = st;

        y = initValue;

        for(itr_i=0; itr_i < length; itr_i++)

        {

                k1 = step*dfunc(t, y); 

                k2 = step*dfunc(t+step/2.0, y+k1/2.0);

                k3 = step*dfunc(t+step/2.0, y+k2/2.0);

                k4 = step*dfunc(t+step, y+k3);



                y += (k1 + 2*(k2 + k3) + k4) / 6.0;



                dest[itr_i] = y;

                t += step;

        }

}

Example


#include <stdio.h>

#include <stdlib.h>



double dfunc(double t, double y){

        return 3*t*t;

}



void main(){

        double st, end, step, initValue;

        double *dest;

        unsigned int i,len;

        st = 0.0;

        end = 1.0;

        step = 0.1;

        len = (end-st)/step;

        dest = new double[len];

        initValue = 0.0;

        RungeKutta4th(dest,st,end,step,

                        initValue,dfunc);

        for(i=0; i < len; i++){

                printf("%f\n",dest[i]);

        }

        delete[] dest;

}