Constant phase element added

CMakeLists.txt

@@ -2,7 +2,14 @@ cmake_minimum_required(VERSION 3.20)
 
 project(ehigeneratorsa)
 
-set(SRC_FILES  main.cpp cap.cpp resistor.cpp model.cpp tokenize.cpp)
+set(SRC_FILES
+	main.cpp
+	cap.cpp
+	resistor.cpp
+	constantphase.cpp
+	model.cpp
+	tokenize.cpp
+)
 
 find_package(PkgConfig REQUIRED)
 pkg_check_modules(SPICE REQUIRED ngspice)
@@ -10,6 +17,6 @@ pkg_check_modules(SPICE REQUIRED ngspice)
 add_executable(${PROJECT_NAME} ${SRC_FILES})
 target_link_libraries(${PROJECT_NAME} ${SPICE_LIBRARIES})
 target_include_directories(${PROJECT_NAME} PUBLIC ${SPICE_CFLAGS} .)
-set_target_properties(${PROJECT_NAME} PROPERTIES COMPILE_FLAGS "-std=c++20 -Wall -O0 -march=native -g" LINK_FLAGS "-flto")
+set_target_properties(${PROJECT_NAME} PROPERTIES COMPILE_FLAGS "-std=c++20 -Wall -O3 -march=native -g" LINK_FLAGS "-flto")
 
 install(TARGETS ${PROJECT_NAME} RUNTIME DESTINATION bin)

cap.cpp

@@ -51,3 +51,8 @@ void Cap::setParam(const std::vector<double>& param)
 
 	_C = param[0];
 }
+
+size_t Cap::paramCount()
+{
+	return 1;
+}

cap.h

@@ -14,4 +14,5 @@ public:
 	virtual std::complex<double> execute(double omega) override;
 	virtual std::vector<double> getParam() override;
 	virtual void setParam(const std::vector<double>& param) override;
+	virtual size_t paramCount() override;
 };

componant.h

@@ -17,4 +17,5 @@ class Componant
 			return std::vector<double>();
 		};
 		virtual void setParam(const std::vector<double>& param){};
+		virtual size_t paramCount(){return 0;}
 };

constantphase.cpp

+#include "constantphase.h"
+#include "tokenize.h"
+#include <cstdlib>
+#include <math.h>
+
+Cpe::Cpe(double q, double alpha): _Q(q), _alpha(alpha)
+{
+
+}
+
+Cpe::Cpe(std::string paramStr)
+{
+	std::vector<std::string> tokens = tokenize(paramStr, ',');
+	if(tokens.size() < 2)
+	{
+		std::cout<<"Warning: to few parameters in "<<__func__<<" parameter string: "<<paramStr<<'\n';
+		_Q = 1e-7;
+		_alpha = 0.9;
+		return;
+	}
+	else
+	{
+		try
+		{
+			_Q = std::stod(tokens[0]);
+			_alpha = std::stod(tokens[1]);
+		}
+		catch(const std::invalid_argument& ia)
+		{
+			std::cout<<"Warning: cant parse parameter in "<<__func__<<" parameter: "<<tokens[0]<<'\n';
+			_Q = 1e-7;
+			_alpha = 0.9;
+		}
+	}
+}
+
+std::complex<double> Cpe::execute(double omega)
+{
+	return std::complex<double>((1.0/(_Q*pow(omega, _alpha)))*cos((M_PI/2)*_alpha),
+	                           0-(1.0/(_Q*pow(omega, _alpha)))*sin((M_PI/2)*_alpha));
+}
+
+std::vector<double> Cpe::getParam()
+{
+	return std::vector<double>({_Q, _alpha});
+}
+
+void Cpe::setParam(const std::vector<double>& param)
+{
+	if(param.size() < 2)
+	{
+		std::cout<<"Warning: invalid parameter list sent to "<<__func__<<'\n';
+		return;
+	}
+
+	_Q = param[0];
+	_alpha = param[1];
+}
+
+size_t Cpe::paramCount()
+{
+	return 2;
+}

constantphase.h

+#pragma once
+#include <complex>
+#include <string>
+#include <vector>
+#include "componant.h"
+
+class Cpe: public Componant
+{
+private:
+	double _Q;
+	double _alpha;
+public:
+	Cpe(std::string paramStr);
+	Cpe(double q = 1e-7, double alpha = 0.9);
+	virtual std::complex<double> execute(double omega) override;
+	virtual std::vector<double> getParam() override;
+	virtual void setParam(const std::vector<double>& param) override;
+	virtual size_t paramCount() override;
+};

