dembo1as.c

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

            -1.e-6,-1.e-5,-1.e-6,-1.e-9,
            -1.e-9,-1.e-3, 0.e0 ,-1.e-3,
            -1.e-5, 0.e0 ,-1.e-4, 0.e0 };
             
    if ( i <= 3 ) cres[i+nh] = cres[i+nh]+1;
    if ( i  > 3 ) goto L700;
    
    switch (i) {
    case 1:
        *gxi = 1.e0-5.36373e-2*x[1]-2.1863746e-2*x[2]-9.7733533e-2*x[3]
                -6.6940803e-3*x[4]*x[5];
        *gxi = *gxi*1.e1;
        
        return;
        
    case 2:
    
        fgeo(x,1.e0,gam2,TRUE,dl2,k2,al2,4,9,12,gxi);
        
        *gxi = *gxi*1.e1;
        
        return;
        
    case 3:

        fgeo(x,1.e0,gam3,TRUE,dl3,k3,al3,3,5,12,gxi);
        
        *gxi = *gxi*1.e2;
        
        return;
    }
    L700:

    *gxi = 1.e2*(x[i-3]-ug[i-3]);
    
    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
    
    void dfgeo(DOUBLE x[],DOUBLE gam[],LOGICAL lin,DOUBLE dl[],void *pk,
               void *pal,INTEGER nlen,INTEGER nanz,DOUBLE g[],INTEGER nx);

    static INTEGER j;
    static DOUBLE  s,p;
    
    static DOUBLE gam2[] = { 0., /* not used : index 0 */
                            -1.e-9,-1.e-9,-1.e-3,-1.e-3,
                            -1.0898645e-1,-1.6108052e-5,-1.0e-23,-1.9304541e-8,
                            -1.e-4 };
                            
    static DOUBLE gam3[] = { 0., /* not used : index 0 */
                            -1.0898645e-1,-1.6108052e-5,-1.e-23,-1.9304541e-8,
                            -1.1184059e-4 };
                            
    static DOUBLE al2[][10] = {
            {0.e0,      0.e0, 0.e0, 0.e0, 0.e0, 0.e0, 0.e0, 0.e0, 0.e0, 0.e0},
             
            {0.e0,      1.e0, 1.e0, 1.e0, 1.e0, 1.e0, 1.e0, 1.e0, 1.e0, 1.e0},
            {0.e0,      1.e0,-1.e0,-1.e0, 1.e0, 1.e0, 1.e0, 1.e0,-1.e0,-1.e0},
            {0.e0,      0.e0, 0.e0, 0.e0, 0.e0, 0.e0,-1.e0, 1.e0, 1.e0, 0.e0},
            {0.e0,      0.e0, 0.e0, 0.e0, 0.e0, 0.e0, 0.e0, 0.e0,-2.e0, 0.e0} };
    
    static INTEGER k2[][10] = {
            {0,      0, 0, 0, 0, 0, 0, 0, 0, 0},
            
            {0,      4, 5, 6, 7, 4, 2, 2, 2,10},
            {0,     12,12,12,12, 5, 5, 4, 4,12},
            {0,      0, 0, 0, 0, 0,12, 5, 4, 0},
            {0,      0, 0, 0, 0, 0, 0, 0,12, 0} };

    static DOUBLE al3[][6] = {
            {0.e0,      0.e0, 0.e0, 0.e0, 0.e0, 0.e0},
             
            {0.e0,      1.e0, 1.e0, 1.e0, 1.e0, 1.e0},
            {0.e0,      1.e0, 1.e0, 1.e0,-1.e0, 1.e0},
            {0.e0,      0.e0, 0.e0, 1.e0, 1.e0, 0.e0} };
     
    static INTEGER k3[][6] = {
            {0,      0, 0, 0, 0, 0},
            
            {0,      4, 2, 2, 2, 1},
            {0,      5, 5, 4, 4, 9},
            {0,      0, 0, 5, 5, 0} };
    
    static DOUBLE dl2[] = { 0., /* not used : index 0 */
            -1.e-6,-1.e-5,-1.e-6,0.e0,0.e0,0.e0,0.e0,0.e0,0.e0,0.e0,0.e0,0.e0 };
            
    static DOUBLE dl3[] = { 0., /* not used : index 0 */
            -1.e-6,-1.e-5,-1.e-6,-1.e-9,
            -1.e-9,-1.e-3, 0.e0 ,-1.e-3,
            -1.e-5, 0.e0 ,-1.e-4, 0.e0 };
             
    /*   *gxi = 1.e0-5.36373e-2*x[1]-2.1863746e-2*x[2]-9.7733533e-2*x[3] */
    /*          -6.6940803e-3*x[4]*x[5]                                  */
    /*   return                                                          */
    
    if ( i > 3 ) return;
    
    cgres[i+nh] = cgres[i+nh]+1;
    
    switch (i) {
    case 1:
        for (j = 1 ; j <= 12 ; j++) {
            gradgi[j] = 0.e0;
        }
        gradgi[1] = -5.36373e-2;
        gradgi[2] = -2.1863746e-2;
        gradgi[3] = -9.7733533e-2;
        gradgi[4] = -x[5]*6.6940803e-3;
        gradgi[5] = -x[4]*6.6940803e-3;
        for (j = 1 ; j <= 12 ; j++) {
            gradgi[j] = gradgi[j]*1.e1;
        }
        return;
        
    case 2:

        dfgeo(x,gam2,TRUE,dl2,k2,al2,4,9,gradgi,12);
        
        for (j = 1 ; j <= 12 ; j++) {
            gradgi[j] = gradgi[j]*1.e1;
        }
        return;
        
    case 3:

        dfgeo(x,gam3,TRUE,dl3,k3,al3,3,5,gradgi,12);
        
        for (j = 1 ; j <= 12 ; j++) {
             gradgi[j] = gradgi[j]*1.e2;
        }
        return;
    }
}

