编写devC++实现用四阶龙格-库塔法求解一阶常微分方程组初值问题y1'=1/(y2-x) y2'=1-1/y1 y1(0)=1 y2(0)=1 x∈[0, 2],h=0.1。 要求: (1) 步长h的值从键盘终端输入; (2) 将每一步的计算结果打印输出
时间: 2023-07-16 22:16:16 浏览: 111
以下是用DevC++编写的四阶龙格-库塔法求解一阶常微分方程组初值问题的代码(注释中有详细的解释):
```c++
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define MAXN 1000
double y1[MAXN+1], y2[MAXN+1]; // 存储y1和y2的数组
double f1(double x, double y1, double y2) { return 1.0 / (y2 - x); } // y1'的函数
double f2(double x, double y1, double y2) { return 1.0 - 1.0 / y1; } // y2'的函数
int main() {
double x0 = 0, y10 = 1, y20 = 1; // 初始值
double x, y1n, y2n, h; // 当前的x、y1、y2和步长h
int n; // 迭代次数
printf("请输入步长h:");
scanf("%lf", &h);
x = x0; y1[0] = y10; y2[0] = y20; // 初始化
printf("n\tx\ty1\ty2\n"); // 输出表头
printf("%d\t%.1lf\t%.6lf\t%.6lf\n", 0, x, y1[0], y2[0]); // 输出初始值
for (n = 1; x <= 2; n++) { // 迭代计算直到x>2
double k11, k12, k21, k22, k31, k32, k41, k42; // 计算k值的中间变量
k11 = h * f1(x, y1[n-1], y2[n-1]);
k12 = h * f2(x, y1[n-1], y2[n-1]);
k21 = h * f1(x + h/2, y1[n-1] + k11/2, y2[n-1] + k12/2);
k22 = h * f2(x + h/2, y1[n-1] + k11/2, y2[n-1] + k12/2);
k31 = h * f1(x + h/2, y1[n-1] + k21/2, y2[n-1] + k22/2);
k32 = h * f2(x + h/2, y1[n-1] + k21/2, y2[n-1] + k22/2);
k41 = h * f1(x + h, y1[n-1] + k31, y2[n-1] + k32);
k42 = h * f2(x + h, y1[n-1] + k31, y2[n-1] + k32);
y1n = y1[n-1] + (k11 + 2*k21 + 2*k31 + k41) / 6;
y2n = y2[n-1] + (k12 + 2*k22 + 2*k32 + k42) / 6;
x += h;
y1[n] = y1n;
y2[n] = y2n;
printf("%d\t%.1lf\t%.6lf\t%.6lf\n", n, x, y1n, y2n); // 输出计算结果
}
return 0;
}
```
注意:这段代码的步长h是要从键盘终端输入的。