main.cpp

@@ -6,18 +6,31 @@
 
 double omegafn(size_t i)
 {
-	return i*10;
+	return i*10000;
 }
 
 int main(int argc, char** argv)
 {
-	std::string modelStr("r{1e3}-(c{1e-6}r{1e3})");
+	std::string modelStr("r{20e3}p{1e-7, 0.9}");
 
 	std::vector<Model::DataPoint> results;
 
 	Model model(modelStr);
+	std::cout<<"Compnants: \n";
+	for(Componant* componant : model.getFlatComponants())
+	{
+		std::cout<<Model::getComponantChar(componant)<<"{";
+		for(size_t i = 0; i < componant->paramCount(); ++i)
+		{
+			std::cout<<componant->getParam()[i];
+			if(i != componant->paramCount()-1)
+				std::cout<<", ";
+		}
+		std::cout<<"}\n";
+	}
+
 	auto start = std::chrono::high_resolution_clock::now();
-	for(size_t i = 0; i < 50; ++i)
+	for(size_t i = 0; i < 101; ++i)
 		results.push_back(model.execute(omegafn(i)));
 	auto end = std::chrono::high_resolution_clock::now();
 	auto duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start);

model.cpp

@@ -3,6 +3,7 @@
 #include "tokenize.h"
 #include "cap.h"
 #include "resistor.h"
+#include "constantphase.h"
 
 Model::Paralell::Paralell(std::vector<Componant*> componantsIn): componants(componantsIn)
 {
@@ -139,6 +140,11 @@ Componant *Model::processBracket(std::string& str)
 					componants.push_back(new Resistor(getParamStr(nodeStr, i)));
 					i = opposingBraket(nodeStr, i, '}');
 					break;
+				case 'p':
+				case 'P':
+					componants.push_back(new Cpe(getParamStr(nodeStr, i)));
+					i = opposingBraket(nodeStr, i, '}');
+					break;
 				case '{':
 					i = opposingBraket(nodeStr, i, '}');
 				case '}':
@@ -193,3 +199,47 @@ Model::DataPoint Model::execute(double omega)
 	}
 	return DataPoint({std::complex<double>(0,0), 0});
 }
+
+void Model::addComponantToFlat(Componant* componant)
+{
+	Paralell* paralell = dynamic_cast<Paralell*>(componant);
+	if(paralell)
+	{
+		for(Componant* element : paralell->componants)
+			addComponantToFlat(element);
+		return;
+	}
+
+	Serial* serial = dynamic_cast<Serial*>(componant);
+	if(serial)
+	{
+		for(Componant* element : serial->componants)
+			addComponantToFlat(element);
+		return;
+	}
+
+	_flatComponants.push_back(componant);
+}
+
+std::vector<Componant*> Model::getFlatComponants()
+{
+	if(!_flatComponants.empty())
+		return _flatComponants;
+
+	addComponantToFlat(_model);
+	return getFlatComponants();
+}
+
+char Model::getComponantChar(Componant* componant)
+{
+	if(dynamic_cast<Resistor*>(componant))
+		return 'r';
+
+	if(dynamic_cast<Cap*>(componant))
+		return 'c';
+
+	if(dynamic_cast<Cpe*>(componant))
+		return 'p';
+
+	return 'x';
+}

model.h

@@ -38,6 +38,7 @@ private:
 	Componant *processBrackets(std::string& str, size_t& bracketCounter);
 	Componant *processBracket(std::string& str);
 	std::string getParamStr(const std::string& str, size_t index);
+	void addComponantToFlat(Componant* componant);
 
 private:
 	Componant *_model = nullptr;
@@ -50,4 +51,5 @@ public:
 	DataPoint execute(double omaga);
 	std::string getModelStr();
 	std::vector<Componant*> getFlatComponants();
+	static char getComponantChar(Componant* componant);
 };

resistor.cpp

@@ -50,3 +50,8 @@ void Resistor::setParam(const std::vector<double>& param)
 
 	_R = param[0];
 }
+
+size_t Resistor::paramCount()
+{
+	return 1;
+}

resistor.h

@@ -13,4 +13,5 @@ public:
 	virtual std::complex<double> execute(double omega)  override;
 	virtual std::vector<double> getParam() override;
 	virtual void setParam(const std::vector<double>& param) override;
+	virtual size_t paramCount() override;
 };