#include "GetdpMaterialWriter.hpp"
#include <sstream>
GetdpMaterialWriter::GetdpMaterialWriter(Material mat_, std::string entity_name_): mat(mat_), entity_name(entity_name_)
{
	isHarmonic = false;
}


std::string GetdpMaterialWriter::writeData()
{

	std::string outString; 
	for (auto& data : mat.getMaterialData())
	{
		if ((data.first == "BHModel"))
		{
			if (mat.getMaterialDataFlag("useTForMagFlux"))
				outString.append(magneticBHWithT(data.second.getValue(), isHarmonic, mat.getMaterialValueByName("TCurie"), mat.getMaterialValueByName("TCurieExponent")));
			else
				outString.append(magneticBH( data.second.getValue(), isHarmonic));
		}
		else if ((data.first == "mu"))
		{
			outString.append("  nu[" + entity_name + "] = 1/(" + data.second.getValue() + "); \n");
			outString.append("  nu_nu0[" + entity_name + "] = 1/(" + data.second.getValue() + " * mu0); \n");
		}
		else if ((data.first == "mu_C"))
		{
			if (mat.getMaterialDataFlag("useTForMagFlux"))
			{
				outString.append(muWithT(data.second.getValue(), mat.getMaterialValueByName("TCurie"), mat.getMaterialValueByName("TCurieExponent")));
			}
			else
			{
				outString.append("  nu_nu0[" + entity_name + "] = 1/(" + data.second.getValue() + " * mu0); \n");
			}
		}
		else if (!data.second.isBoolean())
		{
			outString.append("  " + data.first + "[" + entity_name + "] = " + data.second.getValueWithEnds() + "; \n");
		}
	}
	return outString;
}



std::string GetdpMaterialWriter::magneticBH(std::string data, bool harmonicFlag)
{
	std::string outString;
	std::string matName = "BHModel";

	if (harmonicFlag)
	{
		outString.append("  " + matName + "_" + entity_name + " = {");
		auto BHdata = effectiveBH(data.substr(23));
		for (int i = 0; i < BHdata.size(); i++)
		{
			if (i != 0)
				outString.append(",");
			outString.append(std::to_string(BHdata[i].first) + ",");
			outString.append(std::to_string(BHdata[i].second));
		}
		outString.append("}; \n");
	}
	else
		outString.append("  " + matName + "_" + entity_name + " = " + data.substr(23) + "; \n");

	outString.append("  " + matName + "_" + entity_name + "_b = {}; \n");
	outString.append("  " + matName + "_" + entity_name + "_h = {}; \n");
	outString.append("  For ips In {0: #" + matName + "_" + entity_name + "()/2-1} \n");
	outString.append("  	" + matName + "_" + entity_name + "_h += " + matName + "_" + entity_name + "(2*ips) ; \n");
	outString.append("  	" + matName + "_" + entity_name + "_b += " + matName + "_" + entity_name + "(2*ips+1) ; \n");
	outString.append("  EndFor \n");
	outString.append("  " + matName + "_" + entity_name + "_b2() = " + matName + "_" + entity_name + "_b()^2; \n");
	outString.append("  " + matName + "_" + entity_name + "_nu() = " + matName + "_" + entity_name + "_h()/ " + matName + "_" + entity_name + "_b(); \n");
	outString.append("  " + matName + "_" + entity_name + "_nu(0) = " + matName + "_" + entity_name + "_nu(1); \n");
	outString.append("  " + matName + "_" + entity_name + "_nu_b2() = ListAlt[" + matName + "_" + entity_name + "_b2(), " + matName + "_" + entity_name + "_nu()]; \n");
	outString.append("  nu_nu0[" + entity_name + "] = InterpolationLinear[SquNorm[$1]]{" + matName + "_" + entity_name + "_nu_b2()} ; \n");
	outString.append("  dnudb2[" + entity_name + "] = dInterpolationLinear[SquNorm[$1]]{" + matName + "_" + entity_name + "_nu_b2()} ; \n");
	outString.append("  dhdb_NL[" + entity_name + "] = 2*dnudb2[$1#1] * SquDyadicProduct[$1#1] ; \n");

	return outString;
}


