Skip to content

Commit

Permalink
- support constraint block (optimica)
Browse files Browse the repository at this point in the history 
- add example 



git-svn-id: https://openmodelica.org/svn/OpenModelica/trunk@18716 f25d12d1-65f4-0310-ae8a-bbce733d8d8e
  • Loading branch information
Vitalij Ruge committed Jan 22, 2014
1 parent 9635f1c commit 41253bc
Show file tree
Hide file tree
Showing 5 changed files with 51 additions and 30 deletions.
29 changes: 17 additions & 12 deletions SimulationRuntime/c/optimization/constraints/evalfDiffG.c
View file
Expand Up @@ -68,20 +68,22 @@ Bool evalfDiffG(Index n, double * x, Bool new_x, Index m, Index njac, Index *iRo
/* jac_struc(iRow, iCol, iData->nx, iData->nv, iData->nsi); */
jac_struc(iData, iRow, iCol);

/*
/*
printf("\n m = %i , %i",m ,iData->NRes);
printf("\nk = %i , %i" ,k ,njac);
for(i = 0; i< njac; ++i)
printf("\nJ(%i,%i) = 1; i= %i",iRow[i]+1, iCol[i]+1,i);
assert(0);
*/
*/

}
else
{
ipoptDebuge(iData,x);
long double tmp[3];
int id;
int nng = iData->nx + iData->nc;
double *a1, *a2, *a3;
double *d1, *d2, *d3;

Expand Down Expand Up @@ -118,10 +120,13 @@ Bool evalfDiffG(Index n, double * x, Bool new_x, Index m, Index njac, Index *iRo
break;
}
}
for(;j<iData->nx + iData->nc; ++j)
for(;j<nng; ++j)
{
conJac(iData->J[j], values, iData->nv, &k, j,iData);
}
}
/* printf("\nk = %i , %i" ,k ,njac); */


}

for(; i<iData->nsi; ++i)
Expand All @@ -144,13 +149,12 @@ Bool evalfDiffG(Index n, double * x, Bool new_x, Index m, Index njac, Index *iRo
break;
}
}
for(;j<iData->nx + iData->nc; ++j)
for(;j<nng; ++j)
conJac(iData->J[j], values, iData->nv, &k, j,iData);
}
}
/*assert(k == njac);*/
}

/* assert(k == njac); */
return TRUE;
}

Expand Down Expand Up @@ -481,7 +485,7 @@ static int jac_struc(IPOPT_DATA_ *iData, int *iRow, int *iCol)
{
if(iData->knowedJ[j][l] == 1)
{
iRow[k] = iRow[k-1];
iRow[k] = r + j;
iCol[k++] = c + l;
}
}
Expand Down Expand Up @@ -546,10 +550,11 @@ static int jac_struc(IPOPT_DATA_ *iData, int *iRow, int *iCol)
c += nv;
}

/*
* printf("\n\n%i = %i",iData->njac,k);
* assert(0);
*/
/*
printf("\n\n%i = %i",iData->njac,k);
assert(0);
*/


return 0;
}
Expand Down
6 changes: 5 additions & 1 deletion SimulationRuntime/c/optimization/constraints/evalfG.c
View file
Expand Up @@ -84,15 +84,17 @@ Bool evalfG(Index n, double * v, Bool new_x, int m, Number *g, void * useData)
iData->u2 = x2 + iData->nx; /*6 + 2 = 8*/
iData->u3 = x3 + iData->nx;

if(i > 0)
if(i)
{
/*1*/
functionODE_(x1, iData->u1, iData->time[i*iData->deg + 1], dotx1, iData);
for(j=0; j<iData->nx; ++j)
{
g[k++] = (a1[0]*x0[j] + a1[3]*x3[j] + scaldf[j]*dt[i]*dotx1[j]) - (a1[1]*x1[j] + a1[2]*x2[j]);
}
iData->data->callback->pathConstraints(data,&g[k],&iData->nc);

/*2*/
k += iData->nc;
functionODE_(x2, iData->u2, iData->time[i*iData->deg + 2], dotx2, iData);

Expand All @@ -102,6 +104,7 @@ Bool evalfG(Index n, double * v, Bool new_x, int m, Number *g, void * useData)
}
iData->data->callback->pathConstraints(data,&g[k],&iData->nc);

