User's Guide
for
LSGRG2 Optimization Library

By:

Windward Technologies, Inc.
and
Optimal Methods

 

 

Contents

 

1. Introduction

1.1 LSGRG2 Overview

2. Essential Steps to Using LSGRG2

2.1 Example Problem

2.2 Program for Example Problem

3. How to Use LSGRG2

3.1 LSGRG2 C Function Interface     

3.1.1 Purpose of lgoptG and lsgrg2

3.1.2 Prototypes for lgoptG and lsgrg2

3.1.3 Arguments of lsgrg2

3.1.4 Termination Messages (inform)

3.2 Changing Algorithmic Parameters of LSGRG2

3.2.1 Prototype for lgoptD

3.2.2 Descriptions of LSGRG2 Algorithmic Parameters

3.3 Providing Derivatives to LSGRG2

3.3.1 Providing Derivatives Through lgoptJ

3.3.2 Providing Derivatives Through lgoptP

3.4 Specifying Linear Variables Through lgoptL

4. Things to Try If LSGRG2 Isn't Working

4.1 Errors in the Input Data or in GCOMP

4.2 Incorrect User-Supplied Derivatives

4.3 Inaccurate Numerical Derivatives

4.4 Scaling

4.5 Bounds Which Must Not be Violated

4.6 Local and Global Optima

4.7 Solution Not Accurate Enough

5. LSGRG2 Options and Tolerances

5.1 Tolerances

5.2 Algorithmic Options

6. Unconstrained Minimization and Optimal Control

7. References


1. Introduction

LSGRG2 is a program that solves nonlinear optimization problems of the following form:

Minimize or maximize gp(X),

subject to:

glbi £ gi(X) £ gubi

for i = 0,...,m-1, i ¹ p

xlbj £ xj £ xubj

for j = 0,...,n-1

X is a vector on n variables, x0 ,...,xn-1, and the functions g0 ,...,gm-1 all depend on X.

Any of these functions may be nonlinear. Any of the bounds may be infinite and any of the constraints may be absent. If there are no constraints, the problem is solved as an unconstrained optimization problem. Upper and lower bounds on the variables are optional and, if present, are not treated as additional constraints but are handled separately. The program attempts to solve problems of this form by the Generalized Reduced Gradient Methods (see references 1-4 for a description and further references).

This User’s Guide describes the C interface for LSGRG2. The C version of LSGRG2 is specially designed to be used as part of a larger system. The interface uses features of the C language, such as pointers, to make it easy to build LSGRG2 into an application with a front-end written in C or written in another language (e.g., Visual Basic). These features are explained later in this guide.

LSGRG2 uses first partial derivatives of each function gi with respect to each variable xj. These are automatically computed by finite difference approximation (either forward or central differences) unless the user supplies a function that evaluates them using formulas or other methods. After an initial data entry segment, the program operates in two phases. If the initial values of the variables supplied by the user do not satisfy all gi constraints, a Phase I optimization is started. The Phase I objective function is the sum of the constraint violations plus, optionally, a fraction of the true objective. This optimization terminates either with a message that the problem is infeasible or with a feasible solution. Beware if an infeasibility message is produced, because the program may have become stuck at a local minimum of the Phase I objective (or too large a part of the true objective was incorporated), and the problem may actually have feasible solutions. The suggested remedy, in this case, is to choose different starting values for the variables (or reduce the proportion of the true objective) and try again. These remedies are described in Section 4. "Things to Try If GRG2 Isn’t Working".

Phase II begins with a feasible solution, either found by Phase I or with the user provided starting point if it is feasible, and attempts to optimize the objective function. At the conclusion of Phase II, a full optimization cycle has been completed and summary output is provided.

1.1 LSGRG2 Overview

There are a few essential steps that must be taken in order to present an optimization problem to LSGRG2. First, a C function must be written to compute values of the constraint functions and the objective function for given values of the variables. This C function is referred to as GCOMP (the name of its pointer in LSGRG2) in this manual. Next, the data that describe the optimization problem must be created, organized, and passed to LSGRG2. The data consist of the number of variables, lower and upper bounds on the variables, the number of functions (number of constraint functions plus one for the objective function), the index of the objective function in this set, and lower and upper bounds on the constraints. LSGRG2 also needs to be given the name of a file that it can use to prepare a report on the optimization problem and a title (character string) to identify the problem in the report. Finally, LSGRG2 needs a variable (argument inform) that it can use to return an integer code that reflects the outcome of the optimization process. An example is given below to illustrate how a simple optimization problem is prepared and presented to LSGRG2. If LSGRG2 returns a value of inform=0 or inform=1 this is an indication that the problem has been successfully solved. Other values of inform indicate an outcome that may not be successful. The report should be reviewed in any case. The report contains a summary of starting values and final values for variables, constraints, and the objective function for successful runs. It contains error messages and other information that can be helpful for runs that are not successful.

