add support for sweeping parameters

main.cpp

@@ -4,12 +4,7 @@
 
 #include "model.h"
 
-double omegafn(size_t i)
-{
-	return i*10000;
-}
-
-int main(int argc, char** argv)
+void runSingle()
 {
 	std::string modelStr("r{20e3}p{1e-7, 0.9}");
 
@@ -29,9 +24,10 @@ int main(int argc, char** argv)
 		std::cout<<"}\n";
 	}
 
+	Model::Range omega(0, 1e6, 50);
+
 	auto start = std::chrono::high_resolution_clock::now();
-	for(size_t i = 0; i < 101; ++i)
-		results.push_back(model.execute(omegafn(i)));
+	results = model.sweep(omega);
 	auto end = std::chrono::high_resolution_clock::now();
 	auto duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start);
 
@@ -39,5 +35,41 @@ int main(int argc, char** argv)
 		std::cout<<"omega: "<<res.omega<<" real = "<<res.im.real()<<" im = "<<res.im.imag()<<'\n';
 
 	std::cout<<"time taken: "<<duration.count()<<" us"<<'\n';
+}
+
+void sweepCb(std::vector<Model::DataPoint>& data)
+{
+	static size_t i = 0;
+	++i;
+	if((i & 0x3FF) == 0)
+		std::cout<<'.'<<std::flush;
+}
+
+void runSweep()
+{
+	std::string modelStr("r{20e3}p{1e-7, 0.9}");
+	std::vector<Model::DataPoint> results;
+
+	Model model(modelStr);
+
+	std::vector<Model::Range> parameters;
+	parameters.push_back(Model::Range(1e3, 50e3, 100));
+	parameters.push_back(Model::Range(1e-7, 20e-7, 100));
+	parameters.push_back(Model::Range(0.7, 1.2, 100));
+
+	Model::Range omega(0, 1e6, 25);
+
+	auto start = std::chrono::high_resolution_clock::now();
+	model.sweepParams(parameters, omega, &sweepCb);
+	auto end = std::chrono::high_resolution_clock::now();
+	auto duration = std::chrono::duration_cast<std::chrono::milliseconds>(end - start);
+
+	std::cout<<"\ntime taken: "<<duration.count()<<" ms"<<'\n';
+}
+
+int main(int argc, char** argv)
+{
+	runSingle();
+	runSweep();
 	return 0;
 }

model.cpp

@@ -230,6 +230,101 @@ std::vector<Componant*> Model::getFlatComponants()
 	return getFlatComponants();
 }
 
+std::vector<Model::DataPoint> Model::sweep(const Range& omega)
+{
+	double step = (omega.end - omega.start)/omega.count;
+	double currOmega = omega.start;
+
+	std::vector<DataPoint> results;
+	results.reserve(omega.count);
+	for(size_t i = 0; i < omega.count; ++i)
+	{
+		results.push_back(execute(currOmega));
+		currOmega+=step;
+	}
+	return results;
+}
+
+bool Model::sweepParams(const std::vector<Range>& componantRanges, const Range& omega, std::function<void(std::vector<DataPoint>&)> dataCb)
+{
+	size_t parametersCount = getFlatParametersCount();
+	if(componantRanges.size() != parametersCount)
+	{
+		std::cout<<"Error: a parameter range must be provided for eatch componant parameter\n";
+		return false;
+	}
+
+	for(size_t i = 0; i < parametersCount; ++i)
+	{
+		if(componantRanges[i].count == 0 || (componantRanges[i].count < 2 && componantRanges[i].start != componantRanges[i].end))
+		{
+			std::cout<<"Error: paramter range must specify at least one paramter point if only one paramer point is specified star and end must be the same\n";
+			return false;
+		}
+		else if(componantRanges[i].start > componantRanges[i].end)
+		{
+			std::cout<<"Error: paramter range end-start must be positive\n";
+			return false;
+		}
+	}
+
+	size_t stepsRequired = 1;
+	for(size_t i = 0; i < parametersCount; ++i)
+		stepsRequired*=componantRanges[i].count;
+
+	std::vector<double> currentParam(parametersCount, 0);
+	std::vector<double> stepSize(parametersCount, 0);
+	for(size_t i = 0; i < parametersCount; ++i)
+	{
+		currentParam[i] = componantRanges[i].start;
+		stepSize[i] = (componantRanges[i].end - componantRanges[i].start)/componantRanges[i].count;
+	}
+
+
+	std::cout<<"Executing sweep. Steps requried: "<<stepsRequired<<std::endl;
+
+	for(size_t i = 0; i < stepsRequired; ++i)
+	{
+		for(size_t i = 0; i < parametersCount; ++i)
+		{
+			currentParam[i] += stepSize[i];
+			if(currentParam[i] > componantRanges[i].end)
+				currentParam[i] = componantRanges[i].start;
+			else
+				break;
+		}
+		std::vector<DataPoint> result = sweep(omega);
+		dataCb(result);
+	}
+
+	return true;
+}
+
+bool Model::setFlatParameters(const std::vector<double>& parameters)
+{
+	if(parameters.size() != getFlatParametersCount())
+		return false;
+
+	size_t i = 0;
+	for(Componant* componant : getFlatComponants())
+	{
+		std::vector<double> params;
+		for(size_t j = 0; j < componant->paramCount(); ++j)
+			params.push_back(parameters[i++]);
+		componant->setParam(params);
+	}
+
+	return true;
+}
+
+size_t Model::getFlatParametersCount()
+{
+	size_t count = 0;
+	for(Componant* componant : getFlatComponants())
+		count += componant->paramCount();
+	return count;
+}
+
 char Model::getComponantChar(Componant* componant)
 {
 	if(dynamic_cast<Resistor*>(componant))
@@ -237,9 +332,10 @@ char Model::getComponantChar(Componant* componant)
 
 	if(dynamic_cast<Cap*>(componant))
 		return 'c';
-
 	if(dynamic_cast<Cpe*>(componant))
 		return 'p';
 
 	return 'x';
 }
+
+

model.h

@@ -3,6 +3,7 @@
 #include <complex>
 #include <string>
 #include <vector>
+#include <functional>
 
 #include "componant.h"
 
@@ -14,6 +15,15 @@ public:
 		std::complex<double> im;
 		double omega;
 	};
+	struct Range
+	{
+		double start;
+		double end;
+		size_t count;
+
+		Range(double startI, double endI, size_t countI): start(startI), end(endI), count(countI){}
+		Range() = default;
+	};
 private:
 	class Paralell: public Componant
 	{
@@ -49,7 +59,11 @@ private:
 public:
 	Model(const std::string& str);
 	DataPoint execute(double omaga);
+	std::vector<DataPoint> sweep(const Range& omega);
+	bool sweepParams(const std::vector<Range>& componantRanges, const Range& omega, std::function<void(std::vector<DataPoint>&)> dataCb);
 	std::string getModelStr();
 	std::vector<Componant*> getFlatComponants();
+	size_t getFlatParametersCount();
+	bool setFlatParameters(const std::vector<double>& parameters);
 	static char getComponantChar(Componant* componant);
 };