/* rpa_proj.c: functions for current minlp-project: alan.gms (Manne 1986) 1. double f_function() 2. double g_function() 3. void G_gradient() 4. void G_hessian() 5. void add_minlp() */ int RPA_COUNT=0; int ACTIVATED = FALSE; double f_function (MINLP_ptr *MINLP_PROB,DVECTOR w,DVECTOR dual_x,int k) { /* objective(s) */ int j; double fk; DVECTOR x;DVECTOR d; x = &(w[3]); d = &(w[8]); /* define nonlinear objective functions here */ fk = x[1]*(4*x[1]+3.0*x[2]- x[3])+x[2]*(3.0*x[1]+6.0*x[2]+x[3])+x[3]*(x[2]-x[1]+10.0*x[3]); for(j=0;j<4;j++) fk += 100.0*d[j]; return fk; } /* end of f_function() */ double g_function (MINLP_ptr *MINLP_PROB,DVECTOR w,DVECTOR dual_x,int k) { /* constraint(s) */ double gk; DVECTOR x; x = &(w[3]); gk = f_function (MINLP_PROB,w,dual_x,k); return gk; } /* end of g_function() */ void F_gradient (MINLP_ptr *MINLP_PROB,DVECTOR x,DVECTOR dual_x,DVECTOR grad,int n_g,int k) { minlp_gradient (&f_function,MINLP_PROB,x,dual_x,grad,n_g,k); } /* end of G_gradient() */ void G_gradient (MINLP_ptr *MINLP_PROB,DVECTOR x,DVECTOR dual_x,DVECTOR grad,int n_g,int k) { minlp_gradient (&g_function,MINLP_PROB,x,dual_x,grad,n_g,k); } /* end of G_gradient() */ void G_hessian (MINLP_ptr *MINLP_PROB,DVECTOR x,DVECTOR dual_x,DMATRIX H,int n_w,int k) { minlp_hessian (&g_function,MINLP_PROB,x,dual_x,H,n_w,k); } /* end of G_hessian() */ void add_row (int m_g,MINLP_ptr *MINLP_PROB) { printf ("the alan project is not designed for the sequential mode\r\n"); exit (-1); } /* end of add_row() */ Ralf Östermark 1
void add_minlp (RPA_ptr* ACTIVE_NODE,MINLP_ptr *MINLP_PROB) { FILE *fid; int i,j,n_i,REDUCTION=1,START=1,imax=0,print_flip; double grad[3],H[3][3],buffer; char c_buffer[5],dbg_file[20]; DVECTOR x;DVECTOR LB;DVECTOR UB; ARCHITECTURE_ptr *MESHINFO = GHA_PROB->MESHINFO; if(ACTIVATED EQ FALSE) { ACTIVATED = TRUE; x = &(MINLP_PROB->x[4]); LB = MINLP_PROB->LB; UB = MINLP_PROB->UB; n_i = MINLP_PROB->n_i; H[0][0] = 8.0;H[0][1] = 6.0;H[0][2] = -2; H[1][0] = 6.0;H[1][1] = 12.0;H[1][2] = 2; H[2][0] = -2.0;H[2][1] = 2.0;H[2][2] = 20; grad[0] = 8.0*x[0]+6.0*x[1] -2.0*x[2]; grad[1] = 6.0*x[0]+12.0*x[1]+2.0*x[2]; grad[2] = 2.0*(-x[0]+x[1]) +20.0*x[2]; /* taking a slice of the second order Taylor expansion of f(x) to c: f(x) ~ fk + taylor_2(fk) = fk + grad(fk)(x-x[k]) + (1/2)(x- x[k])'H[k](x-x[k]) near the optimum x ~ x[k] ~ x[k-1] => taylor_2(fk) ~ 0. therefore, also the quadratic term approaches zero, since f(x) is convex by assumption. thus, (x-x[k])'H[k](x-x[k]) ~ (x-x[k])'H[k](x[k]-x[k-1]). since x ~ x[k] ~ x[k-1], x'H[k]x[k] = x[k]'H[k]x ~ 0 (a necessary optimality condition) => we want to minimize c = grad(fk) + x[k]'H[k]x */ for(i=0;i<3;i++) { MINLP_PROB->c[4+i] = grad[i]; for(j=0;j<3;j++) MINLP_PROB->c[4+i] += 0.5*x[j]*H[j][i]; } /* end of i-loop */ buffer = norm (MINLP_PROB->c,7); if(buffer>0.0) {for(j=0;j<3;j++) MINLP_PROB->c[4+j] /= buffer;} /* enter node specific restrictions on LB/UB here */ if(MESHINFO->MESH_SIZE>1) { imax = node_dec2bin (MINLP_PROB,MESHINFO->NODE_INDEX,MESHINFO- >MESH_SIZE,n_i,REDUCTION,START,print_flip=FALSE,"alan"); if(imax<0) { printf ("add_minlp(): imax(%d)<0\n",imax); exit (-1); } /* end of if(imax) */ memcpy (ACTIVE_NODE->LB,LB,n_i*sizeof(double)); memcpy (ACTIVE_NODE->UB,UB,n_i*sizeof(double)); } /* end of if(MESH_SIZE) */ memset (&MINLP_PROB->x[8],'\0',4*sizeof(double)); /* zero the deviations */ memcpy (ACTIVE_NODE->x,MINLP_PROB->x,MINLP_PROB->n*sizeof(double)); Ralf Östermark 2
} /* end of if(ACTIVATED) */ if(MINLP_PROB->debug_print EQ TRUE) { strcpy (dbg_file,"debug_alan"); sprintf (c_buffer,".%d",MESHINFO->NODE_INDEX); strcat (dbg_file,c_buffer); fid= fopen (dbg_file,"w"); fprintf (fid,"add_minlp(DEBUG_OPTLOC=TRUE) G_INCLUSION=%d TREE_COUNT=%3.0f RPA_COUNT=%d:\r\n", MINLP_PROB->G_INCLUSION,MINLP_PROB->TREE_COUNTER,RPA_COUNT); fprintf (fid,"c:\r\n"); for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",MINLP_PROB- >c[j]); fprintf (fid,"\r\n"); fprintf (fid,"A, b, IRTYPE & Lagrange coefficients:\r\n"); for(i=0;i<MINLP_PROB->m;i++) { for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",MINLP_PROB- >A[i][j]); fprintf (fid,"%15.7f %d %9.3f\r\n",MINLP_PROB->b[i],MINLP_PROB- >IRTYPE[i],MINLP_PROB->dual_x[i]); } /* end of i-loop */ for(i=0;i<MINLP_PROB->m;i++) fprintf (fid," %d",MINLP_PROB->IRTYPE[i]); fprintf (fid,"\r\nMINLP_PROB->LB:\r\n"); for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",MINLP_PROB->LB[j]); fprintf (fid,"\r\nMINLP_PROB->UB:\r\n"); for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",MINLP_PROB->UB[j]); fprintf (fid,"\r\nACTIVE_NODE->LB:\r\n"); for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",ACTIVE_NODE->LB[j]); fprintf (fid,"\r\nACTIVE_NODE->UB:\r\n"); for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",ACTIVE_NODE->UB[j]); fprintf (fid,"\r\nMINLP_PROB->x:\r\n"); for(j=0;j<MINLP_PROB->n;j++) fprintf (fid," %6.3f",MINLP_PROB->x[j]); fclose (fid); } /* end of if(DEBUG_OPTLOC) */ } /* end of add_minlp() */ Ralf Östermark 3
Recommend
More recommend