Beyond these basics, LSGRG2 has features to deal with the realities of large optimization problem solving. Some of these features reduce the time needed to solve large problems and other features make it possible to solve difficult problems. LSGRG2 is capable of solving problems with several thousand variables and several thousand constraints.

One important feature that helps LSGRG2 deal with large problems is its ability to represent the matrix of partial derivatives of the constraint functions and the objective function efficiently. LSGRG2 uses a sparse matrix representation (nonzero elements packed by columns) to deal with sparse problems. It is common for large problems to have a sparse matrix of partial derivatives. This happens when each constraint involves only a few variables. A large portion of the time needed to solve large sparse problems is taken up by the automatic computation of partial derivatives. It is possible to replace the automatic derivative computation with user supplied computations. The user-supplied computations can take advantage of the sparsity structure to significantly reduce the computation time. For example, if some constraints are linear functions of some variables then the derivatives of these constraints with respect to the linear variables are constant and don't involve any computation. LSGRG2 provides two formats for user supplied derivatives. These features are described in Section 3.3.1 "Providing Derivatives Through lgoptJ", and Section 3.3.2 "Providing Derivatives Through lgoptP"

LSGRG2 can take advantage of linearity. When each constraint and the objective function are linear in a given variable, derivatives with respect to this variable are constant and it is not necessary to compute them more than once. The linear variables can be specified to LSGRG2 so that it can take advantage of the linearity. See the description of this feature in Section 3.4 "Specifying Linear Variables Through lgoptL".

There are many parameters that control the behavior of the LSGRG2 algorithm. Each parameter has a default value that is appropriate for most problems. The user is not required to take any action in order to use the default parameter values. At times, it may be necessary to set one or more of the parameters to a new value in order to make LSGRG2 more efficient or in order to make it possible to solve a difficult problem. See the description in Section 3.2.2 "Descriptions of LSGRG2 Algorithmic Parameters".

The following table summarizes the user callable functions in the LSGRG2 package and gives the section where they are documented.

 

Function

Description

lgoptG

Required:  provides the name of the function that evaluates constraint functions and the objective function. See Section 2. for an example.

lsgrg2

Required:  main interface, provides required problem specifications. See Section 2. for an example.

lgoptD

Optional:  provides optional parameter values

lgoptJ

Optional:  provides derivative evaluation function, one column format.

lgoptP

Optional:  provides derivative evaluation function, packed column format

lgoptL

Optional:  provides a list of linear variables.

 


2. Essential Steps to Using LSGRG2

This section presents an example problem and a C programs to illustrate the essential steps necessary to prepare a program that uses LSGRG2.

2.1 Example Problem

This problem has 5 variables {X0, X1,..., X4}, bounded above and below, a quadratic objective function {G3}, and three quadratic constraints {G0, G1, G2}, with both upper and lower bounds. It is problem number 11 in the Himmelblau text [7].

Minimize: G3 = 5.3578547X2X2 + 0.8356891X0X4 + 37.293239X0 - 40792.141
 
Subject to:
   0 < G0 = 85.334407 + 0.0056858X1X4 + 0.0006262X0X3 - 0.0022053X2X4 <  92
  90 < G1 = 80.51249  + 0.0071317X1X4 + 0.0029955X0X1 + 0.0021813X2X2 < 110
  20 < G2 = 9.300961  + 0.0047026X2X4 + 0.0012547X0X2 + 0.0019085X2X3 <  25
and,
 
  78 < X0 < 102;  33 < X1 < 45;  27 < X2 < 45;  27 < X3 < 45;  27 < X4 < 45

This problem is solved starting from X = {78, 33, 27, 27, 27} which is infeasible because the third constraint, G2, is not satisfied.

2.2 Program for Example Problem

The main program that calls LSGRG2 appears below. The function h11G evaluates constraint functions and the objective function for the example problem. Note that h11G is passed to LSGRG2 by calling lgoptG prior to calling lsgrg2. The program has an include statement for the header file "lsgrg2.h" which contains prototypes for all LSGRG2 functions.

 
#include "lsgrg2.h"
 
TYPE_stdcall(void) h11G (double g[4], double x[5]);
TYPE_stdcall(void) h11J (double x[5], double g[4],
                  long int jcol, double gradj[4]);
TYPE_stdcall(void) h11P (double x[], double g[],
                  long int igrad[], long int jgrad[],
                  double grad[], long int nzg);
 
