四階龍格庫塔法解振蕩方程c源碼及畫圖

1、vc+6.0 實(shí)現(xiàn)微分方程求解，程序如下：#include <iostream.h>#include <math.h>double f1(double t, double x, double y, double z) return -x + 6 * y + 2 * z;double f2(double t, double x, double y, double z) return -y + 3 * z + 2 * sin(t);double f3(double t, double x, double y, double z) return -z + pow(t,2) *

2、 exp(-t) + cos(t);double rungekutta(double x0, double y0, double z0, double a, double b, int n, double x, double y, double z, double (*f1)(double,double,double,double), double (*f2)(double,double,double,double), double (*f3)(double,double,double,double) double k1, k2, k3, k4 ; double m1, m2, m3, m4

3、; double n1, n2, n3, n4 ; double h; double t;int count;count =0;h = (b - a) / n; x0 = x0; y0 = y0; z0 = z0; for(int i = 1;i<=n;i+) t = a + (i - 1) * h; k1 = f1(t, xi - 1, yi - 1, zi - 1); m1 = f2(t, xi - 1, yi - 1, zi - 1); n1 = f3(t, xi - 1, yi - 1, zi - 1); k2 = f1(t+h/2.0, xi - 1 + k1 * h / 2.

4、0, yi - 1 + m1 * h / 2.0, zi - 1 + n1 * h / 2.0); m2 = f2(t+h/2.0, xi - 1 + k1 * h / 2.0, yi - 1 + m1 * h / 2.0, zi - 1 + n1 * h / 2.0); n2 = f3(t+h/2.0, xi - 1 + k1 * h / 2.0, yi - 1 + m1 * h / 2.0, zi - 1 + n1 * h / 2.0); k3 = f1(t+h/2.0, xi - 1 + k2 * h / 2.0, yi - 1 + m2 * h / 2.0, zi - 1 + n2 *

5、 h / 2.0); m3 = f2(t+h/2.0, xi - 1 + k2 * h / 2.0, yi - 1 + m2 * h / 2.0, zi - 1 + n2 * h / 2.0); n3 = f3(t+h/2.0, xi - 1 + k2 * h / 2.0, yi - 1 + m2 * h / 2.0, zi - 1 + n2 * h / 2.0); k4 = f1(t+h/2.0, xi - 1 + k3 * h / 2.0, yi - 1 + m3 * h / 2.0, zi - 1 + n3 * h / 2.0); m4 = f2(t+h/2.0, xi - 1 + k3

6、 * h / 2.0, yi - 1 + m3 * h / 2.0, zi - 1 + n3 * h / 2.0); n4 = f3(t+h/2.0, xi - 1 + k3 * h / 2.0, yi - 1 + m3 * h / 2.0, zi - 1 + n3 * h / 2.0); xi = xi - 1 + h * (k1 + 2 * k2 + 2 * k3 + k4) / 6.0; yi = yi - 1 + h * (m1 + 2 * m2 + 2 * m3 + m4) / 6.0; zi = zi - 1 + h * (n1 + 2 * n2 + 2 * n3 + n4) / 6.0; +count;cout<<"t= "<<t+0.05<<" "<<"x= "<<xi<<" "<<"y= "<<yi<<" "<<"z= "<<zi<<endl;cout<<"共計(jì)算的次數(shù)為："<<count<<endl;return 0; vo