Skip to content

Commit

Permalink
- num jac for optimization (num gradient missing)
Browse files Browse the repository at this point in the history 
for compare num and sym jac



git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@19087 f25d12d1-65f4-0310-ae8a-bbce733d8d8e
  • Loading branch information
Vitalij Ruge committed Feb 13, 2014
1 parent 465ac1b commit c37f130
Show file tree
Hide file tree
Showing 9 changed files with 118 additions and 18 deletions.
1 change: 1 addition & 0 deletions SimulationRuntime/c/optimization/OptimizationFlags.h
View file
Expand Up @@ -40,5 +40,6 @@
#define OPTIMIZATIONFLAGS_H_

#define DF_STEP(x,s) ( (fmin(fmax(1e-4*fabs(s*x),1e-8),1e-1)))
#define DF_STEP2(x,s) ( (fmin(fmax(1e-2*fabs(s*x),1e-3),1)))

#endif /* OPTIMIZATIONFLAGS_H_ */
99 changes: 89 additions & 10 deletions SimulationRuntime/c/optimization/constraints/evalfG.c
View file
Expand Up @@ -52,6 +52,8 @@ static inline int evalG13(Number *g, IPOPT_DATA_ *iData, double *x0, int i);
static inline int evalG21(Number *g, IPOPT_DATA_ *iData, double *x0, int i);
static inline int evalG22(Number *g, IPOPT_DATA_ *iData, double *x0, int i);
static inline int evalG23(Number *g, IPOPT_DATA_ *iData, double *x0, int i);
static int diff_symColoredODE(double *v, double t, IPOPT_DATA_ *iData, double **J);
static int num_diff_symColoredODE(double *v, double t, IPOPT_DATA_ *iData, double **J);

/*!
* eval s.t.
Expand Down Expand Up @@ -145,15 +147,21 @@ int diff_functionODE(double* v, double t, IPOPT_DATA_ *iData, double **J)
int nJ = iData->nx + iData->nc;
x = v;
u = v + iData->nx;

refreshSimData(x,u,t,iData);
diff_symColoredODE(v,t,iData,J);
for(i = 0;i<iData->nv;++i){
for(j = 0; j <iData->nx; ++j)
J[j][i] *= iData->scalf[j]*iData->vnom[i];
for(; j <nJ; ++j)
J[j][i] *= iData->vnom[i];
}
if(iData->useNumJac){
num_diff_symColoredODE(v,t,iData,J);
for(i = 0;i<iData->nv;++i)
for(j = 0; j <iData->nx; ++j)
iData->numJ[j][i] *= iData->scalf[j];
}else{
refreshSimData(x,u,t,iData);
diff_symColoredODE(v,t,iData,J);
for(i = 0;i<iData->nv;++i){
for(j = 0; j <iData->nx; ++j)
J[j][i] *= iData->scalf[j]*iData->vnom[i];
for(; j <nJ; ++j)
J[j][i] *= iData->vnom[i];
}
}

/*
#ifdef JAC_ADOLC
Expand All @@ -170,6 +178,76 @@ int diff_functionODE(double* v, double t, IPOPT_DATA_ *iData, double **J)
return 0;
}

/*
* function calculates a symbolic colored jacobian matrix by
* author: Willi Braun
*/
int num_diff_symColoredODE(double *v, double t, IPOPT_DATA_ *iData, double **J)
{
DATA * data = iData->data;
const int index = 2;
double*x,*u;
SIMULATION_DATA *sData = (SIMULATION_DATA*)iData->data->localData[0];
int i,j,l,ii,nx;
int *cC,*lindex;

x = v;
u = x + iData->nx;

nx = data->simulationInfo.analyticJacobians[index].sizeCols;
cC = (int*)data->simulationInfo.analyticJacobians[index].sparsePattern.colorCols;
lindex = (int*)data->simulationInfo.analyticJacobians[index].sparsePattern.leadindex;

memcpy(iData->lhs, sData->realVars + iData->nx, sizeof(double)*iData->nx);
memcpy(iData->vsave, v, sizeof(double)*nx);

for(ii = 0; ii<nx; ++ii){
iData->eps[ii] = DF_STEP(v[ii], iData->vnom[ii]);;
}

for(i = 1; i < data->simulationInfo.analyticJacobians[index].sparsePattern.maxColors + 1; ++i){

for(ii = 0; ii<nx; ++ii)
if(cC[ii] == i){
v[ii] = iData->vsave[ii] + iData->eps[ii];
//printf("\nlv[%i] = %g\t eps[%i] = %g",ii,v[ii], ii,iData->eps[ii]);
}

functionODE_(x, u, t, iData->lhs, iData);
if(iData->nc > 0)
iData->data->callback->pathConstraints(iData->data,iData->lhs + iData->nx,&iData->nc);

for(ii = 0; ii<nx; ++ii)
if(cC[ii] == i)
{
v[ii] = iData->vsave[ii] - iData->eps[ii];
//printf("\nrv[%i] = %g\t eps[%i] = %g",ii,v[ii], ii,iData->eps[ii]);
}

functionODE_(x, u, t, iData->rhs, iData);
if(iData->nc > 0)
iData->data->callback->pathConstraints(iData->data,iData->rhs + iData->nx,&iData->nc);

memcpy(v, iData->vsave, sizeof(double)*nx);

for(ii = 0; ii < nx; ii++){
if(cC[ii] == i){
if(ii == 0) j = 0;
else j = lindex[ii-1];

for(; j<lindex[ii]; ++j){
l = data->simulationInfo.analyticJacobians[index].sparsePattern.index[j];
J[l][ii] = (iData->lhs[l] - iData->rhs[l])/(2.0*iData->eps[ii]);
}
}
}
}


return 0;
}


