Added ISystemInitialization interface to cpp runtime

git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@11133 f25d12d1-65f4-0310-ae8a-bbce733d8d8e

Compiler/susan_codegen/SimCode/SimCodeCpp.tpl

@@ -222,6 +222,12 @@ case SIMCODE(modelInfo = MODELINFO(__)) then
    <%giveZeroFunc1(zeroCrossings,simCode)%>
    <%isODE(simCode)%>
    <%DimZeroFunc(simCode)%>
+   <%DimInitEquations(simCode)%>
+   <%DimUnfixedStates(simCode)%>
+   <%DimUnfixedParameters(simCode)%>
+   <%DimIntialResiduals(simCode)%>
+   <%GetIntialStatus(simCode)%>
+   <%SetIntialStatus(simCode)%>
    <%checkConditions(zeroCrossings,whenClauses,simCode)%>
    <%handleSystemEvents(zeroCrossings,whenClauses,simCode)%>
    <%saveall(modelInfo,simCode)%>
@@ -1742,7 +1748,7 @@ template generateClassDeclarationCode(SimCode simCode)
 match simCode
 case SIMCODE(modelInfo = MODELINFO(__)) then
   <<
-  class <%lastIdentOfPath(modelInfo.name)%>: public IDAESystem ,public IContinous ,public IEvent ,public ISystemProperties <%if modelInfo.labels then ',public IReduceDAE '%>,public SystemDefaultImplementation
+  class <%lastIdentOfPath(modelInfo.name)%>: public IDAESystem ,public IContinous ,public IEvent ,public ISystemProperties, public  ISystemInitialization <%if modelInfo.labels then ',public IReduceDAE '%>,public SystemDefaultImplementation
   { 
   public: 
       <%lastIdentOfPath(modelInfo.name)%>(IGlobalSettings& globalSettings); 
@@ -1752,8 +1758,7 @@ case SIMCODE(modelInfo = MODELINFO(__)) then
   
   private:    
     //Methods:
-	//Saves all variables before an event is handled, is needed for the pre, edge and change operator
-     void saveAll();
+	
 
 	 void resetHelpVar(const int index);
 	 
@@ -2039,9 +2044,7 @@ case SIMCODE(modelInfo = MODELINFO(__)) then
 	virtual int getDimVars(const INDEX index = ALL_VARS) const;
 	//Provide number (dimension) of right hand sides (equations and/or residuals) according to the index 
 	virtual int getDimRHS(const INDEX index = ALL_VARS)const;
-    //(Re-) initialize the system of equations
-	virtual void init(double ts,double te);
-	//Resets all time events
+  	//Resets all time events
 	virtual void resetTimeEvents();
 	//Set current integration time 
     virtual void setTime(const double& t);
@@ -2079,6 +2082,8 @@ case SIMCODE(modelInfo = MODELINFO(__)) then
 	virtual void handleEvent(unsigned long index);
 	//Checks if a discrete variable has changed and triggers an event 
 	virtual bool checkForDiscreteEvents();
+	 //Saves all variables before an event is handled, is needed for the pre, edge and change operator
+	virtual void saveAll();
 	//Returns the vector with all time events 
 	virtual event_times_type getTimeEvents();
 	// No input
@@ -2097,6 +2102,16 @@ case SIMCODE(modelInfo = MODELINFO(__)) then
 	virtual bool hasStateDependentMass(); 
 	// System is able to provide the Jacobian symbolically 
 	virtual bool provideSymbolicJacobian();
+	//Systeminitialization methods
+	virtual unsigned int getDimInitEquations() ; 
+	virtual unsigned int getDimUnfixedStates() ; 
+	virtual unsigned int getDimUnfixedParameters(); 
+	virtual unsigned int getDimIntialResiduals() ;
+	virtual bool initial();
+	virtual void setInitial(bool);
+	/// (Re-) initialize the system of equations and bounded parameters
+	virtual void init(double ts,double te);
+	
 	<%if modelInfo.labels then 
 	<<
 	// Returns the history object to query simulation results 
@@ -4880,6 +4895,88 @@ end DimZeroFunc;
 
 
 
+
+template DimInitEquations(SimCode simCode)
+::=
+match simCode
+case SIMCODE(modelInfo = MODELINFO(__)) then
+match modelInfo
+case MODELINFO(varInfo=VARINFO(__)) then
+  <<
+   unsigned int <%lastIdentOfPath(modelInfo.name)%>::getDimInitEquations()
+    {
+      return <%varInfo.numInitEquations%>;
+    }
+  >>
+  end match
+end DimInitEquations;
+
+template DimUnfixedStates(SimCode simCode)
+::=
+match simCode
+case SIMCODE(modelInfo = MODELINFO(__)) then
+  <<
+   unsigned int <%lastIdentOfPath(modelInfo.name)%>::getDimUnfixedStates()
+    {
+      return 0;
+    }
+  >>
+end DimUnfixedStates;
+
+
+template DimUnfixedParameters(SimCode simCode)
+::=
+match simCode
+case SIMCODE(modelInfo = MODELINFO(__)) then
+match modelInfo
+case MODELINFO(varInfo=VARINFO(__)) then
+  <<
+   unsigned int <%lastIdentOfPath(modelInfo.name)%>::getDimUnfixedParameters()
+    {
+      return <%varInfo.numResiduals%>;
+    }
+  >>
+  end match
+end DimUnfixedParameters;
+
+
+template DimIntialResiduals(SimCode simCode)
+::=
+match simCode
+case SIMCODE(modelInfo = MODELINFO(__)) then
+  <<
+   unsigned int <%lastIdentOfPath(modelInfo.name)%>::getDimIntialResiduals()
+    {
+      return 0;
+    }
+  >>
+end DimIntialResiduals;
+
+template SetIntialStatus(SimCode simCode)
+::=
+match simCode
+case SIMCODE(modelInfo = MODELINFO(__)) then
+  <<
+   void <%lastIdentOfPath(modelInfo.name)%>::setInitial(bool status)
+    {
+      _initial = status;
+    }
+  >>
+end SetIntialStatus;
+
+template GetIntialStatus(SimCode simCode)
+::=
+match simCode
+case SIMCODE(modelInfo = MODELINFO(__)) then
+  <<
+   bool <%lastIdentOfPath(modelInfo.name)%>::initial()
+    {
+      return _initial;
+    }
+  >>
+end GetIntialStatus;
+
+
 template daeExpRelation(Operator op, Integer index,Exp exp1, Exp exp2, Context context, Text &preExp,Text &varDecls,SimCode simCode) ::=
   let e1 = daeExp(exp1, context, &preExp, &varDecls,simCode)
   let e2 = daeExp(exp2, context, &preExp, &varDecls,simCode)
@@ -5629,8 +5726,7 @@ template handleSystemEvents(list<ZeroCrossing> zeroCrossings,list<SimWhenClause>
     int iter=0;
     while(restart && !(iter++ > 10))
     {
-    	//save all variables for pre and edge operators
-	saveAll();
+    
     <%zeroCrossingsCode%>
     double h[<%helpvarlength(simCode)%>];
     <%helpvarvector(whenClauses,simCode)%>
@@ -5751,7 +5847,7 @@ end conditionvarZero1;
 template saveconditionvar(list<ZeroCrossing> zeroCrossings,SimCode simCode)
 ::=
   (zeroCrossings |> ZERO_CROSSING(__) hasindex i0 =>
-    conditionvarZero1(i0, relation_, simCode)
+    saveconditionvar1(i0, relation_, simCode)
   ;separator="\n")
 end saveconditionvar;
 

SimulationRuntime/cpp/Source/SettingsFactory/Implementation/GlobalSettings.cpp

@@ -5,8 +5,8 @@
 
 GlobalSettings::GlobalSettings()
 : _startTime			(0.0)
-, _endTime			(1)
-, _hOutput			(1e-4)
+, _endTime			(3)
+, _hOutput			(1e-3)
 , _resultsOutput		(true)	
 , _infoOutput		(true)
 , selected_solver("Euler")

SimulationRuntime/cpp/Source/Solver/Implementation/CMakeLists.txt

@@ -2,7 +2,7 @@ cmake_minimum_required (VERSION 2.6)
 
 project(Solver)
 # add the solver default implementation library
-add_library(Solver SHARED SolverDefaultImplementation.cpp SolverSettings.cpp  )
+add_library(Solver SHARED SolverDefaultImplementation.cpp SolverSettings.cpp Initialization.cpp )
 
 install (TARGETS Solver DESTINATION bin)
 

SimulationRuntime/cpp/Source/Solver/Implementation/Initialization.cpp

@@ -0,0 +1,38 @@
+#include "stdafx.h"
+#include "Initialization.h"
+
+
+Initialization::Initialization(ISystemInitialization* system_initialization)
+:_system(system_initialization)
+{
+}
+
+
+Initialization::~Initialization(void)
+{
+}
+void Initialization::initializeSystem(double start_time, double end_time)
+{
+	IContinous* continous_system = dynamic_cast<IContinous*>(_system);
+	IEvent* event_system = dynamic_cast<IEvent*>(_system);
+
+	//Initialization of continous equations and bounded parameters
+	_system->init(start_time,end_time);
+
+	//Intialization of discrete equations 
+	_system->setInitial(true);
+   bool restart=true;
+   int iter=0;
+
+   while(restart && !(iter++ > 10))
+   {
+	 event_system->saveAll();
+     continous_system->update(IContinous::ALL);
+	  restart = event_system->checkForDiscreteEvents();
+   }
+
+   event_system->saveAll();
+
+
+ _system->setInitial(false);
+}

SimulationRuntime/cpp/Source/Solver/Implementation/Initialization.h

@@ -0,0 +1,16 @@
+#pragma once
+#include "System/Interfaces/ISystemInitialization.h"
+#include "System/Interfaces/IEvent.h"
+#include "System/Interfaces/IContinous.h"
+class Initialization
+{
+public:
+	Initialization(ISystemInitialization* system_initialization);
+	~Initialization(void);
+	void initializeSystem(double start_time, double end_time);
+private:
+
+	ISystemInitialization* _system;
+
+};
+

SimulationRuntime/cpp/Source/Solver/Implementation/SolverDefaultImplementation.cpp

@@ -40,6 +40,7 @@ SolverDefaultImplementation::SolverDefaultImplementation(IDAESystem* system, ISo
 , _outputCommand		(IDAESystem::WRITE)
 
 {
+	_initialization = new Initialization(dynamic_cast<ISystemInitialization*>(_system));
 }
 SolverDefaultImplementation::~SolverDefaultImplementation()
 {
@@ -51,6 +52,7 @@ SolverDefaultImplementation::~SolverDefaultImplementation()
 		delete [] _zeroValLastSuccess;
 	if(_events)
 		delete [] _events;
+	delete _initialization;
 }
 
 		void SolverDefaultImplementation::setStartTime(const double& t)
@@ -78,11 +80,13 @@ void SolverDefaultImplementation::init()
 {
 	IContinous* continous_system = dynamic_cast<IContinous*>(_system);
 	IEvent* event_system =  dynamic_cast<IEvent*>(_system);
+	ISystemInitialization* init_system = dynamic_cast<ISystemInitialization*>(_system);
 	// Set current start time to the system
 	continous_system->setTime(_tCurrent);
 
 	// Assemble the system 
-	continous_system->init(_settings->getGlobalSettings()->getStartTime(),_settings->getGlobalSettings()->getEndTime());
+	//init_system->init(_settings->getGlobalSettings()->getStartTime(),_settings->getGlobalSettings()->getEndTime());
+	_initialization->initializeSystem(_settings->getGlobalSettings()->getStartTime(),_settings->getGlobalSettings()->getEndTime());
 
 
 	//// Write out head line

SimulationRuntime/cpp/Source/Solver/Implementation/SolverDefaultImplementation.h

@@ -3,9 +3,10 @@
 #include "System/Interfaces/IDAESystem.h"		// System interface
 #include "System/Interfaces/IContinous.h"		
 #include "System/Interfaces/IEvent.h"
+#include "System/Interfaces/ISystemInitialization.h"
 #include "Solver/Interfaces/IDAESolver.h"				// Solver interface
 #include "Solver/Interfaces/ISolverSettings.h"			// SolverSettings interface
-
+#include "Initialization.h"
 /// typedef to hand over (callback) functions to fortran routines
 typedef int (*U_fp)(...);
 
@@ -115,7 +116,7 @@ class BOOST_EXTENSION_SOLVER_DECL SolverDefaultImplementation
 	IDAESystem::OUTPUT	
 		_outputCommand;					///< Controls the output
 
-
+	Initialization* _initialization;
 
 private:
 	/// Definition of signum function

SimulationRuntime/cpp/Source/Solver/Implementation/SolverSettings.cpp

@@ -9,7 +9,7 @@ SolverSettings::SolverSettings( IGlobalSettings* globalSettings)
 	, _hLowerLimit	(10*UROUND)
 	, _endTimeTol	(1e-6)
 	, _zeroTol		(1e-6)
-	, _zeroTimeTol	(1e-6)
+	, _zeroTimeTol	(1e-10)
 	,_zeroRatio(1.0)
 
 {	

SimulationRuntime/cpp/Source/System/Implementation/AlgLoopSolverFactory.cpp

@@ -17,7 +17,33 @@ AlgLoopSolverFactory::~AlgLoopSolverFactory()
 {
 	if(algLoop->getDimVars(IAlgLoop::REAL) > 0)
 	{
-		std::string newton_name(NEWTON_LIB );
+		//std::string newton_name("Kinsol.dll" );
+		//type_map types;
+		//if(!load_single_library(types, newton_name))
+		//	throw std::invalid_argument(" Newton library could not be loaded");
+		//std::map<std::string, factory<IAlgLoopSolver,IAlgLoop*, IKinsolSettings*> >::iterator iter;
+		//std::map<std::string, factory<IAlgLoopSolver,IAlgLoop*, IKinsolSettings*> >& Newtonfactory(types.get());
+		//std::map<std::string, factory<IKinsolSettings> >::iterator iter2;
+		//std::map<std::string, factory<IKinsolSettings> >& Newtonsettingsfactory(types.get());
+		//iter2 = Newtonsettingsfactory.find("KinsolSettings");
+		//if (iter2 ==Newtonsettingsfactory.end()) 
+		//{
+		//	throw std::invalid_argument("No such Newton Settings");
+		//}
+		//boost::shared_ptr<IKinsolSettings> algsolversetting= boost::shared_ptr<IKinsolSettings>(iter2->second.create());
+		//_algsolversettings.push_back(algsolversetting);
+		////Todo load or configure settings
+		////_algsolversettings->load("config/Newtonsettings.xml");
+		//iter = Newtonfactory.find("KinsolCall");
+		//if (iter ==Newtonfactory.end()) 
+		//{
+		//	throw std::invalid_argument("No such Newton Solver");
+		//}
+
+		//boost::shared_ptr<IAlgLoopSolver> algsolver= boost::shared_ptr<IAlgLoopSolver>(iter->second.create(algLoop,algsolversetting.get()));
+		//_algsolvers.push_back(algsolver);
+		//return algsolver;
+		std::string newton_name("Newton.dll" );
 		type_map types;
 		if(!load_single_library(types, newton_name))
 			throw std::invalid_argument(" Newton library could not be loaded");

SimulationRuntime/cpp/Source/System/Implementation/AlgLoopSolverFactory.h

@@ -3,6 +3,7 @@
 #include "System/Interfaces/IAlgLoop.h"				// Interface for algebraic loops
 #include "System/Interfaces/IAlgLoopSolver.h"		// Interface for algebraic loops
 #include "System/Newton/Interfaces/INewtonSettings.h"
+#include "System/KinSol/Interfaces/IKinSolSettings.h"
 #include "System/Interfaces/IAlgLoopSolverFactory.h"	
 /*****************************************************************************/
 /**
@@ -28,6 +29,7 @@ class AlgLoopSolverFactory : public IAlgLoopSolverFactory
 	virtual boost::shared_ptr<IAlgLoopSolver> createAlgLoopSolver(IAlgLoop* algLoop);
 
 private:
+	//std::vector<boost::shared_ptr<IKinsolSettings> > _algsolversettings;
 	std::vector<boost::shared_ptr<INewtonSettings> > _algsolversettings;
 	std::vector<boost::shared_ptr<IAlgLoopSolver> > _algsolvers;
 };

SimulationRuntime/cpp/Source/System/Implementation/SystemDefaultImplementation.cpp

@@ -30,11 +30,7 @@ void SystemDefaultImplementation::Terminate(string msg)
 {
 	throw std::runtime_error(msg);
 }
-bool SystemDefaultImplementation::initial()
-{
-	return _initial;
 
-}
 int SystemDefaultImplementation::getDimVars(const IContinous::INDEX index) const	
 {
 	int i=0;