Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- #include <iostream>
- #include <cmath>
- using namespace std;
- struct Point {
- double x;
- double y;
- Point() {x = y = 0;}
- Point(double _x, double _y) : x(_x), y(_y) {};
- };
- double df(Point x) {return 0.2 * x.x + pow(x.y, 2);}
- void print(Point* x, int n) {
- cout << endl << "n = " << n << endl << endl;
- cout << "===== points =====" << endl;
- for (int i = 0; i < n + 1; i++) {
- cout << "[" << i << "]: " << "x = " << x[i].x << " y = " << x[i].y << endl;
- cout << endl;
- }
- }
- Point* rkSolve(Point s, double a, double b, double df(Point x), int n) {
- double k1, k2, k3, k4;
- double h = (b - a) / n;
- Point* res = new Point[n + 1];
- for (int i = 0; i < n + 1; i++) {
- res[i] = s;
- k1 = df(s);
- k2 = df(Point(s.x + 0.5 * h, s.y + 0.5 * h * k1));
- k3 = df(Point(s.x + 0.5 * h, s.y + 0.5 * h * k2));
- k4 = df(Point(s.x + h, s.y + h * k3));
- s.y = s.y + (h / 6) * (k1 + 2 * k2 + 2 * k3 + k4);
- s.x += h;
- }
- return res;
- }
- double checkEps(Point* s1, Point* s2, int n) {
- double eps = -1, curEps;
- for (int i = 1; i < n + 1; i++) {
- curEps = fabs(s1[i].y - s2[i * 2].y) / 15;
- eps = max(eps, curEps);
- }
- return eps;
- }
- void solve(Point s, double a, double b, double df(Point x), double eps) {
- int n = 10;
- Point *s1, *s2;
- do {
- n *= 2;
- s1 = rkSolve(s, a, b, df, n);
- s2 = rkSolve(s, a, b, df, n * 2);
- } while (checkEps(s1, s2, n) >= eps);
- print(s1, n);
- }
- int main() {
- Point s;
- double a, b, eps;
- cout << "Enter start point: ";
- cin >> s.x >> s.y;
- cout << "Enter a & b: ";
- cin >> a >> b;
- cout << "Enter eps: ";
- cin >> eps;
- solve(s, a, b, df, eps);
- return 0;
- }
Add Comment
Please, Sign In to add comment