void main() {
 
    long int    nvars, nfuns, nobj, inform;
    double      xlb[5], xub[5], glb[4], gub[4], xx[5];
    char        *title, *report;
 
    nvars = 5 ;  /*  number of variables  */
 
    /*  set bounds on variables  */
    xlb[0] = 78.0;    xub[0] = 102.0;
    xlb[1] = 33.0;    xub[1] = 45.0;
    xlb[2] = 27.0;    xub[2] = 45.0;
    xlb[3] = 27.0;    xub[3] = 45.0;
    xlb[4] = 27.0;    xub[4] = 45.0;
 
    nfuns =  4 ;  /*  number of functions  */
    nobj  =  3 ;  /*  minimize g[3]        */
 
    /*  set bounds on constraints  */
    glb[0] = 0.0;     gub[0] = 92.0;
    glb[1] = 90.0;    gub[1] = 110.0;
    glb[2] = 20.0;    gub[2] = 25.0;
 
    /*  title for report  */
    title = "LSGRG2 C Example:   Himmelblau Problem 11.";
 
    /*  file name for LSGRG2 report  */
    report = "H11.txt";
 
    /*  starting values for variables  */
    xx[0] = 78; xx[1] = 33; xx[2] = 27; xx[3] = 27; xx[4] = 27;
 
/*  if neither of the following two lines is active
    LSGRG2 will compute derivatives automatically
    by finite difference approximations             */
/*  lgoptJ (h11J);  use h11J to compute derivatives */ 
/*  lgoptP (h11P);  use h11P to compute derivatives */
    lgoptG (h11G);  /*  GCOMP pointer  */
 
    lsgrg2 (nvars, xlb, xub, nfuns, -nobj, glb, gub,
            title, report, xx, &inform);
 
}  /*  end h11main  */
 
TYPE_stdcall(void) h11G (double  g[4], double  x[5]) {
/*  Himmelblau Problem 11; nvars=5; nfuns=4; nobj=3;  */
 
    g[0] = 85.334407 + 0.0056858*x[1]*x[4] + 0.0006262*x[0]*x[3]
                     - 0.0022053*x[2]*x[4];
    g[1] = 80.51249  + 0.0071317*x[1]*x[4] + 0.0029955*x[0]*x[1]
                     + 0.0021813*x[2]*x[2];
    g[2] = 9.300961  + 0.0047026*x[2]*x[4] + 0.0012547*x[0]*x[2]
                     + 0.0019085*x[2]*x[3];
    g[3] = 5.3578547*x[2]*x[2] + 0.8356891*x[0]*x[4] + 37.293239*x[0]
                     - 40792.141;
 
}  /*  end h11G  */

3. How to Use LSGRG2

This section describes how to use the LSGRG2 C function interface. The LSGRG2 package includes lgoptG, lsgrg2, and other functions that select options. All required information for a given problem is passed through lgoptG and lsgrg2. The user supplied C function that evaluates the constraint functions and the objective function is passed through lgoptG. All other required information, such as variable bounds and constraint bounds, is passed through lsgrg2. The purpose of this design is to provide flexibility when LSGRG2 is used as part of a larger system. Several GCOMP functions, representing different problems or variations of one problem, can be packaged together and passed to LSGRG2 one at a time through lgoptG. This design makes it possible to pass GCOMP as a component of a separate system that may be written in a language other than C or to call lsgrg2 from one system and provide GCOMP through another system.

The remainder of Section 3 describes lsgrg2 and its support functions lgoptG, lgoptD, lgoptJ, lgoptP and lgoptL. Most problems can be solved with only lgoptG and lsgrg2. Although, for very large problems, much computation time can be saved by providing an efficient C function to compute partial derivatives. The support functions lgoptJ, and lgoptP are used for this purpose. The other support functions, lgoptD and lgoptL, can used for difficult problems or to modify the behavior of LSGRG2 to take advantage of special, problem specific, conditions.

3.1 LSGRG2 C Function Interface

This section describes lsgrg2, the C interface for LSGRG2, and its support function lgoptG. Most problems can be solved with only lgoptG and lsgrg2. The prototypes that are given in this section provide specifications for the arguments of these functions. Note that the file lsgrg2.h, part of the LSGRG2 package, contains prototypes and should be used with programs that call lsgrg2. The remaining material in this section describes the arguments of lgoptG and lsgrg2.

3.1.1 Purpose of lgoptG and lsgrg2

The problem solved by LSGRG2 is restated here in terms of its C variables.

 
Minimize     g[nobj]
 
Subject to:  glb[i] <= g[i] <= gub[i], i = 0,..., nfuns-1
                                       i != nobj
             xlb[j] <= x[j] <= xub[j], j = 0,..., nvars-1

where x[j] are variables and g[i] are nonlinear functions that specify problem constraints.

3.1.2 Prototypes for lgoptG and lsgrg2

The prototypes for all functions in the LSGRG2 package are included in the header file "lsgrg2.h" and can be accessed by including this file in your program.

    TYPE_stdcall(void) lgoptG (
                 void (STDCALL *GCOMP) (
                 double g[], double x[]));
    TYPE_stdcall(void) lsgrg2 (
                 long int    nvars,
                 double      *xlb,
                 double      *xub,
                 long int    nfuns,
                 long int    nobj,
                 double      *glb,
                 double      *gub,
                 char        *title,
                 char        *report,
                 double      *xx,
                 long int    *inform);

Notes: The macros TYPE_stdcall(void) and STDCALL are defined through the header file lsgrg2.h for the target system.