/*3*/
k += iData->nc;
functionODE_(x3, iData->u3, iData->time[i*iData->deg + 3], dotx3, iData);

Expand All @@ -110,6 +113,7 @@ Bool evalfG(Index n, double * v, Bool new_x, int m, Number *g, void * useData)
g[k++] = (a3[0]*x0[j] + a3[2]*x2[j] + scaldf[j]*dt[i]*dotx3[j]) - (a3[1]*x1[j] + a3[3]*x3[j]);
}
iData->data->callback->pathConstraints(data,&g[k],&iData->nc);
k += iData->nc;
}
else
{
Expand Down
32 changes: 21 additions & 11 deletions SimulationRuntime/c/optimization/initialOptimizer/allocate_ipopt.c
View file
Expand Up @@ -43,7 +43,7 @@

#include "../localFunction.h"

static int local_jac_struct(IPOPT_DATA_ *iData);
static int local_jac_struct(IPOPT_DATA_ *iData, int *nng);
static int check_nominal(IPOPT_DATA_ *iData, double min, double max, double nominal, short set, int i, double x0);

/*!
Expand All @@ -55,6 +55,7 @@ int allocateIpoptData(IPOPT_DATA_ *iData)
int deg1, deg2, i, j;
int nJ = iData->nc + iData->nx;
int ng = iData->NRes+iData->nc*iData->deg*iData->nsi;
int nng;
deg1 = iData->deg + 1;
deg2 = deg1 + 1;

Expand Down Expand Up @@ -129,11 +130,19 @@ int allocateIpoptData(IPOPT_DATA_ *iData)
for(i = 0; i < iData->nv; i++)
iData->mH[i] = (long double*) calloc(iData->nv, sizeof(long double));

if(iData->nc > 0)
for(i = iData->nx; i<ng-iData->nc; i+=nJ)
for(j=0;j<iData->nc;++j)
iData->gmin[i+j] = -1e32;
nng = ng-iData->nc;
if((int)iData->nc > (int)0){
for(i = iData->nx; i<nng; i+=nJ)
for(j=0;j<(int)iData->nc;++j)
{
iData->gmin[i+j] = -1e320;
}
}

/*
for(i = 0; i<ng;++i)
printf("gmin = %g \t gmax = %g\n",iData->gmin[i],iData->gmax[i]);
*/
return 0;
}

Expand Down Expand Up @@ -247,8 +256,8 @@ int loadDAEmodel(DATA *data, IPOPT_DATA_ *iData)
iData->endN = iData->NV - iData->nv;

/* iData->njac = iData->nX*iData->nsi*(iData->nv + iData->deg) + iData->nX*(iData->nv-1); */
local_jac_struct(iData);
iData->njac = iData->deg*(iData->nlocalJac-iData->nx+iData->nsi*iData->nlocalJac+iData->deg*iData->nsi*iData->nx) - iData->nc*iData->nv;
local_jac_struct(iData, &id);
iData->njac = iData->deg*(iData->nlocalJac-iData->nx+iData->nsi*iData->nlocalJac+iData->deg*iData->nsi*iData->nx)-iData->nc*iData->deg*id;
iData->nhess = 0.5*iData->nv*(iData->nv + 1)*(1+iData->deg*iData->nsi);

allocateIpoptData(iData);
Expand Down Expand Up @@ -438,16 +447,16 @@ int move_grid(IPOPT_DATA_ *iData)
* function calculates a jacobian matrix struct
* author: Willi
*/
static int local_jac_struct(IPOPT_DATA_ *iData)
static int local_jac_struct(IPOPT_DATA_ *iData, int * nng)
{
DATA * data = iData->data;
const int index = 2;
int **J;
short **Hg;

int i,j,l,ii,nx;
int i,j,l,ii,nx, id;
int *cC,*lindex;
int nJ = (iData->nx+iData->nc);
id = 0;

iData->nlocalJac = 0;
iData->knowedJ = (int**) malloc( nJ* sizeof(int*));
Expand Down Expand Up @@ -489,6 +498,7 @@ static int local_jac_struct(IPOPT_DATA_ *iData)
l = data->simulationInfo.analyticJacobians[index].sparsePattern.index[j];
J[l][ii] = 1;
++iData->nlocalJac;
if(l>= iData->nx) id++;
}
}
}
Expand Down Expand Up @@ -537,7 +547,7 @@ static int local_jac_struct(IPOPT_DATA_ *iData)
printf("\n");
}
}

