Skip to content

Commit

Permalink
Dynamically allocated sparse structure and generic jacobian evaluation
Browse files Browse the repository at this point in the history 
  - Changed `SPARSE_PATTERN` to be dynamically allocated
    - Updated code generation
    - Updated C runtime solvers
  - Added generic function to evaluate jacobian for ida and dassl
  - Added jacobianSymbolical.c and jacobianSymbolical.h
  • Loading branch information
@AnHeuermann @adrpo
AnHeuermann authored and adrpo committed Jul 22, 2019
1 parent 03cf86a commit 6ecd8e3
Show file tree
Hide file tree
Showing 25 changed files with 342 additions and 252 deletions.
47 changes: 27 additions & 20 deletions OMCompiler/Compiler/Template/CodegenC.tpl
View file
Expand Up @@ -2122,7 +2122,9 @@ template functionSetupLinearSystemsTemp(list<SimEqSystem> linearSystems, String
DATA *data = (DATA*) ((void**)dataIn[0]);
threadData_t *threadData = (threadData_t*) ((void**)dataIn[1]);
const int equationIndexes[2] = {1,<%ls.index%>};
<% if ls.partOfJac then 'ANALYTIC_JACOBIAN* parentJacobian = data->simulationInfo->linearSystemData[<%ls.indexLinearSystem%>].parentJacobian;'%>
<% if ls.partOfJac then
'ANALYTIC_JACOBIAN* parentJacobian = data->simulationInfo->linearSystemData[<%ls.indexLinearSystem%>].parentJacobian;'
%>
ANALYTIC_JACOBIAN* jacobian = NULL;
<%varDeclsRes%>
<% if profileAll() then 'SIM_PROF_TICK_EQ(<%ls.index%>);' %>
Expand Down Expand Up @@ -2844,7 +2846,7 @@ template generateStaticSparseData(String indexName, String systemType, list<tupl
let sizeleadindex = listLength(sparsepattern)
let colPtr = genSPCRSPtr(listLength(sparsepattern), sparsepattern, "colPtrIndex")
let rowIndex = genSPCRSRows(lengthListElements(unzipSecond(sparsepattern)), sparsepattern, "rowIndex")
let colorString = genSPColors(colorList, "inSysData->sparsePattern.colorCols")
let colorString = genSPColors(colorList, "inSysData->sparsePattern->colorCols")
<<
void initializeSparsePattern<%indexName%>(<%systemType%>* inSysData)
{
Expand All @@ -2853,20 +2855,21 @@ template generateStaticSparseData(String indexName, String systemType, list<tupl
<%rowIndex%>
/* sparsity pattern available */
inSysData->isPatternAvailable = 'T';
inSysData->sparsePattern.leadindex = (unsigned int*) malloc((<%sizeleadindex%>+1)*sizeof(int));
inSysData->sparsePattern.index = (unsigned int*) malloc(<%sp_size_index%>*sizeof(int));
inSysData->sparsePattern.numberOfNoneZeros = <%sp_size_index%>;
inSysData->sparsePattern.colorCols = (unsigned int*) malloc(<%sizeleadindex%>*sizeof(int));
inSysData->sparsePattern.maxColors = <%maxColor%>;
inSysData->sparsePattern = (SPARSE_PATTERN*) malloc(sizeof(SPARSE_PATTERN));
inSysData->sparsePattern->leadindex = (unsigned int*) malloc((<%sizeleadindex%>+1)*sizeof(int));
inSysData->sparsePattern->index = (unsigned int*) malloc(<%sp_size_index%>*sizeof(int));
inSysData->sparsePattern->numberOfNoneZeros = <%sp_size_index%>;
inSysData->sparsePattern->colorCols = (unsigned int*) malloc(<%sizeleadindex%>*sizeof(int));
inSysData->sparsePattern->maxColors = <%maxColor%>;
/* write lead index of compressed sparse column */
memcpy(inSysData->sparsePattern.leadindex, colPtrIndex, (<%sizeleadindex%>+1)*sizeof(int));
memcpy(inSysData->sparsePattern->leadindex, colPtrIndex, (<%sizeleadindex%>+1)*sizeof(int));
for(i=2;i<<%sizeleadindex%>+1;++i)
inSysData->sparsePattern.leadindex[i] += inSysData->sparsePattern.leadindex[i-1];
inSysData->sparsePattern->leadindex[i] += inSysData->sparsePattern->leadindex[i-1];
/* call sparse index */
memcpy(inSysData->sparsePattern.index, rowIndex, <%sp_size_index%>*sizeof(int));
memcpy(inSysData->sparsePattern->index, rowIndex, <%sp_size_index%>*sizeof(int));
/* write color array */
<%colorString%>
Expand Down Expand Up @@ -4757,7 +4760,7 @@ match sparsepattern
let sizeleadindex = listLength(sparsepattern)
let colPtr = genSPCRSPtr(listLength(sparsepattern), sparsepattern, "colPtrIndex")
let rowIndex = genSPCRSRows(lengthListElements(unzipSecond(sparsepattern)), sparsepattern, "rowIndex")
let colorString = genSPColors(colorList, "jacobian->sparsePattern.colorCols")
let colorString = genSPColors(colorList, "jacobian->sparsePattern->colorCols")
let indexColumn = (jacobianColumn |> JAC_COLUMN(numberOfResultVars=n) => '<%n%>';separator="\n")
let tmpvarsSize = (jacobianColumn |> JAC_COLUMN(columnVars=vars) => listLength(vars);separator="\n")
let index_ = listLength(seedVars)
Expand All @@ -4777,20 +4780,21 @@ match sparsepattern
jacobian->seedVars = (modelica_real*) calloc(<%index_%>,sizeof(modelica_real));
jacobian->resultVars = (modelica_real*) calloc(<%indexColumn%>,sizeof(modelica_real));
jacobian->tmpVars = (modelica_real*) calloc(<%tmpvarsSize%>,sizeof(modelica_real));
jacobian->sparsePattern.leadindex = (unsigned int*) malloc((<%sizeleadindex%>+1)*sizeof(int));
jacobian->sparsePattern.index = (unsigned int*) malloc(<%sp_size_index%>*sizeof(int));
jacobian->sparsePattern.numberOfNoneZeros = <%sp_size_index%>;
jacobian->sparsePattern.colorCols = (unsigned int*) malloc(<%index_%>*sizeof(int));
jacobian->sparsePattern.maxColors = <%maxColor%>;
jacobian->sparsePattern = (SPARSE_PATTERN*) malloc(sizeof(SPARSE_PATTERN));
jacobian->sparsePattern->leadindex = (unsigned int*) malloc((<%sizeleadindex%>+1)*sizeof(int));
jacobian->sparsePattern->index = (unsigned int*) malloc(<%sp_size_index%>*sizeof(int));
jacobian->sparsePattern->numberOfNoneZeros = <%sp_size_index%>;
jacobian->sparsePattern->colorCols = (unsigned int*) malloc(<%index_%>*sizeof(int));
jacobian->sparsePattern->maxColors = <%maxColor%>;
/* write lead index of compressed sparse column */
memcpy(jacobian->sparsePattern.leadindex, colPtrIndex, (<%sizeleadindex%>+1)*sizeof(int));
memcpy(jacobian->sparsePattern->leadindex, colPtrIndex, (<%sizeleadindex%>+1)*sizeof(int));
for(i=2;i<<%sizeleadindex%>+1;++i)
jacobian->sparsePattern.leadindex[i] += jacobian->sparsePattern.leadindex[i-1];
jacobian->sparsePattern->leadindex[i] += jacobian->sparsePattern->leadindex[i-1];
/* call sparse index */
memcpy(jacobian->sparsePattern.index, rowIndex, <%sp_size_index%>*sizeof(int));
memcpy(jacobian->sparsePattern->index, rowIndex, <%sp_size_index%>*sizeof(int));
/* write color array */
<%colorString%>
Expand Down Expand Up @@ -5394,7 +5398,10 @@ case e as SES_LINEAR(lSystem=ls as LINEARSYSTEM(__), alternativeTearing = at) th
messageClose(LOG_DT);
}
<% if profileSome() then 'SIM_PROF_TICK_EQ(modelInfoGetEquation(&data->modelData->modelDataXml,<%ls.index%>).profileBlockIndex);' %>
<% if ls.partOfJac then 'data->simulationInfo->linearSystemData[<%ls.indexLinearSystem%>].parentJacobian = jacobian;'%>
<% if ls.partOfJac then
'data->simulationInfo->linearSystemData[<%ls.indexLinearSystem%>].parentJacobian = jacobian;'
%>

retValue = solve_linear_system(data, threadData, <%ls.indexLinearSystem%>, &aux_x[0]);

/* check if solution process was successful */
Expand Down
4 changes: 2 additions & 2 deletions OMCompiler/SimulationRuntime/c/Makefile.objs
View file
Expand Up @@ -69,11 +69,11 @@ else
SOLVER_OBJS_MINIMAL=$(SOLVER_OBJS_FMU)
endif
ifeq ($(OMC_MINIMAL_RUNTIME),)
SOLVER_OBJS=$(SOLVER_OBJS_MINIMAL) kinsolSolver$(OBJ_EXT) linearSolverKlu$(OBJ_EXT) linearSolverLis$(OBJ_EXT) linearSolverUmfpack$(OBJ_EXT) dassl$(OBJ_EXT) radau$(OBJ_EXT) sym_solver_ssc$(OBJ_EXT) nonlinearSolverNewton$(OBJ_EXT) newtonIteration$(OBJ_EXT) ida_solver$(OBJ_EXT) irksco$(OBJ_EXT) dae_mode$(OBJ_EXT)
SOLVER_OBJS=$(SOLVER_OBJS_MINIMAL) kinsolSolver$(OBJ_EXT) linearSolverKlu$(OBJ_EXT) linearSolverLis$(OBJ_EXT) linearSolverUmfpack$(OBJ_EXT) dassl$(OBJ_EXT) radau$(OBJ_EXT) sym_solver_ssc$(OBJ_EXT) nonlinearSolverNewton$(OBJ_EXT) newtonIteration$(OBJ_EXT) ida_solver$(OBJ_EXT) irksco$(OBJ_EXT) dae_mode$(OBJ_EXT) jacobianSymbolical$(OBJ_EXT)
else
SOLVER_OBJS=$(SOLVER_OBJS_MINIMAL)
endif
SOLVER_HFILES = dassl.h dae_mode.h delay.h epsilon.h events.h external_input.h fmi_events.h ida_solver.h linearSystem.h mixedSystem.h model_help.h nonlinearSystem.h nonlinearValuesList.h radau.h sym_solver_ssc.h solver_main.h stateset.h
SOLVER_HFILES = dassl.h dae_mode.h delay.h epsilon.h events.h external_input.h fmi_events.h ida_solver.h linearSystem.h mixedSystem.h model_help.h nonlinearSystem.h nonlinearValuesList.h radau.h sym_solver_ssc.h solver_main.h stateset.h jacobianSymbolical.h

INITIALIZATION_OBJS = initialization$(OBJ_EXT)
INITIALIZATION_HFILES = initialization.h
Expand Down
8 changes: 4 additions & 4 deletions OMCompiler/SimulationRuntime/c/optimization/DataManagement/DerStructure.c
View file
Expand Up @@ -225,10 +225,10 @@ static inline void local_jac_struct(DATA * data, OptDataDim * dim, OptDataStruct
h_index = s->indexABCD[index];
/******************************/
sizeCols = data->simulationInfo->analyticJacobians[h_index].sizeCols;
maxColors = data->simulationInfo->analyticJacobians[h_index].sparsePattern.maxColors + 1;
cC = (unsigned int*) data->simulationInfo->analyticJacobians[h_index].sparsePattern.colorCols;
lindex = (unsigned int*) data->simulationInfo->analyticJacobians[h_index].sparsePattern.leadindex;
pindex = data->simulationInfo->analyticJacobians[h_index].sparsePattern.index;
maxColors = data->simulationInfo->analyticJacobians[h_index].sparsePattern->maxColors + 1;
cC = (unsigned int*) data->simulationInfo->analyticJacobians[h_index].sparsePattern->colorCols;
lindex = (unsigned int*) data->simulationInfo->analyticJacobians[h_index].sparsePattern->leadindex;
pindex = data->simulationInfo->analyticJacobians[h_index].sparsePattern->index;

s->seedVec[index] = (modelica_real **)malloc((maxColors)*sizeof(modelica_real*));
/**********************/
Expand Down
16 changes: 8 additions & 8 deletions OMCompiler/SimulationRuntime/c/optimization/DataManagement/MoveData.c
View file
Expand Up @@ -820,14 +820,14 @@ void diffSynColoredOptimizerSystem(OptData *optData, modelica_real **J, const in
const int h_index = optData->s.indexABCD[index];
ANALYTIC_JACOBIAN* jacobian = &(data->simulationInfo->analyticJacobians[h_index]);
const long double * scaldt = optData->bounds.scaldt[m];
const unsigned int * const cC = jacobian->sparsePattern.colorCols;
const unsigned int * const lindex = jacobian->sparsePattern.leadindex;
const unsigned int * const cC = jacobian->sparsePattern->colorCols;
const unsigned int * const lindex = jacobian->sparsePattern->leadindex;
const int nx = jacobian->sizeCols;
const int Cmax = jacobian->sparsePattern.maxColors + 1;
const int Cmax = jacobian->sparsePattern->maxColors + 1;
const int dnx = optData->dim.nx;
const int dnxnc = optData->dim.nJ;
const modelica_real * const resultVars = jacobian->resultVars;
const unsigned int * const sPindex = jacobian->sparsePattern.index;
const unsigned int * const sPindex = jacobian->sparsePattern->index;
long double scalb = optData->bounds.scalb[m][n];

const int * index_J = (index == 3)? optData->s.indexJ3 : optData->s.indexJ2;
Expand Down Expand Up @@ -881,12 +881,12 @@ void diffSynColoredOptimizerSystemF(OptData *optData, modelica_real **J){
const int index = 4;
const int h_index = optData->s.indexABCD[index];
ANALYTIC_JACOBIAN* jacobian = &(data->simulationInfo->analyticJacobians[h_index]);
const unsigned int * const cC = jacobian->sparsePattern.colorCols;
const unsigned int * const lindex = jacobian->sparsePattern.leadindex;
const unsigned int * const cC = jacobian->sparsePattern->colorCols;
const unsigned int * const lindex = jacobian->sparsePattern->leadindex;
const int nx = jacobian->sizeCols;
const int Cmax = jacobian->sparsePattern.maxColors + 1;
const int Cmax = jacobian->sparsePattern->maxColors + 1;
const modelica_real * const resultVars = jacobian->resultVars;
const unsigned int * const sPindex = jacobian->sparsePattern.index;
const unsigned int * const sPindex = jacobian->sparsePattern->index;

modelica_real **sV = optData->s.seedVec[index];

Expand Down
1 change: 0 additions & 1 deletion OMCompiler/SimulationRuntime/c/simulation/simulation_runtime.cpp
View file
Expand Up @@ -50,7 +50,6 @@
#include <regex.h>
#endif


/* ppriv - NO_INTERACTIVE_DEPENDENCY - for simpler debugging in Visual Studio
*
*/
Expand Down

0 comments on commit 6ecd8e3

Please sign in to comment.