std::string GetdpMaterialWriter::magneticBHWithT(std::string data, bool harmonicFlag, std::string TCurie, std::string Texp)
{
	std::string outString;
	std::string matName = "BHModel";

	outString.append("  nu0= 1/mu0; \n");
	outString.append("  TFuncBelowTc_" + entity_name + "[] = (1-($1/" + TCurie + ")^" + Texp + "); \n");
	outString.append("  TFuncAboveTc_" + entity_name + "[] = 0; \n");
	outString.append("  TFunc_" + entity_name + "[] = $1 < " + TCurie + " ? TFuncBelowTc_" + entity_name + "[$1] : TFuncAboveTc_" + entity_name + "[]; \n");

	if (harmonicFlag)
	{
		outString.append("  " + matName + "_" + entity_name + " = {");
		auto BHdata = effectiveBH(data.substr(23));
		for (int i = 0; i < BHdata.size(); i++)
		{
			if (i != 0)
				outString.append(",");
			outString.append(std::to_string(BHdata[i].first) + ",");
			outString.append(std::to_string(BHdata[i].second));
		}
		outString.append("}; \n");
	}
	else
		outString.append("  " + matName + "_" + entity_name + " = " + data.substr(23) + "; \n");

	outString.append("  " + matName + "_" + entity_name + "_b = {}; \n");
	outString.append("  " + matName + "_" + entity_name + "_h = {}; \n");
	outString.append("  For ips In {0: #" + matName + "_" + entity_name + "()/2-1} \n");
	outString.append("  	" + matName + "_" + entity_name + "_h += " + matName + "_" + entity_name + "(2*ips) ; \n");
	outString.append("  	" + matName + "_" + entity_name + "_b += " + matName + "_" + entity_name + "(2*ips+1) ; \n");
	outString.append("  EndFor \n");

	outString.append("  " + matName + "_" + entity_name + "_b2() = " + matName + "_" + entity_name + "_b()^2; \n");
	outString.append("  " + matName + "_" + entity_name + "_nu() = " + matName + "_" + entity_name + "_h()/ " + matName + "_" + entity_name + "_b(); \n");
	outString.append("  " + matName + "_" + entity_name + "_nu(0) = " + matName + "_" + entity_name + "_nu(1); \n");
	outString.append("  " + matName + "_" + entity_name + "_mu() = " + matName + "_" + entity_name + "_b()/ " + matName + "_" + entity_name + "_h(); \n");
	outString.append("  " + matName + "_" + entity_name + "_mu(0) = " + matName + "_" + entity_name + "_mu(1); \n");
	outString.append("  " + matName + "_" + entity_name + "_nu_b2() = ListAlt[" + matName + "_" + entity_name + "_b2(), " + matName + "_" + entity_name + "_nu()]; \n");
	outString.append("  " + matName + "_" + entity_name + "_mu_b() = ListAlt[" + matName + "_" + entity_name + "_b(), " + matName + "_" + entity_name + "_mu()]; \n");
	outString.append("  " + matName + "_" + entity_name + "_mu_b2() = ListAlt[" + matName + "_" + entity_name + "_b2(), " + matName + "_" + entity_name + "_mu()]; \n");
	outString.append("  nu_" + entity_name + "_0[] = InterpolationLinear[SquNorm[$1]]{" + matName + "_" + entity_name + "_nu_b2()} ; \n");
	outString.append("  nu_nu0[" + entity_name + "] = 1/(TFunc_" + entity_name + "[$2] * (1/nu_" + entity_name + "_0[$1] - mu0) + mu0) ; \n");
	outString.append("  dnudb2[" + entity_name + "] = dInterpolationLinear[SquNorm[$1]]{" + matName + "_" + entity_name + "_nu_b2()}* TFunc_" + entity_name + "[$2] * nu_nu0[$1,$2]^2 / nu_" + entity_name + "_0[$1]^2; \n");


	// INVERSE nu(T) formulation. Do not use!!!
	//	outString.append("  nu[" + entity_name + "] = nu0 + (nu_" + entity_name + "_0[$1] - nu0)*TFunc[$2] ; \n");
	//	outString.append("  dnudb2["+ entity_name + "] = - dInterpolationLinear[SquNorm[$1]]{BHModel_tube_mu_b2()}* TFunc[$2]; \n");

	outString.append("  dhdb_NL[" + entity_name + "] = 2*dnudb2[$1,$2] * SquDyadicProduct[$1] ; \n");
	return outString;
}



std::vector<std::pair<float, float>> GetdpMaterialWriter::effectiveBH(std::string data)
{
	int method = 0;
	/* method = 1 - active energy method, else - simple energy method */

	std::vector<std::pair<float, float> > outVec;
	std::vector<float> bVec;
	std::vector<float> bVec_new;
	std::vector<float> hVec;
	std::vector<float> dataVec;
	std::stringstream bhStream(data.substr(1, data.size() - 2));

	float dp;
	while (bhStream >> dp)
	{
		dataVec.push_back(dp);

		if (bhStream.peek() == ',')
			bhStream.ignore();
	}


	for (size_t i = 0; i < dataVec.size() / 2; i++)
	{
		hVec.push_back(dataVec[i * 2]);
		bVec.push_back(dataVec[i * 2 + 1]);
	}

	if ((bVec[0] != 0) and (hVec[0] != 0))
	{
		bVec.insert(bVec.begin(), 0.0);
		hVec.insert(hVec.begin(), 0.0);
	}

	if (bVec.size() == 1)
		return { std::make_pair(hVec[0], bVec[0]) };

	bVec_new.push_back(bVec[0]);

	if (method != 1)
	{
		for (size_t i = 1; i < bVec.size(); i++)
		{
			float locSum = 0;
			for (size_t j = 1; j <= i; j++)
			{
				locSum = locSum + (hVec[j] - hVec[j - 1]) * (bVec[j - 1] + bVec[j]) / 2;
			}
			bVec_new.push_back((2 / hVec[i]) * locSum);
		}
	}
	else
	{
		//active energy method not implemented
	}

	for (size_t i = 0; i < bVec.size(); i++)
	{
		outVec.push_back(std::make_pair(hVec[i], bVec_new[i]));
	}
	return outVec;
}


void GetdpMaterialWriter::setHarmonic(bool b )
{
	isHarmonic = b;
}

std::string GetdpMaterialWriter::muWithT(std::string data, std::string TCurie, std::string Texp)
{
	std::string outString;

	outString.append("  nu0= 1/mu0; \n");

	outString.append("  TFuncBelowTc_" + entity_name + "[] = (1-($1/" + TCurie + ")^" + Texp + "); \n");
	outString.append("  TFuncAboveTc_" + entity_name + "[] = 0; \n");
	outString.append("  TFunc_" + entity_name + "[] = $1 < " + TCurie + " ? TFuncBelowTc_" + entity_name + "[$1] : TFuncAboveTc_" + entity_name + "[]; \n");
	outString.append("  nu_nu0[" + entity_name + "] = 1/( (1+(" + data + "-1)*TFunc_" + entity_name + "[$2])*mu0); \n");

	return outString;
}