#include "C_DGLSolver.h"
using namespace std;

C_DGLSolver::C_DGLSolver() {
}

// Konstruktor für DGL 1. Ordnung
C_DGLSolver::C_DGLSolver(CMyVektor(*f_DGL_System)(CMyVektor y, double x)) {
	isDGLSystem = true;
	fDglSystemZeiger = f_DGL_System;
}

// Konstruktor für DGL n-ter Ordnung
C_DGLSolver::C_DGLSolver(double(*f_DGL_nterOrdnung)(CMyVektor y, double x)) {
	isDGLSystem = false;
	fDglNterOrdnungZeiger = f_DGL_nterOrdnung;
}

C_DGLSolver::~C_DGLSolver() {
}

CMyVektor C_DGLSolver::euler(CMyVektor y, double xStart, double XEnd, double schritte) {
	double schrittweite = (XEnd - xStart) / schritte;
	CMyVektor yCopy(y.getDimension());
	CMyVektor yAbl(y.getDimension());
	yCopy = y;

	double x = xStart;
	int s = 0;
	cout << "h = " << schrittweite << endl;

	for (; s < schritte; s++, x += schrittweite) {
		yAbl = ableitungen(yCopy, x);

		cout << "Schritt " << s << ": " << endl
			<< "\t" << " 𝑥 = " << x << endl
			<< "\t" << " 𝑦 = " << yCopy << endl
			<< "\t" << "𝑦' = " << yAbl << endl << endl;

		for (int i = 1; i <= yCopy.getDimension(); i++)
			yCopy.setElement(i, yCopy.getElement(i) + schrittweite * yAbl.getElement(i));
	}

	cout << "Ende bei " << endl;
	cout << "\t" << " 𝑥 = " << x << endl;
	cout << "\t" << " 𝑦 = " << yCopy << endl;
	return yCopy;
}

CMyVektor C_DGLSolver::heun(CMyVektor y, double xStart, double XEnd, double schritte) {
	double schrittweite = (XEnd - xStart) / schritte;
	CMyVektor yCopy(y.getDimension());
	CMyVektor yCopyTest(y.getDimension());
	yCopy = y;

	CMyVektor yAbl(y.getDimension());
	CMyVektor yAblTest(y.getDimension());
	CMyVektor fx_mStg(y.getDimension());

	double x = xStart;
	int s = 0;
	cout << "h = " << schrittweite << endl;

	for (; s < schritte; x += schrittweite, s++) {
		yAbl = ableitungen(yCopy, x);

		cout << "Schritt " << s << ": " << endl;
		cout << "\t" << " 𝑥  = " << x << endl;
		cout << "\t" << " 𝑦  = " << yCopy << endl;

		cout << "\t" << "𝑦ₒ' = " << yAbl << endl << endl;


		for (int i = 1; i <= yCopy.getDimension(); i++)
			yCopyTest.setElement(i, yCopy.getElement(i) + schrittweite * yAbl.getElement(i));

		yAblTest = ableitungen(yCopyTest, x + schrittweite);

		for (int i = 1; i <= fx_mStg.getDimension(); i++)
			fx_mStg.setElement(i, 0.5 * (yAbl.getElement(i) + yAblTest.getElement(i)));

		for (int i = 1; i <= yCopy.getDimension(); i++)
			yCopy.setElement(i, yCopy.getElement(i) + schrittweite * fx_mStg.getElement(i));

		cout << "\t" << " 𝑦ₜ = " << yCopyTest << endl;
		cout << "\t" << "𝑦ₜ' = " << yAblTest << endl << endl;
		cout << "\t" << "𝑦ₘ' = " << fx_mStg << endl << endl;
	}

	cout << "Ende bei " << endl
		<< "\t" << " 𝑥  = " << x << endl
		<< "\t" << " 𝑦  = " << yCopy << endl;

	return yCopy;
}

CMyVektor C_DGLSolver::ableitungen(CMyVektor y, double x) {
	CMyVektor result(y.getDimension());

	if (isDGLSystem) {
		result = fDglSystemZeiger(y, x);
	}
	else {
		for (int i = 1; i < y.getDimension(); i++)
			result.setElement(i, y.getElement(i + 1));
			
		result.setElement(y.getDimension(), fDglNterOrdnungZeiger(y, x));
	}
	return result;
}