/* **************************************************************************** */
/* evalution of a function of a geometric programming problem described by      */
/* the model                                                                    */
/* fx  =  con + sum{i=1,nx} x[i]*dl[i]                                          */
/*            + sum{i=1,nanz} (gam[i]*(prod{j=1,nlen}pow(x[k[j][i]],al[j][i]))  */
/* if lin  =  FALSE dl may be undefined                                         */
/* **************************************************************************** */
void fgeo(DOUBLE x[],DOUBLE con,DOUBLE gam[],LOGICAL lin,DOUBLE dl[],void *pk,
          void *pal,INTEGER nlen,INTEGER nanz,INTEGER nx,DOUBLE *fx) {
          
    static INTEGER i,j,il,*k;
    static DOUBLE  s,p,expo,*al;
    
    k  = pk;
    al = pal;
    
    s = con;
    
    if ( ! lin ) goto L200;
    
    for (i = 1 ; i <= nx ; i++) {
        s = s+dl[i]*x[i];
    }
    L200 :
  
    for (i = 1 ; i <= nanz ; i++) {
    
        if ( gam[i] == 0.e0 ) goto L600;
        
        p = 1.e0;
        for (j = 1 ; j <= nlen ; j++) {
            il = k[j*(nanz+1)+i];           /* il = k[j][i]; */
    
            if ( il == 0 ) goto L500;
    
            expo = al[j*(nanz+1)+i];        /* expo = al[j][i]; */
    
            if ( expo == 0.e0 ) goto L500;
    
            p = p*exp(expo*log(fabs(x[il])));
    
            L500:;
        }

        s = s+gam[i]*p;
    
        L600:;
    }

    *fx = s;
    
    return;
}

/* **************************************************************************** */
/*          computation of the gradient of a function given by fgeo above       */
/* **************************************************************************** */
void dfgeo(DOUBLE x[],DOUBLE gam[],LOGICAL lin,DOUBLE dl[],void *pk,
           void *pal,INTEGER nlen,INTEGER nanz,DOUBLE g[],INTEGER nx) {
           
    static INTEGER i,j,il,l,*k;
    static DOUBLE  s,p,fc,expo,fij,*al;
    
    k  = pk;
    al = pal;
    
    for (l = 1 ; l <= nx ; l++) {
        s = 0.e0;
        for (i = 1 ; i <= nanz ; i++) {
        
            if ( gam[i] == 0.e0 ) goto L400;
            
            p  = 1.e0;
            fc = 0.e0;
            for (j = 1 ; j <= nlen ; j++) {
                il = k[j*(nanz+1)+i];           /* il = k[j][i]; */
                if( il == 0 ) goto L300;
                if( il != l ) goto L100;
                fc = 1.e0;
                
                L100:

                expo = al[j*(nanz+1)+i];        /* expo = al[j][i]; */
                if ( expo == 0.e0 ) goto L300;
                fij = 1.e0;
                if ( il != l ) goto L200;
                fij  = expo;
                expo = expo-1.e0;
                
                L200:

                p = p*fij*exp(expo*log(fabs(x[il])));
                
                L300:;
            }
            if ( fc != 0.e0 ) s = s+p*gam[i];
                
            L400:;
        }
        if ( lin ) s = s+dl[l];
        g[l] = s;
    }
    return;
}

/* **************************************************************************** */
/*                        user functions (if bloc == TRUE)                      */
/* **************************************************************************** */
void eval_extern(INTEGER mode) {
    #define  X extern
    #include "o8comm.h"
    #include "o8fint.h"
    #undef   X
    #include "o8cons.h"

    return;
}