To use LSGRG2: prepare a function, consistent with GCOMP above, and pass it to LSGRG2 by calling lgoptG; establish values for all LSGRG2 arguments, then call lsgrg2. The function GCOMP must evaluate g[i], i=0 to nfuns-1, for given values of variables x[j], j=0 to nvars-1. One of these functions, with index nobj, is the objective function; the others are constraints.

    lgoptG (gcomp);
    lsgrg2 (nvars, xlb, xub, nfuns, nobj, glb, gub, 
            title, report, xx, &inform);

3.1.3 Arguments of lsgrg2

 
Argument        Description
 
    nvars       Number of variables. Variables are numbered 0 to
                nvars-1.
 
    xlb         Array of length nvars containing lower bounds of
                the variables.  If variable j has no lower bound,
                set xlb[j] = -1.0e30.
 
    xub         Array of length nvars containing upper bounds of
                the variables.  If variable j has no upper bound,
                set xub[j] = 1.0e30.
 
                Note:  If you wish to fix a variable at a certain
                value, say V, and have LSGRG2 leave it unchanged,
                set xlb[j]=V and xub[j]=V.
 
    nfuns       Number of functions, including the objective
                function. Functions are numbered 0 to nfun-1.
 
    nobj        Index of the objective function.  If nobj > 0,
                the objective function is maximized.  If nobj <
                0, the objective function is minimized and -nobj
                is used as the index of the objective function.
 
    glb         Array of length nfuns containing lower bounds of
                the constraint functions. If function i has no
                lower bound, set glb[i] = -1.0e30.
 
    gub         Array of length nfuns containing upper bounds of
                the constraint functions. If function i has no
                upper bound, set glb[i] = 1.0e30.
 
                Notes:  It does not matter what value is used for
                the bounds of the objective function in glb and
                gub.  Set these values to 0.0.
 
                Equality constraints are indicated by setting
                glb[j]=B and gub[i]=B where B is the value of the
                constraint.
 
                If a constraint g[i] is to be ignored by LSGRG2,
                set glb[i] = -1.0e30 and gub[i] = 1.0e30
 
    title       Character string which identifies the problem in
                the LSGRG2 report.
 
    report      Name of the file for the LSGRG2 report.  If this
                file does not exist, it is created.  If it does
                exists, it will be overwritten.
 
    xx          Array of length nvars containing initial values
                of the variables in xx[0] to xx[nvars-1].
 
                On exit, xx is set to the final values of the
                optimal variable settings determined by LSGRG2.
 
    inform      Termination message output from LSGRG2.
 
                inform      Message
 
                   0        Kuhn-Tucker conditions satisfied.
                            This is the best possible indicator
                            that an optimal point has been found.
 
                   1        Fractional change in objective less
                            than EPSTOP for NSTOP consecutive
                            iterations.  This is not as good as
                            inform=0, but still indicates the
                            likelihood that an optimal point has
                            been found.
 
                   2        All remedies have failed to find a
                            better point.  User should check
                            functions and bounds for consistency
                            and, perhaps, try other starting
                            values.
 
                   3        Number of completed one dimensional
                            searches exceeded LIMSER.  User
                            should check functions and bounds for
                            consistency and, perhaps, try other
                            starting values.  It might help to
                            increase limser.  Use
                            lgoptD ("limser", larger_value)  to
                            do this.
 
                   4        Objective function is unbounded.
                            LSGRG2 has observed dramatic change
                            in the objective function over
                            several steps.  This is a good
                            indication that the objective
                            function is unbounded.  If this is
                            not the case, the user should check
                            functions and bounds for consistency.
 
                   5        Feasible point not found.  LSGRG2 was
                            not able to find a feasible point.
                            If the problem is believed to be
                            feasible, the user should check
                            functions and bounds for consistency
                            and, perhaps, try other starting
                            values.
 
                   6        Degeneracy has been encountered.  The
                            point returned may be close to
                            optimal.  The user should check
                            functions and bounds for consistency
                            and, perhaps, try other starting
                            values.
 
                  -1        Fatal Error.  Some condition, such as
                            nvars < 0, was encountered.  LSGRG2
                            documented the condition in the
                            report and terminated.  In this case,
                            the user needs to correct the input
                            and rerun LSGRG2.
 
                  -2        Fatal Error while opening report file.
                            LSGRG2 terminated because the report file could
                            not be opened.  Check the file name given in
                            LSGRG2 argument report.
 
                  -3        Fatal Error.
                            Same as inform = -1.  In this case, LSGRG2 did
                            not open the report file.
                            Check input arguments and/or rerun with the
                            report option turned on so GRG2 can document the
                            condition in the report.

3.1.4 Termination Messages (inform)

LSGRG2 returns the outcome of the optimization process through argument inform. All possible values of inform are described in the previous section. When the results of LSGRG2 are not as expected, the user must determine what to do. This section provides some information on this topic and more information is given in Section 4. "Things to Try If GRG2 Isn't Working".