/*
* function calculates a symbolic colored jacobian matrix by
* author: Willi Braun
Expand Down Expand Up @@ -203,7 +281,7 @@ int diff_symColoredODE(double *v, double t, IPOPT_DATA_ *iData, double **J)
if(cC[ii] == i){
if(ii == 0) j = 0;
else j = lindex[ii-1];

for(; j<lindex[ii]; ++j){
l = data->simulationInfo.analyticJacobians[index].sparsePattern.index[j];
J[l][ii] = data->simulationInfo.analyticJacobians[index].resultVars[l];
Expand All @@ -222,6 +300,7 @@ int diff_symColoredODE(double *v, double t, IPOPT_DATA_ *iData, double **J)
return 0;
}


/*!
* helper evalfG
* author: Vitalij Ruge
Expand Down
7 changes: 5 additions & 2 deletions SimulationRuntime/c/optimization/initialOptimizer/allocate_ipopt.c
View file
Expand Up @@ -90,8 +90,8 @@ int allocateIpoptData(IPOPT_DATA_ *iData)
iData->d1 = (double*)malloc(deg2*sizeof(double));
iData->d2 = (double*)malloc(deg2*sizeof(double));
iData->d3 = (double*)malloc(deg2*sizeof(double));
iData->lhs = (double*)malloc(iData->nx*sizeof(double));
iData->rhs = (double*)malloc(iData->nx*sizeof(double));
iData->lhs = (double*)malloc((int)iData->nJ*sizeof(double));
iData->rhs = (double*)malloc((int)iData->nJ*sizeof(double));
iData->dotx0 = (double*)malloc(iData->nx*sizeof(double));
iData->dotx1 = (double*)malloc(iData->nx*sizeof(double));
iData->dotx2 = (double*)malloc(iData->nx*sizeof(double));
Expand Down Expand Up @@ -123,6 +123,9 @@ int allocateIpoptData(IPOPT_DATA_ *iData)
iData->sv = (double*)malloc(iData->nv*sizeof(double));
iData->sh = (double*)malloc(iData->nJ*sizeof(double));

iData->vsave = (double*)malloc(iData->nv*sizeof(double));
iData->eps = (double*)malloc(iData->nv*sizeof(double));

iData->J0 = (double**) malloc(iData->nJ * sizeof(double*));
for(i = 0; i < iData->nJ; i++)
iData->J0[i] = (double*) calloc(iData->nv, sizeof(double));
Expand Down
2 changes: 1 addition & 1 deletion SimulationRuntime/c/optimization/initialOptimizer/initial_guess.c
View file
Expand Up @@ -157,7 +157,7 @@ static int initial_guess_ipopt_sim(IPOPT_DATA_ *iData,SOLVER_INFO* solverInfo)
v += iData->nv;
/* updateContinuousSystem(iData->data); */
rotateRingBuffer(iData->data->simulationData, 1, (void**) iData->data->localData);
}
}
}

