#include <iostream>
#include <functional>
#include <cmath>
#include <fstream>
#include <string>
#include <vector>

#define Nstep 3200
#define A -1.
#define B 1.
#define MAXIT 100000000
#define TOL 1e-9
#define M_PI		3.14159265358979323846

std::vector<double> Set1DGrid(double a,double b, std::function<double(double)> f);
void SetBC(double lower, double upper, std::vector<double> & solution);
void IterativeSolver(std::vector<double> & u,std::vector<double> & rhs, double tol);
void WriteFile(std::string filename, std::vector<double> f, std::vector<double> grid);
double dx = (B-A)/(Nstep-1);

std::function<double(double)> Grid = [](double x) { return x; };
std::function<double(double)> RHS = [](double x) { return -M_PI*exp(x)*(M_PI*exp(x)*cos(M_PI*exp(x))+sin(M_PI*exp(x))); };
std::function<double(double)> FirstGuess = [](double x) { return 0; };
std::function<double(double)> ExactResult = [](double x) { return cos(M_PI*exp(x)); };

int main (int argc, char *argv[]){
	std::vector<double> MyGrid = Set1DGrid(A,B, Grid);
	std::vector<double> rhs = Set1DGrid(A,B, RHS);
	std::vector<double> u = Set1DGrid(A,B, FirstGuess);
	std::vector<double> exact = Set1DGrid(A,B, ExactResult);
	SetBC(cos(M_PI*exp(-1)),cos(M_PI*exp(1)),u);
	IterativeSolver(u,rhs,TOL);
	WriteFile("solution.txt",u,MyGrid);
	WriteFile("exact.txt",exact,MyGrid);
	return 0;
}
std::vector<double> Set1DGrid(double a,double b, std::function<double(double)> f){
	std::vector<double> grid(Nstep);
	for (int i=0;i<Nstep;++i)
		grid[i]=f(a+i*(b-a)/(Nstep-1));
	return grid;
}
void SetBC(double lower, double upper,std::vector<double> & solution){
	solution[0]=lower;
	solution[Nstep-1]=upper;
}
void IterativeSolver(std::vector<double> & u,std::vector<double> & rhs, double tol){
	int iter=0;
	std::vector<double> prevU = u;
	double max;
	while (iter<=MAXIT){
		max=0;
		for (int i=1;i<Nstep-1;++i){
			u[i]=(prevU[i+1]+prevU[i-1]-dx*dx*rhs[i])*0.5;
			max = fabs(u[i]-prevU[i]) > max? fabs(u[i]-prevU[i]):max;
		}
		if (max < tol){
			std::cout<<"#Solution converged in "<<iter<<" iterations with error "<<
			TOL<<std::endl;
			break;
		}
		prevU=u;
		++iter;
	}
}
void WriteFile(std::string filename,std::vector<double> f,std::vector<double> grid){
	std::ofstream outfile;
    outfile.open(filename.c_str());
    if (!outfile.is_open())
    {
        std::cout<<"Cannot open file "<<filename<<"."<<std::endl;
        return;
    }
    for (int i=0;i<Nstep;++i)
		outfile<<grid[i]<<"  "<<f[i]<<std::endl;
	outfile.close();
}