Лабы / 5lab 10 var

#include <iostream>
#include <cmath>
#include <vector>
#include <iomanip>
#include <Windows.h>
using namespace std;

double solution(double t) {
	return (0.00346020761)*exp(-pow(t, 2))*(2905*exp(4*t) - (15 + 68*t)*cos(t) + (8 + 17*t)*sin(t));
}

double f(double t, double x) {
	return -2 * (t - 2) * x + t*exp(-pow(t,2)) * cos(t);
}

double RK(double t, double h, double x) {
	double k1, k2, k3, k4;
	k1 = h * f(t, x);
	k2 = h * f(t + h / 3, x + k1/3);
	k3 = h * f(t + 2*h / 3, x - k1/3 + k2 );
	k4 = h * f(t + h, x + k1 - k2 + k3);
	return (k1 + 3 * k2 + 3 * k3 + k4) / 8;
}

int main() {
	SetConsoleCP(1251);
	SetConsoleOutputCP(1251);
	//Условия
	double t0 = 0, x0 = 10, h = 0.5, b = 6;
	int s = 4;
	double eps = 0.0001;

	vector<double> X, T;
	X.push_back(x0); T.push_back(t0);

	double x11, x21, x22; 
	bool flag = true;

	for (double t = t0; t < b; t += h) {
		/*x11 = X.back() + RK(t, h, X.back());
		x21 = X.back() + RK(t, h / 2, X.back());
		x22 = x21 + RK(t + h / 2, h / 2, x21);
		double D = abs(x11 - x22) / (pow(2.0, s) - 1);
		if (D >= eps)
			h = h / 2;*/

		while (flag) {

			x11 = X.back() + RK(t, h, X.back());
			x21 = X.back() + RK(t, h / 2, X.back());
			x22 = x21 + RK(t + h / 2, h / 2, x21);
			double D = abs(x11 - x22) / (pow(2.0, s) - 1);
			
			/*int s = 1;
			if (D < eps / (pow(2, s + 1)))
				h = h * 2;
			s++;*/

			if (D >= eps)
				h = h / 2;
			else
				flag = false;
		}
		X.push_back(x22); T.push_back(t + h); flag = true;
	}


	while (!X.empty()) 
	{
		/*cout << "X[";
		printf("%.6lf", T.back());
		cout << "] = ";
		printf("%.6lf", X.back());
		cout << '\t';
		printf("%.6lf", solution(T.back()));
		cout << " (";
		printf("%.6lf", abs(X.back() - solution(T.back())));
		cout<<")"<<endl;*/
		cout << "X[" << T.back() << "] = " << X.back() << '\t' << solution(T.back()) << " (" << abs(X.back() - solution(T.back())) << ")" << endl;
		X.pop_back(); T.pop_back();
	}
	system("pause");
}