- #1

cppIStough

- 22

- 2

C++:

```
#include <iostream>
#include <vector>
class RK4 {
private:
double x0, y0, x, h;
public:
RK4(double x0_, double y0_, double x_, double h_)
{
x0 = x0_;
y0 = y0_;
x = x_;
h = h_;
std::vector<double> vec_x, vec_y;
int n = (x - x0)/h;
for (int i = 0; i < n; ++i)
{
vec_x[i] = i*h;
vec_y[i] = y0;
}
}
double f(double x_, double y_) {
return y_;
}
void rkf(std::vector<double> & vec_x, std::vector<double> & vec_y) {
for(int i = 0; i < vec_x.size(); ++i) {
double k1 = f(vec_x[i], vec_y[i]);
double k2 = f(vec_x[i] + h / 2, vec_y[i] + h * k1 / 2);
double k3 = f(vec_x[i] + h / 2, vec_y[i] + h * k2 / 2);
double k4 = f(vec_x[i] + h, vec_y[i] + h * k3);
vec_y[i+1] = vec_y[i] + h / 6 * (k1 + 2 * k2 + 2 * k3 + k4);
}
}
};
int main() {
double h = 0.0001;
double x0 = 0;
double y0 = 1;
double x = 1;
RK4 ob1(x0, y0, h, x);
return 0;
}
```