Messages inform=0, 1, and 2 are by far the most common. Message inform=0 implies the highest level of confidence that at least a local optimum has been found, message inform=1 less confidence, message inform=2 even less. In message inform=0, the Kuhn-Tucker conditions are first-order necessary conditions that hold if the current point is at least a local optimum and all functions have continuous first partial derivatives. For further explanation, see references 1, 2, and 7. In message inform=2, the following sequence of events has occurred: (1) No improved point was located along the last search direction, (2) a change of basis was attempted (if one had not already been done), (3) if the search direction was not the negative reduced gradient, this direction is tried, (4) if any variables with values at a bound have reduced gradient components indicating that releasing them from that bound could improve the objective, one such variable is allowed to leave its bound. In other words, LSGRG2 tried all known remedies, and none of these remedies improved the objective function, so the program terminated.

Regardless of which of termination messages inform=0, 1, and 2 is returned, the current point may be (nearly) optimal. Message inform=0 may fail to appear because the variables or constraints of the problem are poorly scaled; see Section 4.4 "Scaling".

Message inform=5 is returned when Phase I (see Section 1. "Introduction" for an explanation) terminates and the final point is not feasible. In this case, there may be no point satisfying all problem constraints.

For things to try if you are unsatisfied with the solution found by LSGRG2, see Section 4. "Things to Try If GRG2 Isn't Working".

3.2 Changing Algorithmic Parameters of LSGRG2

The behavior of LSGRG2 can be modified by specifying values of one or more of the parameters that control the behavior of the algorithm. This is done by calling lgoptD.

All parameters have default settings that are appropriate for most problems. Parameter names are case insensitive. "EPNEWT", "Epnewt", and "epnewt" are equivalent. lgoptD must be called before calling lsgrg2.

3.2.1 Prototype for lgoptD

The prototype for lgoptD is:

    TYPE_stdcall(void) lgoptD(char *name, double Dvalue);

Parameters are set with a statement of the form,

    lgoptD("parameter", value);

before the call to lsgrg2.

