lp_problem.c

sqp程序包。用sqp算法实现非线性约束的优化求解

/* ****************************************************************************
*/
/*                                 user functions
*/
/* ****************************************************************************
*/
#include "o8para.h"

main() {
    void donlp2(void);

    donlp2();

    exit(0);
}
/*  lp-problem               */
/*  solution is (1,2,3,0,0) degenerate! */

/* **************************************************************************** */
/*                              donlp2 standard setup                           */
/* **************************************************************************** */
void setup0(void) {
    #define  X extern
    #include "o8comm.h"
    #undef   X

    static INTEGER i,j,nuser = 5;
    static DOUBLE  xst0[] = {0., /* not used : index 0 */
                             0.,0.,0.,0.,0.};
/* this is a bad initial point since we have 6 "active" constraints */
/* name is ident of the example/user and can be set at users will           */
/* the first static char [] must be alphabetic.  40 static char [s] maximum */

    strcpy(name,"lp-problem");

    /* x is initial guess and also holds the current solution */
    /* problem dimension n = dim(x), nh = dim(h), ng = dim(g) */
    
    n      = nuser;
    nh     = 1;
    ng     = 10;
/* that means nres = 11 . coding in gunit and gconst : 
   [0] = objective function, [1]..[nh] = equality constraints ,
   [nh+1] .... [nres] = inequality constraints */

    analyt = TRUE;
    epsdif = 1.e-16;
    del0   = 1.e0;
    tau0   = 1.e4;

    for (i = 1 ; i <= n ; i++) {
        x[i] = xst0[i];
    }
    /* gunit-array, see donlpdoc.ps */

    gunit[1][0] = -1;
    gunit[2][0] = 0;
    gunit[3][0] = 0;
    gunit[1][1] = -1;
    gunit[2][1] = 0; 
    gunit[3][1] = 0;     

    for (j = 2 ; j <= 11 ; j++)
      gunit[1][j] = 1;
	
    for (j = 2 ; j <= 6 ; j++)
      {
	gunit[2][j] = j-1;   /* 1,2,3,4,5 */
	gunit[3][j] = 1;     /* lower bound */
      }
    for (j = 7 ; j <= 11 ; j++)
      {
	gunit[2][j] = j-6;  /* 1,2,3,4,5*/
	gunit[3][j] = -1;   /* upper bound */
      }
    return;
}

/* **************************************************************************** */
/*                                 special setup                                */
/* **************************************************************************** */
void setup(void) {
    #define  X extern
    #include "o8comm.h"
    #undef   X
    
    return;
}

/* **************************************************************************** */
/*  the user may add additional computations using the computed solution here   */
/* **************************************************************************** */
void solchk(void) {
    #define  X extern
    #include "o8comm.h"
    #undef   X
    #include "o8cons.h"

    return;
}

/* **************************************************************************** */
/*                               objective function                             */
/* **************************************************************************** */
void ef(DOUBLE x[],DOUBLE *fx) {
    #define  X extern
    #include "o8fuco.h"
    #undef   X

    icf = icf+1;
/* f = 2*x_1 + 3*x_2 + 4*x_3-10*x_4-5*x_5    maximize! */
    *fx = -(2.0*x[1]+3.0*x[2]+4.0*x[3]-10.0*x[4]-5.0*x[5]);
    
    return;
}

/* **************************************************************************** */
/*                          gradient of objective function                      */
/* **************************************************************************** */
void egradf(DOUBLE x[],DOUBLE gradf[]) {
    #define  X extern
    #include "o8fuco.h"
    #undef   X


    icgf = icgf+1;
/* f = 2*x_1 + 3*x_2 + 4*x_3-10*x_4-5*x_5 maximize!*/    
    gradf[1] = -2.0;
    gradf[2] = -3.0;
    gradf[3] = -4.0;
    gradf[4]=  10.0;
    gradf[5]=   5.0;

            
    return;
}

/* **************************************************************************** */
/*                compute the i-th equality constaint, value is hxi             */
/* **************************************************************************** */
void eh(INTEGER i,DOUBLE x[],DOUBLE *hxi) {
    #define  X extern
    #include "o8fuco.h"
    #undef   X
    *hxi = x[1]+x[2]+x[3]+x[4]+x[5]-6.0;
    return;

}

/* **************************************************************************** */
/*              compute the gradient of the i-th equality constraint            */
/* **************************************************************************** */
void egradh(INTEGER i,DOUBLE x[],DOUBLE gradhi[]) {
    #define  X extern
    #include "o8fuco.h"
    #undef   X
    INTEGER j;
      for ( j =1 ; j <= 5 ; j++ ) gradhi[j]=1.0;
    return;

}

/* **************************************************************************** */
/*              compute the i-th inequality constaint, bounds included          */
/* **************************************************************************** */
void eg(INTEGER i,DOUBLE x[],DOUBLE *gxi) {
    #define  X extern
    #include "o8fuco.h"
    #undef   X
    
    if ( gunit[1][i+nh] == -1 ) cres[i+nh] = cres[i+nh]+1;

    switch (i) {
    case 1:
      *gxi = x[1];
      return;
    case 2:
      *gxi = x[2];
      return;
    case 3:
      *gxi = x[3];
      return;
    case 4:
      *gxi = x[4];
      return;
    case 5:
      *gxi = x[5];
      return;
    case 6:
      *gxi = - x[1] + 1.0;
      return;
    case 7:
      *gxi = - x[2] + 2.0;
      return;
    case 8:
      *gxi = - x[3] + 3.0;
      return;
    case 9:
      *gxi = - x[4] + 4.0;
      return;
    case 10:
      *gxi = - x[5] + 5.0;
      return;
    }    


}

/* **************************************************************************** */
/*              compute the gradient of the i-th inequality constraint          */
/*          not necessary for bounds, but constant gradients must be set        */
/*                      here e.g. using dcopy from a data-field                 */
/* **************************************************************************** */
void egradg(INTEGER i,DOUBLE x[],DOUBLE gradgi[]) {
    #define  X extern
    #include "o8fuco.h"
    #undef   X
    return;
/*     only  bounds */

}
void eval_extern(INTEGER mode) {
/*  no external evalution necessary here */

    #define  X extern
    #include "o8comm.h"
    #include "o8fint.h"
    #undef   X
    #include "o8cons.h"

    return;
}