*nng = id;
return 0;
}

Expand Down
8 changes: 5 additions & 3 deletions SimulationRuntime/c/optimization/lagrangeFun/ipopt_hessian.c
View file
Expand Up @@ -135,6 +135,8 @@ Bool ipopt_h(int n, double *v, Bool new_x, double obj_factor, int m, double *lam
iData->sh[i] = ll[i+2*iData->nx];
}
}
for(; i< (int)iData->nx + iData->nc; ++i)
iData->sh[i] = ll[i];
num_hessian(x, iData->time[p], iData, iData->sh ,iData->lagrange,mayer_yes,obj_factor);
}

Expand Down Expand Up @@ -189,8 +191,8 @@ static int sumLagrange(IPOPT_DATA_ *iData, double * erg,int ii, int i, int j, in

sum = 0.0;
if(iData->Hg[j][i])
for(l = 0; l<iData->nx; ++l)
sum += iData->H[l][i][j];
for(l = 0; l<iData->nx; ++l)
sum += iData->H[l][i][j];

if(iData->lagrange){
if(ii)
Expand Down Expand Up @@ -238,7 +240,7 @@ static int num_hessian(double *v, double t, IPOPT_DATA_ *iData, double *lambda,
for(j = i; j < iData->nv; ++j)
{
if(iData->Hg[i][j]){
for(l = 0; l< iData->nx; ++l)
for(l = 0; l< (int)iData->nx + iData->nc; ++l)
{
if(iData->knowedJ[l][j] + iData->knowedJ[l][i] >= 2)
iData->H[l][i][j] = (long double)lambda[l]*(iData->J[l][j] - iData->J0[l][j])/h;
Expand Down
6 changes: 3 additions & 3 deletions SimulationRuntime/c/optimization/mainOptimizer/ipoptODE.c
View file
Expand Up @@ -103,7 +103,7 @@ int startIpopt(DATA* data, SOLVER_INFO* solverInfo, int flag)
if(flag == 5)
{
nlp = CreateIpoptProblem(iData->NV, iData->Vmin, iData->Vmax,
iData->NRes, iData->gmin, iData->gmax, iData->njac, iData->nhess, 0, &evalfF,
iData->NRes +iData->nc*iData->deg*iData->nsi, iData->gmin, iData->gmax, iData->njac, iData->nhess, 0, &evalfF,
&evalfG, &evalfDiffF, &evalfDiffG, &ipopt_h);

AddIpoptNumOption(nlp, "tol", iData->data->simulationInfo.tolerance);
Expand Down Expand Up @@ -137,13 +137,13 @@ int startIpopt(DATA* data, SOLVER_INFO* solverInfo, int flag)
/* AddIpoptNumOption(nlp,"derivative_test_perturbation",1e-6); */
AddIpoptIntOption(nlp, "max_iter", 5000);

res = IpoptSolve(nlp, (*iData).v, NULL, &obj, (*iData).mult_g, (*iData).mult_x_L, (*iData).mult_x_U, (void*)iData);
res = IpoptSolve(nlp, iData->v, NULL, &obj, iData->mult_g, iData->mult_x_L, iData->mult_x_U, (void*)iData);

FreeIpoptProblem(nlp);

if(ACTIVE_STREAM(LOG_IPOPT))
{
for(i =0; i<iData->nv;i++)
for(i =0; i<iData->nv;++i)
if(iData->pFile[i])
fclose(iData->pFile[i]);
if(iData->pFile)
Expand Down

0 comments on commit 41253bc

Please sign in to comment.