3.2.2 Descriptions of LSGRG2 Algorithmic Parameters

    Parameter   Description and default value
 
    epnewt      A constraint is assumed to be binding if it is
                within epnewt of one of its bounds.
                Default: lgoptD ("epnewt", 1.0e-4);
 
    epinit      If it is desired to run the problem with epnewt
                initially set fairly large and then tighten at
                the end of the optimization;  this is
                accomplished by assigning epinit the initial
                tolerance and epnewt the final one.
                Default: lgoptD ("epinit", 1.0e-4);
 
    epstop      This specifies the LSGRG2 convergence criteria.
                If the fractional change in the objective
                function is less than epstop for nstop
                consecutive iterations, the program will accept
                the current point as optimal. LSGRG2 will accept
                the current point as optimal if the Kuhn-Tucker
                optimality conditions are satisfied to within
                epstop.
                Default: lgoptD ("epstop", 1.0e-4);
 
    epspiv      If, in constructing the basis inverse, the
                absolute value of a prospective pivot element is
                less than epspiv, the pivot will be rejected and
                another pivot element will be sought.
                Default: lgoptD ("epspiv", 1.0e-6);
 
    ph1eps      If ph1eps is nonzero, the phase 1 objective is
                augmented by a multiple of the true objective.
                The multiple is selected so that, at the initial
                point, the ratio of the true objective and the
                sum of the infeasibilities isph1eps.
                Default: lgoptD ("ph1eps", 0.0);
 
    pstep       This is the step size used in parshf and parshc
                for estimating partial derivatives of functions
                with respect to the variables.
                Default: lgoptD ("pstep", sqrt(machine eps));
 
    nstop       If the fractional change in the objective
                function is less than epstop for nstop
                consecutive iterations, LSGRG2 will accept the
                current point as optimal.
                Default: lgoptD ("nstop", 3);
 
    itlim       If the Newton procedure takes itlim iterations
                without converging, the iterations are stopped
                and corrective action taken.
                Default: lgoptD ("itlim", 10);
 
    limser      If the number of completed one dimensional
                searches exceeds limser,  LSGRG2 terminates and
                returns inform = 3.
                Default: lgoptD ("limser", 10000);
 
    limser_p1   If the number of completed one dimensional
                searches exceeds limser in phase 1,  LSGRG2 
                terminates and returns inform = 3.
                Default: limser_p1 = limser;
 
    ipr         Print level for LSGRG2 report.
                   1  Print one summary line for each
                      one dimensional search.
                   2  Prints more detailed information of
                      the optimization process.
                Values of ipr >= 2 and <= 6 are permitted, but
                require knowledge of the internal workings of
                LSGRG2 and are not recommended for general use.
 
                Default: lgoptD ("ipr", 1);
 
    iquad       Method for initial estimates of basic variables
                for each one dimensional search
                   0    Tangent vectors and linear extrapolation
                   1    Quadratic extrapolation
 
                Default: lgoptD ("iquad", 0);
 
    kderiv      Method for computing estimates of partial
                derivatives of functions with respect to
                variables.
                   0    Forward difference approximations
                   1    Central difference approximations
 
                Default: lgoptD ("kderiv", 0);
 
    ckgrad      Check user supplied derivatives (see lgoptJ and
                lgoptP descriptions) with numerical
                approximations.  Derivatives are checked at the 
                starting point and any differences are documented
                in the report file.
                   0    Don't check user supplied derivatives
                   1    check user supplied derivatives
                   2    check user supplied derivatives and
                terminate if differences are detected
                Default: lgoptD ("ckgrad", 0);
 
    ckzero      Check zero derivatives.  The structure of the Jacobian
                is determined by evaluating derivatives at the starting
                point.  Only nonzero derivative values are stored.  If
                a derivative is zero at the starting point, but nonzero 
                elsewhere, the structure of the Jacobian will not be correct.
                LSGRG2 uses extra function evaluations to make sure that
                all nonzero derivatives are detected.  This process can take
                as many as nvars evaluations and may not be necessary.  To
                turn off the check on zero derivatives, use
                lgoptD ("ckzero", 0);
                Default: lgoptD ("ckgrad", 1);
 
    modcg       modcg and maxhes (see maxhes below) control the
                use of a conjugate gradient (CG) method. If the
                number of superbasic variables exceeds maxhes,
                the CG method indicated by modcg is used. When a
                CG method is selected, the Hessian approximation
                is not needed so maxhes is set to 0.  This
                reduces memory requirements by
                ((nvars)*(nvars+1))/2
                   1    Fletcher-Reeves.
                   2    Polak-Ribiere.
                   3    Perry.
                   4    1 step DFP.
                   5    1 step BFS.
                   6    Limited memory BFGS method.
 
                Default: lgoptD ("modcg", 6);
 
    maxcg       Length of memory (nvars*maxcg) for the limited
                memory BFGS method.  maxcg must be > 1.
                Increasing maxcg may improve performance for some
                problems.
                Default: lgoptD ("maxcg", 3);
 
    hscale      Automatic Hessian scaling for modcg == 6.
                   1, initial Hessian is not scaled.
                   0, initial Hessian is scaled.
                Default: lgoptD ("hscale", 1);
 
    maxhes      Maximum number of rows for the Hessian
                approximation for the BFGS algorithm.  LSGRG2
                allocates memory to accommodate 100 superbasics.
                Smaller settings of maxhes can be used to reduce
                memory requirements and to force a CG method to
                be used (see modcg description above).  Larger
                settings allow the BFGS algorithm to be used when
                the number of superbasics exceeds 100.  In
                particular, maxhes = 0 forces a CG method to be
                used at each iteration. Other values of maxhes
                depend on the number of superbasic variables to
                determine when to use a CG method.
                Default: lgoptD ("maxhes", max (nvars, 100) );
 
    iscale      Scaling.
                   0    No scaling.
                   1    The problem is scaled so that the maximum
                value of any row or column of the initial
                gradient array is less than or equal to 1.0
                Default: lgoptD ("iscale", 0);
 
    isclag      Rescale the Jacobian after every isclag line
                searches are performed.  if isclag = 0, then no
                rescaling is performed if isclag > 0, then rescale 
                after isclag line searches.
                Default: lgoptD ("isclag", 0);
 
    aijtol      Zero tolerance for partial derivatives.  Partial
                derivatives are set to zero if they are less than
                aijtol in absolute value.
                Default: lgoptD ("aijtol", 1.0e-10);
 
    fullj       Treat all partial derivatives as if they are nonzero.
                When this option is used by calling
                lgoptD ("fullj", 1);  aijtol is ignored.
                Default: lgoptD ("fullj", 0);
 
    pivpct      when choosing pivots for each row, all columns
                with updated entries within pivpct of the maximum
                will be considered.
                Default: lgoptD ("pivpct", 0.1);
 
    mxtabu      Maximum length of the tabu list in chuzq.
                Default: lgoptD ("mxtabu", 25);
 
    funpr       Estimate of the accuracy in the constraints when
                evaluated in gcomp.
                Default: lgoptD ("funpr", 1.0e-8);
 
    cndtol      Largest condition number which is acceptable.
                Blocks whose condition exceeds cdntol are
                considered to be ill-conditioned.
                Default: lgoptD ("cndtol", 1.0e8);
 
    idglim      maximum number of consecutive degenerate steps
                before the LSGRG2 terminates.
                Default: lgoptD ("idglim", 25);
 
    epboun      A variable is considered to be at its bound if it
                is within epboun of the bound.
                Default: lgoptD ("epboun", 1.0e-6);
 
    epdeg       If degeneracy is encountered, the variable bounds
                are perturbed by this relative tolerance, and the
                problem re-solved.
                Default: lgoptD ("epdeg", 1.0e-4);
 
    multsb      If multsb == 1, superbasics are projected on to
                their bounds during line search.
                If multsb == 0, the line search will end if a
                superbasic variable hits a bound.
                Default: lgoptD ("multsb", 1);
 
    ibvblm      If a basic variable is in the basis at a bound
                for "ibvblm" consecutive iterations, then consbs
                does not try to replace it.
                Default: lgoptD ("ibvblm", 2);
 
    hrdbnd      If hrdbnd == 0, basic variables are allowed to
                violate their bounds during Newton's method. If
                hrdbnd == 1, basic variables will never violate
                their bounds during Newton's method.
                Default: lgoptD ("hrdbnd", 0);
 
    fixpiv      If fixpiv == 1, dynamic adjustment of the pivot
                tolerance "pivpct" will be made. If fixpiv != 1,
                the pivot tolerance "pivpct" is fixed at the
                specified level.
                Default: lgoptD ("fixpiv", 0);
 
    lenz        Length of the work array.
                Default: lgoptD ("lenz", 100000);
 
    minimize    Minimize the objective function.
                Default: lgoptD ("minimize", 1);
 
    maximize    Maximize the objective function.
                Default: lgoptD ("minimize", 1);
 
    report      Normally, a report file is produced.  The report file
                can be turned off by, lgoptD ("report", 0);
                Default: lgoptD ("report", 1);
 
    flush       This option causes the buffer for the report file to be
                flushed after each iteration.  Useful for debugging.
                To envoke this option, lgoptD ("flush", 1);
                Default: lgoptD ("flush", 0);
 
    pr_start    Print starting values of functions and variables
                in the report file.  These values are included in
                the report by default when nvars <= 30 and 
                nfuns <= 31.  Use lgoptD ("pr_start", 1); to
                include starting values of functions in the report.  
                Use lgoptD ("pr_start", 2); to include starting values 
                of variables in the report. Use lgoptD ("pr_start", 3); 
                to include starting values of functions and variables
                in the report. Use lgoptD ("pr_start", 0); when you do 
                not want starting values of functions or variables in the 
                report.
 
    pr_final    Print final values of functions and variables
                in the report file.  These values are included in
                the report by default when nvars <= 30 and 
                nfuns <= 31.  Use lgoptD ("pr_final", 1); to
                include final values of functions in the report.  
                Use lgoptD ("pr_final", 2); to include final values 
                of variables in the report. Use lgoptD ("pr_final", 3); 
                to include final values of functions and variables
                in the report.   Use lgoptD ("pr_final", 0); when you do 
                not want final values of functions or variables in the 
                report.
 
    pr_detail   Print detailed algorithm statistics in the report.
                Default: lgoptD ("pr_detail", 0);
 
    default     This returns all lgoptD options to their default
                setting. Use lgoptD ("default", 1);