for(i = 0, id=0; i<iData->NV;i++,++id){
Expand Down
4 changes: 4 additions & 0 deletions SimulationRuntime/c/optimization/ipoptODEstruct.h
View file
Expand Up @@ -196,6 +196,10 @@ typedef struct IPOPT_DATA_
int matrixC;
int matrixD;

short useNumJac;
double *vsave;
double *eps;

/*Debuger*/
FILE **pFile;
long index_debug_iter;
Expand Down
3 changes: 0 additions & 3 deletions SimulationRuntime/c/optimization/localFunction.h
View file
Expand Up @@ -69,11 +69,8 @@ Bool ipopt_h(int n, double *x, Bool new_x, double obj_factor, int m, double *lam

/*sym JAC*/
int diff_functionODE(double *v, double t, IPOPT_DATA_ *iData, double **J);
int diff_symColoredODE(double *v, double t, IPOPT_DATA_ *iData, double **J);
int diff_symColoredObject(IPOPT_DATA_ *iData, double *gradF, int this_it);

/*num JAC */

/*model*/
int functionODE_(double * x, double *u, double t, double * dotx, IPOPT_DATA_ *iData);
int refreshSimData(double *x, double *u, double t, IPOPT_DATA_ *iData);
Expand Down
13 changes: 12 additions & 1 deletion SimulationRuntime/c/optimization/mainOptimizer/ipoptODE.c
View file
Expand Up @@ -226,7 +226,7 @@ static int res2file(IPOPT_DATA_ *iData,SOLVER_INFO* solverInfo)
* set optimizer options
* author: Vitalij Ruge
**/
static int set_optimizer_flags(IPOPT_DATA_ *iData,IpoptProblem *nlp)
static int set_optimizer_flags(IPOPT_DATA_ *iData, IpoptProblem *nlp)
{
char *cflags;
AddIpoptNumOption(*nlp, "tol", iData->data->simulationInfo.tolerance);
Expand All @@ -253,6 +253,17 @@ static int set_optimizer_flags(IPOPT_DATA_ *iData,IpoptProblem *nlp)
warningStreamPrint(LOG_STDOUT, 1, "not support ipopt_hesse=%s",cflags);
}

iData->useNumJac = 0;
cflags = (char*)omc_flagValue[FLAG_IPOPT_JAC];
if(cflags){
if(!strcmp(cflags,"NUM"))
iData->useNumJac = 1;
else if(!strcmp(cflags,"SYM") || !strcmp(cflags,"sym"))
iData->useNumJac = 0;
else
warningStreamPrint(LOG_STDOUT, 1, "not support ipopt_hesse=%s",cflags);
}

cflags = (char*)omc_flagValue[FLAG_LS_IPOPT];
if(cflags)
AddIpoptStrOption(*nlp, "linear_solver", cflags);
Expand Down
6 changes: 5 additions & 1 deletion SimulationRuntime/c/util/simulation_options.c
View file
Expand Up @@ -44,6 +44,7 @@ const char *FLAG_NAME[FLAG_MAX+1] = {
/* FLAG_INTERACTIVE */ "interactive",
/* FLAG_IOM */ "iom",
/* FLAG_IPOPT_HESSE*/ "ipopt_hesse",
/* FLAG_IPOPT_JAC*/ "ipopt_jac",
/* FLAG_IPOPT_INIT*/ "ipopt_init",
/* FLAG_L */ "l",
/* FLAG_LOG_FORMAT */ "logFormat",
Expand Down Expand Up @@ -78,7 +79,8 @@ const char *FLAG_DESC[FLAG_MAX+1] = {
/* FLAG_ILS */ "[int] default: 1",
/* FLAG_INTERACTIVE */ "specify interactive simulation",
/* FLAG_IOM */ "value specifies the initialization optimization method",
/* FLAG_IPOPT_HESSE */ "value specifies the hessematrix for Ipopt",
/* FLAG_IPOPT_HESSE */ "value specifies the hessian for Ipopt",
/* FLAG_IPOPT_JAC */ "value specifies the jacobian for Ipopt",
/* FLAG_IPOPT_INIT */ "value specifies the initial guess for optimization",
/* FLAG_L */ "value specifies a time where the linearization of the model should be performed",
/* FLAG_LOG_FORMAT */ "value specifies the log format of the executable. -logFormat=text (default) or -logFormat=xml",
Expand Down Expand Up @@ -114,6 +116,7 @@ const char *FLAG_DETAILED_DESC[FLAG_MAX+1] = {
/* FLAG_INTERACTIVE */ "specify interactive simulation",
/* FLAG_IOM */ "value specifies the initialization optimization method",
/* FLAG_IPOPT_HESSE */ "value specifies the hessematrix for Ipopt(OMC, BFGS, const)",
/* FLAG_IPOPT_JAC */ "value specifies the jacobian for Ipopt(SYM, NUM)",
/* FLAG_IPOPT_INIT */ "value specifies the initial guess for optimization (sim, const)",
/* FLAG_L */ "value specifies a time where the linearization of the model should be performed",
/* FLAG_LOG_FORMAT */ "value specifies the log format of the executable. -logFormat=text (default) or -logFormat=xml",
Expand Down Expand Up @@ -149,6 +152,7 @@ const int FLAG_TYPE[FLAG_MAX] = {
/* FLAG_INTERACTIVE */ FLAG_TYPE_FLAG,
/* FLAG_IOM */ FLAG_TYPE_OPTION,
/* FLAG_IPOPT_HESSE */ FLAG_TYPE_OPTION,
/* FLAG_IPOPT_JAC */ FLAG_TYPE_OPTION,
/* FLAG_IPOPT_INIT */ FLAG_TYPE_OPTION,
/* FLAG_L */ FLAG_TYPE_OPTION,
/* FLAG_LOG_FORMAT */ FLAG_TYPE_OPTION,
Expand Down
1 change: 1 addition & 0 deletions SimulationRuntime/c/util/simulation_options.h
View file
Expand Up @@ -51,6 +51,7 @@ enum _FLAG
FLAG_INTERACTIVE,
FLAG_IOM,
FLAG_IPOPT_HESSE,
FLAG_IPOPT_JAC,
FLAG_IPOPT_INIT,
FLAG_L,
FLAG_LOG_FORMAT,
Expand Down

0 comments on commit c37f130

Please sign in to comment.