//#include"Ctridiag.h"
//#include<iostream>
//#include<math.h>
//
////Ctridiag::Ctridiag() {};
//Ctridiag::Ctridiag(double& dx):sub(dx),sup(dx),diag(dx) {}
//Ctridiag::Ctridiag(int n, double& dx) : sub(dx), sup(dx), diag(dx) {
//	this->init(n,dx);
//}
//
//void Ctridiag::init(int n, double& dx) {
//	this->diag = Cvec(n, dx);
//	this->sub = Cvec(n - 1, dx);
//	this->sup = Cvec(n - 1, dx);
//}
//
//Ctridiag::Ctridiag(const Ctridiag& in, double& dx) :sub(dx), sup(dx), diag(dx) {
//	*this = in;
//}
//
//Ctridiag::Ctridiag(const int n, double& dx, Creal& sub_in, Creal& diag_in, Creal& sup_in) : sub(dx), sup(dx), diag(dx) {
//	this->init(n,dx);
//	this->sub = sub_in;
//	this->diag = diag_in;
//	this->sup = sup_in;
//}
//
//Ctridiag::Ctridiag(const int n, double& dx, double sub_in, double diag_in, double sup_in) :sub(dx), sup(dx), diag(dx) {
//	this->init(n, dx);
//	this->sub = sub_in;
//	this->diag = diag_in;
//	this->sup = sup_in;
//}
//
//Ctridiag::~Ctridiag() {
//
//}
//
//void Ctridiag::print() {
//	std::cout << "[(-,-) , " << diag[0] << " , " << sup[0] << "]" << std::endl;
//	for (int i = 1; i < diag.size - 1; i++) {
//		std::cout << "[" << sub[i - 1] << " , " << diag[i] << " , " << sup[1] << "]" << std::endl;
//	}
//	std::cout << "[" << sub[diag.size - 2] << " , " << diag[diag.size - 1] << " , (-,-)]" << std::endl;
//}
//
//Ctridiag Ctridiag::operator=(double& diag_in) {
//	if (sub.size == 0) {
//		sub = 0;
//		sup = 0;
//	}
//	this->diag = diag_in;
//	return *this;
//}
//
//Ctridiag Ctridiag::operator=(Cvec& in) {
//	if (sub.size == 0) {
//		sub.init(in.size - 1);
//		sub = 0;
//		sup.init(in.size - 1);
//		sup = 0;
//	}
//	diag = in;
//	return *this;
//}
//
//Ctridiag Ctridiag::operator+(const Ctridiag& in) {
//	Ctridiag out(in.diag.dx);
//	out.sub = sub + in.sub;
//	out.diag = diag + in.diag;
//	out.sup = sup + in.sup;
//
//	return out;
//}
//Ctridiag Ctridiag::operator-(const Ctridiag& in) {
//	Ctridiag out(in.diag.dx);
//	out.sub = sub - in.sub;
//	out.diag = diag - in.diag;
//	out.sup = sup - in.sup;
//
//	return out;
//}
//
//Ctridiag Ctridiag::operator/(double& in) {
//	Ctridiag out(diag.dx);
//	out.sub = sub / in;
//	out.diag = diag / in;
//	out.sup = sup / in;
//	return out;
//}
//Ctridiag Ctridiag::operator*(const double& in) {
//	Ctridiag out(diag.dx);
//	out.sub = sub * in;
//	out.diag = diag * in;
//	out.sup = sup * in;
//	return out;
//}
//
//Cvec Ctridiag::operator*(const Cvec& invec){
//	Cvec res(invec.size, invec.dx);
//	long i, n;
//	if (invec.size != this->diag.size) {
//		std::cout << "ERROR: Rows of Matrix must agree with size of vector" << std::endl;
//		exit(0);
//	}
//	n = this->diag.size;
//
//	// first row
//	res.vec[0] = this->diag[0] * invec.vec[0] + this->sup[0] * invec.vec[1];
//
//	// inner region
//	for (i = 1; i < (n - 1); i++) res.vec[i] = this->sub[i - 1] * invec.vec[i - 1] + this->diag[i] * invec.vec[i] + this->sup[i] * invec.vec[i + 1];
//
//	// last row
//	res.vec[n - 1] = this->sub[i - 1] * invec.vec[i - 1] + this->diag[i] * invec.vec[i];
//	
//	return res;
//}
//
//Cvec Ctridiag::solve(const Cvec& b)const{
//	Cvec x(b.size, b.dx);
//	int i, n;
//	n = this->diag.size;
//
//	for (i = 0; i < n; i++) x.vec[i] = b.vec[i];    // copy inhomogeneity to x vector
//	Ctridiag a(n, b.dx);   // make temporary matrix
//	a = *this; // copy matrix into temorary one
//
//	// begin solve
//	// forward sweep to eliminate sub diag elements (dont need to update sub diag elements, as it is clear that they are set to zero)
//	for (i = 0; i < (n - 1); i++)
//	{
//		a.diag[i + 1] += a.sup[i] * (-a.sub[i] / a.diag[i]);
//		x.vec[i + 1] += x.vec[i] * (-a.sub[i] / a.diag[i]);
//	}
//
//	// backward loop to eliminate super diag (again, dont need to update super diag elements, as it is clear that they are set to zero)
//	for (i = (n - 2); i >= 0; i--)
//	{
//		x.vec[i + 1] /= a.diag[i + 1];
//		x.vec[i] += x.vec[i + 1] * (-a.sup[i]);
//	}
//	x.vec[0] /= a.diag[0];
//
//	return x;
//}