3.3 Providing Derivatives to LSGRG2

The LSGRG2 package contains two functions that compute derivatives numerically. One of these functions, parshf, computes derivatives by forward differences requiring nvars calls to GCOMP, and the other one, parshc, computes derivatives by central differences using 2*nvars calls to GCOMP. The forward difference function is the default and its performance is sufficient for most problems. The user may select the central difference function when higher accuracy is needed in the derivatives by using lgoptD("kderiv", 1). If LSGRG2 fails to reach an optimal solution with these numerical methods and if it is suspected that inaccurate derivatives are at fault, the user may prepare a function to compute the derivatives. There are other reasons for providing derivatives to LSGRG2. When the number of variables is in the thousands, it takes thousands of calls to GCOMP to compute the necessary derivatives. There may be structural aspects of the problem, such as linearity, that can be used to advantage to substantially reduce the time needed to compute derivatives. LSGRG2 provides two ways to do this.

The first way of providing derivatives is through lgoptJ. In this case, the user provides a C function that computes the derivatives with respect to one variable at a time and returns an array of length nfuns containing the derivatives to LSGRG2. This way of providing derivatives can reduce computation time substantially. This format for providing derivatives is described in Section 3.3.1 "Providing Derivatives Through lgoptJ".

The second way of providing derivatives is through lgoptP. This way results in the greatest efficiency gain, but is also the most difficult. It can have a good payoff for very large problems. The user provides a C function that computes the derivatives and returns them in an array stored in packed column format. The packed column format can result in a 25% timesaving over the single column format. This format for providing derivatives is described in Section 3.3.2 "Providing Derivatives Through lgoptP".

The most common problem with user provided derivatives is that the C function that performs the computations contains errors and returns incorrect derivative values. The behavior of LSGRG2 can be very perplexing under these circumstances and care must be taken to prevent this from happening. See Section 4.2 "Incorrect User-Supplied Derivatives".

3.3.1 Providing Derivatives Through lgoptJ

The prototype for lgoptJ is:

    TYPE_stdcall(void) lgoptJ(
        void (STDCALL *Pname)(
            double    x[],
            double    g[],
            long int  jcol,
            double    gradj[] ));

The user provides a C function to LSGRG2 by inserting a statement of the form,

    lgoptJ(parshj);

before the call to lsgrg2. The C function parshj must compute derivatives of the constraint functions and the objective function (nfuns functions) with respect to variable x[jcol] and return them in the array grad, elements grad[0], ..., grad[nfuns-1]. It must conform to the prototype given for lgoptJ so that it can be called with a statement of the form,

    parshj(x, g, jcol, grad);

Example for lgoptJ

The function to compute derivatives with respect to one variable at a time for the example problem (discussed in Section 2.1 "Example Problem"), named h11J, is given below. It can be passed to LSGRG2 by calling lgoptJ(h11J) prior to calling lsgrg2.

TYPE_stdcall(void) h11J (double x[5], double g[4],
                  long int jcol, double gradj[4]) {
/*  Himmelblau Problem 11; nvars=5; nfuns=4; nobj=3;  */
 
    if(jcol == 0) {
        gradj[0] =  0.0006262*x[3];
        gradj[1] =  0.0029955*x[1];
        gradj[2] =  0.0012547*x[2];
        gradj[3] =  0.8356891*x[4] + 37.293239;
        }
    else if(jcol == 1) {
        gradj[0] =  0.0056858*x[4];
        gradj[1] =  0.0071317*x[4] + 0.0029955*x[0];
        gradj[2] =  0.0;
        gradj[3] =  0.0;
        }
    else if(jcol == 2) {
        gradj[0] = -0.0022053*x[4];
        gradj[1] =  0.0021813*x[2]*2.0;
        gradj[2] =  0.0047026*x[4] + 0.0019085*x[3] + 0.0012547*x[0];
        gradj[3] =  5.3578547*x[2]*2.0;
        }
    else if(jcol == 3) {
        gradj[0] =  0.0006262*x[0];
        gradj[1] =  0.0;
        gradj[2] =  0.0019085*x[2];
        gradj[3] =  0.0;
        }
    else if(jcol == 4) {
        gradj[0] =  0.0056858*x[1] - 0.0022053*x[2];
        gradj[1] =  0.0071317*x[1];
        gradj[2] =  0.0047026*x[2];
        gradj[3] =  0.8356891*x[0];
    }
 
    return;
 
}   /* end h11J */

3.3.2 Providing Derivatives Through lgoptP

The prototype for lgoptP is:

TYPE_stdcall(void) lgoptP(
    void (STDCALL *Pname)(
        double    x[],
        double    g[],
        long int  igrad[],
        long int  jgrad[],
        double    gradj[],
        long int  nzg ));

The user provides a C function to LSGRG2 by inserting a statement of the form,

    lgoptP(parshp);

before the call to lsgrg2. The C function parshp must conform to the prototype given for lgoptP so that it can be called with a statement of the form,

    parshp(x, g, igrad, jgrad, grad, ngz);

It must compute the derivatives at the point x and return them in array grad stored in packed column format.

Arguments igrad and jgrad are the row and column index arrays for grad. Array jgrad has nvars+1 elements and array igrad has nz elements. jgrad[j]is the index of the beginning of column j for j=0, ..., nvars-1. The elements grad[k], k = jgrad[j],..., jgrad[j+1]-1, are all the nonzero elements of column j. igrad[k] is the row index of grad[k]. jgrad[nvars] is the index of the first unused element of grad and igrad. parshp must return the index arrays igrad and jgrad the first time it is called; indicated by ngz>0. On the first call, ngz also provides the number of elements available for grad and igrad. These arrays must not use more that nzg elements. If all derivatives cannot be returned in the array grad, elements grad[0] to grad[nzg], parshp must set jgrad[nvars]=0 and return. This is an error condition that occurs when working storage is not sufficient. It can usually be corrected by using lgoptD("lenz", larger_number_of_units) to increase the amount of working storage.

When parshp is called, after the first call, nzg is 0. parshp must compute values for derivatives at the point x. Derivatives that are constant need not be computed and the index arrays igrad and jgrad must not be changed.

Example for lgoptP

This problem has 4 variables {X0, X1,..., X3}, bounded above and below, a quadratic objective function {G3}, and three linear constraints {G0, G1, G2}, with both lower and upper bounds.

Minimize: X0X0 + X1X1 + X2X2 + X3X3

Subject to:

            0 < 2X0 + X1 < 10
            0 < X0 + 2X1 + X2 < 10
            0 < X1 + 2X2 + X3 < 10

and,

            2 < Xi < 100, i = 0,1,2,3

This problem is solved starting from X = {5, 5, 5, 5} which is infeasible because the third constraint, G2, is not satisfied. The array of partial derivatives is,

2

1

0

0

1

2

1

0

0

1

2

1

2X0

2X1

2X2

2X3

The representation of this array in packed column form is illustrated in the following table. The rows labeled index and column are not part of the representation, but are provided for clarity.

index

0