User's Guide
|
Contents
2.2 C Program for Example Problem
2.3 Fortran Program for Example Problem
2.4 MATLAB Program for Example Problem
3.1.1 C Prototypes for goptG and grg2
3.1.2 Fortran Interface for goptgf and grg2f
3.1.3 MATLAB Interface for grg2m
3.1.4 Description of GRG2 Arguments
3.1.5 Termination Messages (inform)
3.2 Changing GRG2 Algorithmic Parameters
3.2.2 Fortran Interface for goptdf
3.2.3 MATLAB Interface for goptdm
3.2.4 Descriptions of GRG2 Algorithmic Parameters
3.3 Providing Derivatives to GRG2
3.3.2 Fortran Interface for goptpf
3.3.3 MATLAB Interface for goptpm
3.4 Using goptL to Set Variable and Function Labels
3.5 Understanding the GRG2 Report
4. Things to Try If GRG2 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.5 Bounds Which Must Not be Violated
4.7 Solution Not Accurate Enough
5. GRG2 Options and Tolerances
6. Unconstrained Minimization and Optimal Control
GRG2 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=1,...,m, i ¹ p |
|
xlbj £ xj £ xubj |
for j=1,...,n. |
X is a vector on n variables, x1 ,...,xn, and the functions g1 ,...,gm 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 solves problems of this form by the Generalized Reduced Gradient Methods (see references 1-4 for a description and further references).
GRG2 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 is 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.
This User’s Guide describes the PC Windows interface for GRG2. The PC Windows DLL version is specially designed for the MS Windows operating system. The DLL can be used by C programs, by Fortran programs, by MATLAB, and by any Windows application capable of referencing a DLL. This manual describes the C interface and gives examples using C, Fortran, and MATLAB. The Visual Basic interface is described in a README file on the GRG2 product disk.
There are a few essential steps that must be taken in order to present an optimization problem to GRG2. First, a function/subroutine must be written to compute values of the constraint functions and the objective function for given values of the variables. This function/subroutine is referred to as GCOMP in this manual. Next, the data that describe the optimization problem must be created, organized, and passed to GRG2. 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. GRG2 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, GRG2 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 GRG2. If GRG2 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.
It is possible to replace the automatic derivative computation with user-supplied computations. The user-supplied computations can often take advantage of the problem 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. This feature is described in Section 3.3 "Providing Derivatives to GRG2"
There are many parameters that control the behavior of the GRG2 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 GRG2 more efficient or in order to make it possible to solve a difficult problem. See the description in Section 3.2 "Changing Algorithmic Parameters of GRG2".
The following table summarizes the user callable functions in the GRG2 program and gives the section where they are documented.
|
Function |
Description |
|
goptG, goptgf, & goptgm |
Required: provides the name of the function that evaluates constraint functions and the objective function. See Section 3.1 for discussion and examples. |
|
grg2, grg2f, & grg2m |
Required: main interface, provides required problem specifications. See Section 3.1 for discussion and examples. |
|
goptD, goptdf, & goptdm |
Optional: optional parameter values. See Section 3.2 for discussion and examples. |
|
goptP, goptpf, & goptpm |
Optional: derivative evaluation function. See Section 3.3 for discussion and examples. |
|
goptL |
Optional: labels for variables and functions. See Section 3.4 for discussion and examples. |
This section presents an example problem with C, Fortran, and MATLAB programs to illustrate the essential steps necessary to prepare a program to use GRG2.
This problem has 5 variables {X1, X2, ..., X5}, bounded above and below, a quadratic objective function {G4}, and three quadratic constraints {G1, G2, G3}, with both upper and lower bounds. It is problem number 11 in the Himmelblau text [7].
Minimize: G4 = 5.3578547X3X3 + 0.8356891X1X5 + 37.293239X1 - 40792.141
Subject to:
0 < G1 = 85.334407 + 0.0056858X2X5 + 0.0006262X1X4 - 0.0022053X3X5 < 92 90 < G2 = 80.51249 + 0.0071317X2X5 + 0.0029955X1X2 + 0.0021813X3X3 < 110 20 < G3 = 9.300961 + 0.0047026X3X5 + 0.0012547X1X3 + 0.0019085X3X4 < 25
and,
78 < X1 < 102; 33 < X2 < 45; 27 < X3 < 45; 27 < X4 < 45; 27 < X5 < 45
This problem is solved starting from X = {78, 33, 27, 27, 27} which is infeasible because the third constraint, G3, is not satisfied at this point.
This section presents a C program that uses GRG2 to solve the example problem. The C function h11G evaluates constraint functions and the objective function and h11P evaluates derivatives. The program is contained in file Ex1.c on the GRG2 product disk.
#include "grg2.h"
TYPE_export(void) h11G(double g[5], double x[6]);
TYPE_export(void) h11P (double g[5], double x[6],
int nfuns, int nvars, float *grad);
void main () {
int nvars, nfuns, nobj, inform;
double xx[25], xlb[25], xub[25], glb[25], gub[25];
char *title, *report;
nvars = 5; /* number of variables */
nfuns = 4; /* number of functions */
nobj = 4; /* use -nobj in the call to grg2 to minimize g[4] */
/* bounds on variables */
xlb[1] = 78.0; xub[1] = 102.0;
xlb[2] = 33.0; xub[2] = 45.0;
xlb[3] = 27.0; xub[3] = 45.0;
xlb[4] = 27.0; xub[4] = 45.0;
xlb[5] = 27.0; xub[5] = 45.0;
/* bounds on constraint functions */
glb[1] = 0.0; gub[1] = 92.0;
glb[2] = 90.0; gub[2] = 110.0;
glb[3] = 20.0; gub[3] = 25.0;
/* title for report */
title = "GRG2 C Example: Himmelblau Problem 11";
/* file name for GRG2 report */
report = "H11.txt";
/* starting values for variables */
xx[1] = 78; xx[2] = 33; xx[3] = 27; xx[4] = 27; xx[5] = 27;
/* GRG2 will use h11P to evaluate derivatives */
goptP (h11P);
goptG(h11G); /* h11G will be used to evaluate constraints and
the objective function */
grg2 (nvars, xlb, xub, nfuns, -nobj, glb, gub, title, report,
xx, &inform);
} /* end h11main */
TYPE_export(void) h11G (double g[5], double x[6]) {
/* Himmelblau Problem 11; nvars = 5; nfuns = 4; nobj = 4; */
g[1] = 85.334407 + 0.0056858*x[2]*x[5] + 0.0006262*x[1]*x[4]
- 0.0022053*x[3]*x[5];
g[2] = 80.51249 + 0.0071317*x[2]*x[5] + 0.0029955*x[1]*x[2]
+ 0.0021813*x[3]*x[3];
g[3] = 9.300961 + 0.0047026*x[3]*x[5] + 0.0012547*x[1]*x[3]
+ 0.0019085*x[3]*x[4];
g[4] = 5.3578547*x[3]*x[3] + 0.8356891*x[1]*x[5] + 37.293239*x[1]
- 40792.141;
} /* end h11G */
TYPE_export(void) h11P (double g[5], double x[6],
int nfuns, int nvars, float *grad) {
/* Himmelblau Problem 11; nvars = 5; nfuns = 4; nobj = 4; */
/* h11P computes partial derivatives
grad[i*(nvars+1)+j] = Dg[i]/Dx[j] */
grad[1*(nvars+1)+1] = (float) (0.0006262*x[4]);
grad[1*(nvars+1)+2] = (float) (0.0056858*x[5]);
grad[1*(nvars+1)+3] = (float) (-0.0022053*x[5]);
grad[1*(nvars+1)+4] = (float) (0.0006262*x[1]);
grad[1*(nvars+1)+5] = (float) (0.0056858*x[2] - 0.0022053*x[3]);
grad[2*(nvars+1)+1] = (float) (0.0029955*x[2]);
grad[2*(nvars+1)+2] = (float) (0.0071317*x[5] + 0.0029955*x[1]);
grad[2*(nvars+1)+3] = (float) (0.0021813*x[3]*2.0);
grad[2*(nvars+1)+4] = (float) (0.0);
grad[2*(nvars+1)+5] = (float) (0.0071317*x[2]);
grad[3*(nvars+1)+1] = (float) (0.0012547*x[3]);
grad[3*(nvars+1)+2] = (float) (0.0);
grad[3*(nvars+1)+3] = (float) (0.0047026*x[5] + 0.0019085*x[4]
+ 0.0012547*x[1]);
grad[3*(nvars+1)+4] = (float) (0.0019085*x[3]);
grad[3*(nvars+1)+5] = (float) (0.0047026*x[3]);
grad[4*(nvars+1)+1] = (float) (0.8356891*x[5] + 37.293239);
grad[4*(nvars+1)+2] = (float) (0.0);
grad[4*(nvars+1)+3] = (float) (5.3578547*x[3]*2.0);
grad[4*(nvars+1)+4] = (float) (0.0);
grad[4*(nvars+1)+5] = (float) (0.8356891*x[1]);
} /* end h11P */
This section presents a Fortran program that uses GRG2 to solve the example problem. The Fortran subroutine h11gf evaluates constraint functions and the objective function and h11pf evaluates derivatives. The program is contained in file Ex1.for on the GRG2 product disk.
program h11fmain
implicit none
! GRG2F arguments
integer nvars, nfuns, nobj, inform
double precision xx(25), xlb(25), xub(25), glb(25), gub(25)
character title*50, report*50
! externals
external h11gf, h11pf
nvars = 5 ! number of variables
nfuns = 4 ! number of functions
nobj = 4 ! use -nobj in the call to grg2f to minimize g(4)
! bounds on variables
xlb(1) = 78.0; xub(1) = 102.0;
xlb(2) = 33.0; xub(2) = 45.0;
xlb(3) = 27.0; xub(3) = 45.0;
xlb(4) = 27.0; xub(4) = 45.0;
xlb(5) = 27.0; xub(5) = 45.0;
! bounds on constraint functions
glb(1) = 0.0; gub(1) = 92.0;
glb(2) = 90.0; gub(2) = 110.0;
glb(3) = 20.0; gub(3) = 25.0;
glb(4) = 0.0; gub(4) = 0.0;
! variable starting values
xx(1) = 78; xx(2) = 33; xx(3) = 27; xx(4) = 27; xx(5) = 27
! report filename and title
report = 'H11f.txt'
title = 'GRG2 Fortran Example: Himmelblau Problem 11'
call goptpf (h11pf) ! h11pf will be used to evaluate derivatives
call goptgf (h11gf) ! h11gf will be used to evaluate constraints
! and the objective function
call grg2f (nvars, xlb, xub, nfuns, -nobj, glb, gub,
& title, report, xx, inform)
end ! h11fmain
subroutine h11gf (g, x)
! Himmelblau Problem 11; nvars = 5; nfuns = 4; nobj = 4;
implicit none
! arguments
double precision g(*), x(*)
g(1) = 85.334407 + 0.0056858*x(2)*x(5) + 0.0006262*x(1)*x(4)
& -0.0022053*x(3)*x(5)
g(2) = 80.51249 + 0.0071317*x(2)*x(5) + 0.0029955*x(1)*x(2)
& + 0.0021813*x(3)*x(3)
g(3) = 9.300961 + 0.0047026*x(3)*x(5) + 0.0012547*x(1)*x(3)
& + 0.0019085*x(3)*x(4)
g(4) = 5.3578547*x(3)*x(3) + 0.8356891*x(1)*x(5) + 37.293239*x(1)
& - 40792.141
return
end ! h11gf
subroutine h11pf (g, x, nfuns, nvars, grad)
! Himmelblau Problem 11; nvars = 5; nfuns = 4; nobj = 4;
implicit none
! arguments
integer nfuns, nvars
double precision g(nfuns), x(nvars)
real grad(nfuns, nvars)
! h11pf computes partial derivatives
! grad(i,j) = Dg(i)/Dx(j)
grad(1, 1) = 0.0006262*x(4)
grad(1, 2) = 0.0056858*x(5)
grad(1, 3) = -0.0022053*x(5)
grad(1, 4) = 0.0006262*x(1)
grad(1, 5) = 0.0056858*x(2) - 0.0022053*x(3)
grad(2, 1) = 0.0029955*x(2)
grad(2, 2) = 0.0071317*x(5) + 0.0029955*x(1)
grad(2, 3) = 0.0021813*x(3)*2.0
grad(2, 4) = 0.0
grad(2, 5) = 0.0071317*x(2)
grad(3, 1) = 0.0012547*x(3)
grad(3, 2) = 0.0
grad(3, 3) = 0.0047026*x(5) + 0.0019085*x(4)
& + 0.0012547*x(1)
grad(3, 4) = 0.0019085*x(3)
grad(3, 5) = 0.0047026*x(3)
grad(4, 1) = 0.8356891*x(5) + 37.293239
grad(4, 2) = 0.0
grad(4, 3) = 5.3578547*x(3)*2.0
grad(4, 4) = 0.0
grad(4, 5) = 0.8356891*x(1)
return
end ! h11pf
This section presents a MATLAB program that uses GRG2 to solve the example problem. The MATLAB function hmbl11 evaluates constraint functions and the objective function and hmbl11p evaluates derivatives. The program is contained in files Ex1p.m, Hmbl11.m, and Hmbl11p.m on the GRG2 product disk.
% Example problem for grg2m
% Problem Himmelblau 11
% 5 variables and 3 constraints
% Constraints and objective defined in hmbl11.m
% Partial derivatives evaluated in hmbl11p.m
xbnd = [78 102; 33 45; 27 45; 27 45; 27 45];
gbnd = [0 92; 90 110; 20 25; 0 0];
nobj = -4;
gcomp = 'hmbl11';
goptpm ('hmbl11p');
xstrt = [78.62 33.44 31.07 44.18 35.22]';
title = 'Matlab: Himmelblau 11';
report = 'hmbl11p.txt';
[xopt, inform] = grg2m (xbnd, gbnd, nobj, gcomp, xstrt, title, report);
' Minimized objective for hmbl11 problem is'
g = hmbl11(xopt);
g(abs(nobj))
function g = hmbl11 (x)
% Himmelblau Problem 11, 5 variables and 4 functions */
g(1) = 85.334407 + 0.0056858*x(2)*x(5) + 0.0006262*x(1)*x(4) ...
- 0.0022053*x(3)*x(5);
g(2) = 80.51249 + 0.0071317*x(2)*x(5) + 0.0029955*x(1)*x(2) ...
+ 0.0021813*x(3)*x(3);
g(3) = 9.300961 + 0.0047026*x(3)*x(5) + 0.0012547*x(1)*x(3) ...
+ 0.0019085*x(3)*x(4);
g(4) = 5.3578547*x(3)*x(3) + 0.8356891*x(1)*x(5) + 37.293239*x(1) ...
- 40792.141;
% end hmbl11
function grad = hmbl11p (g, x, nfuns, nvars)
% Himmelblau Problem 11, 5 variables and 4 functions */
% compute partial derivatives
% Dg(i)/Dx(j) to grad(i, j)
grad(1, 1) = 0.0006262*x(4);
grad(1, 2) = 0.0056858*x(5);
grad(1, 3) = -0.0022053*x(5);
grad(1, 4) = 0.0006262*x(1);
grad(1, 5) = 0.0056858*x(2) - 0.0022053*x(3);
grad(2, 1) = 0.0029955*x(2);
grad(2, 2) = 0.0071317*x(5) + 0.0029955*x(1);
grad(2, 3) = 0.0021813*x(3)*2.0;
grad(2, 4) = 0.0;
grad(2, 5) = 0.0071317*x(2);
grad(3, 1) = 0.0012547*x(3);
grad(3, 2) = 0.0;
grad(3, 3) = 0.0047026*x(5) + 0.0019085*x(4) + 0.0012547*x(1);
grad(3, 4) = 0.0019085*x(3);
grad(3, 5) = 0.0047026*x(3);
grad(4, 1) = 0.8356891*x(5) + 37.293239;
grad(4, 2) = 0.0;
grad(4, 3) = 5.3578547*x(3)*2.0;
grad(4, 4) = 0.0;
grad(4, 5) = 0.8356891*x(1);
% end hmbl11p
This section describes how to use the GRG2 interface for C, Fortran, and MATLAB. There are two required steps to use GRG2. Step 1 is to provide GRG2 with the user-supplied function/subroutine that evaluates the constraint functions and the objective function. Step 2 is to provide all other required information, such as variable bounds and constraint bounds.
Section 3.1 describes the GRG2 interface for these essential steps. The following table gives a quick reference of the language specific interfaces for carrying out Steps 1 and Step 2. Refer to Section 2 for complete example programs.
Language |
Interface |
goptG(gcomp);
grg2 (nvars, xlb, xub, nfuns, nobj, glb, gub,
title, report, xx, &inform); |
|
call goptgf(gcomp)
call grg2f (nvars, xlb, xub, nfuns, nobj, glb, gub,
title, report, xx, inform) |
|
[xopt, inform] = grg2m (xbnd, gbnd, nobj, xstrt, 'gcomp'
title, report); |
Section 3.2 describes how to use other GRG2 support functions to handle difficult problems or to modify the behavior of GRG2 to take advantage of special, problem specific, conditions.
Section 3.3 describes how to provide derivatives. Most problems can be solved without providing derivatives (GRG2 computes then automatically) but, for some problems, much computation time can be saved by providing an efficient function/subroutine to compute partial derivatives.
The discussion in the remainder of this Section uses the following problem statement:
Minimize gnobj
Subject to:
|
glbi £ gi £ gubi, |
i = 1,..., nfuns |
|
|
i ¹ nobj |
|
xlbj £ xj £ xubj, |
j = 1,..., nvars |
where xj are variables and gi are nonlinear functions that specify problem constraints.
This section describes the GRG2 interface for providing the user-supplied function/subroutine that evaluates the constraint functions and the objective function and the interface for providing all other required information, such as variable bounds and constraint bounds
The prototypes for all C functions in the GRG2 program are included in the header file grg2.h and can be accessed by including this file in your C program.
TYPE_export(void) goptG (
void STDCALL *GCOMP) (
double g[], double x[]));
TYPE_export(void) grg2 (
int nvars,
double *xlb,
double *xub,
int nfuns,
int nobj,
double *glb,
double *gub,
char *title,
char *report,
double *xx,
int *inform);
Notes: The macros TYPE_export(void) and STDCALL are defined through the header file grg2.h for the target system.
The interfaces for Fortran subroutines goptgf and grg2f in the GRG2 program are:
subroutine goptgf (gcomp)
external gcomp
subroutine grg2f (nvars, xlb, xub, nfuns, nobj, glb, gub,
& title, report, xx, inform)
integer nvars, nfuns, nobj, inform
double precision xlb(nvars), xub(nvars), xx(nvars),
& glb(nfuns), gub(nfuns)
character title*(*), report*(*)
The arguments for grg2m are documented in the MATLAB help file grg2m.m in the GRG2 product disk.
[xopt, inform] = grg2m (xbnd, gbnd, nobj, 'gcomp', xstrt, 'title', 'report');
Input:
xbnd - nvars by 2 matrix.
xbnd(i,1) is the lower bound for variable i
xbnd(i,2) is the upper bound for variable i
gbnd - nfuns by 2 matrix.
xbnd(i,1) is the lower bound for variable i
xbnd(i,2) is the upper bound for variable I
nobj - Scalar
gcomp - MATLAB function that evaluates the problem constraint
functions and the objective function. g = gcomp(x).
xstrt - Vector of length nvars
xstrt(1) to xstrt(nvars) initial values for variables.
title - Title written to the GRG2 report file
report - The GRG2 report is written to the file named 'report'
Output:
inform - Scalar
xopt - Vector of length nvars.
xopt(1) to xopt(nvars) are the final values of optimal variable
settings determined by GRG2
The arguments for grg2m are described further in the next section. The MATLAB interface differs in the follow way: argument xbnd of grg2m combines xlb and xub; argument gbnd combines glb and gub. The number of variables, nvars, is determined by the length of xbnd and the number of functions, nfuns, is determined by the length of gbnd. See Section 2.4 MATLAB Program for Example Problem
Argument Description
nvars Number of variables. Variables are numbered 1 to nvars.
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 and
have GRG2 leave it unchanged, set its entries in xlb
and xub to that value.
nfuns Number of functions, including the objective function.
Functions are numbered 1 to nfun.
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] = gub[i] = B where B is the value of the
constraint.
If a constraint g[i] is to be ignored by GRG2,
set glb[i] = -1.0e30 and gub[i] = 1.0e30
title character string which identifies the problem
in the GRG2 report.
report Name of the file for the GRG2 report. If this file does not
exist it will be created. If it exists, it will be
overwritten.
xx Array of length nvars containing initial values of
the variables in xx[1] to xx[nvars].
On exit, xx is set to the final values of the optimal
variable settings determined by GRG2.
inform Termination message output from GRG2.
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.
4 Objective function is unbounded.
GRG2 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.
GRG2 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. Degeneracy occurs rarely and
you can usually get around it by making a small
change in the starting point or making a small
change in one or more of the bounds on functions
or variables.
7 Noisy and nonsmooth function values.
Possible singularity or error in the function
evaluations.
8 Optimization process terminated by user request.
9 Maximum number of function evaluations exceeded.
-1 Fatal Error.
Some condition, such as nvars < 0, was encountered.
GRG2 documented the condition in the report and
terminated. In this case, the user needs to
correct the input and rerun GRG2.
-2 Fatal Error while opening report file.
GRG2 terminated because the report file could not
be opened. Check the file name given in GRG2
argument report.
-3 Fatal Error.
Same as inform = -1. In this case, GRG2 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.
GRG2 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 GRG2 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, GRG2 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 terminates and the final point is not feasible (see Section 1. "Introduction" for an explanation) . 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 GRG2, see Section 4. "Things to Try If GRG2 Isn't Working".
The behavior of GRG2 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 goptD in C programs, calling goptdf in Fortran programs, and calling goptdm in MATLAB programs. All GRG2 parameters have default settings that are 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 GRG2 more efficient or in order to make it possible to solve a difficult problem. The GRG2 parameter descriptions in Section 3.2.4 are described for use in C programs. The translation needed to use them with Fortran or MATLAB is straightforward.
Parameter names are case insensitive. "EPNEWT", "Epnewt", and "epnewt" are equivalent. goptD must be called before calling grg2.
The prototype for goptD is:
TYPE_export(void) goptD (char *name, double Dvalue);
Parameters are set by inserting a statement of the form,
goptD("parameter", value);
before the call to grg2. Note that argument Dvalue in gotpD is type double. Type int values can be used and they are coerced to type double as long as the grg2.h header file is used.
In a C program, to set the GRG2 option that checks user-supplied derivatives, use
goptD("ckgrad", 1);
before the call to grg2.
The Fortran interface for subroutine goptdf is:
subroutine goptdf (name, dvalue)
character name*(*)
double precision dvalue
Fortran requires that goptdf argument dvalue be provided as a double precision value (numbers of the form, 1d0, 1.0d-5, ...).
In a Fortran program, to set the GRG2 option that checks user-supplied derivatives, use
call goptdf('ckgrad', 1d0)
before the call to grg2.
The MATLAB interface for goptdm is:
goptdm ('name', value)
In a MATLAB program, to set the GRG2 option that checks user-supplied derivatives, use
goptdm('ckgrad', 1)
before the call to grg2.
Parameter Description and default value
epnewt A constraint is assumed to be
binding if it is within epnewt
of one of its bounds.
Default: goptD ("epnewt", 1.0e-6);
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: goptD ("epinit", 1.0e-6);
epstop This specifies the GRG2 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.
GRG2 will accept the current point as optimal
if the Kuhn-Tucker optimality conditions are
satisfied to within epstop.
Default: goptD ("epstop", 1.0e-4);
epskt Convergence, inform = 1, requires that the K-T
factor is less than epskt.
Default: goptD ("epskt", 0.01);
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
elememt will be sought.
Default: goptD ("epspiv", 1.0e-4);
.
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 is
ph1eps.
Default: goptD ("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: goptD ("pstep", 1.0e-8);
nstop If the fractional change in the objective
function is less than epstop for nstop
consecutive iterations, GRG2 will accept the
current point as optimal.
Default: goptD ("nstop", 3);
itlim If the Newton procedure takes itlim
iterations without converging, the iterations
are stopped and corrective action taken.
Default: goptD ("itlim", 10);
limser If the number of completed one dimensional
searches exceeds limser, GRG2 terminates
and returns inform = 3.
Default: goptD ("limser", 10000);
ipr Print level for GRG2 report.
1 Print one summary line for each one
dimensional search.
2 Provides more detailed information of the
optimization process and progress.
Values of ipr >= 2 and <= 6 are permitted, but
require knowledge of the internal workings of
GRG2 and are not recommended for general use.
Default: goptD ("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: goptD ("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: goptD ("kderiv", 0);
ckgrad Check user supplied derivatives (see goptP descriptions)
with numerical approximations. Any differences are written
to 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);
modcg modcg and maxr (see maxr below) control the use
of a conjugate gradient (CG) method.
If the number of superbasic variables
exceeds maxr, the CG method selected by modcg
is used.
The CG methods do not need the Hessian approximation.
When maxr is set to 0 the memory requirements of GRG2
are reduced by ((nvars+1)*(nvars+2))/2 and the CG method
selected by modcg is used. When the number of superbasics
is <= maxr, the BFGS variable metric method is used.
1 Fletcher-Reeves.
2 Polak-Ribiere.
3 Perry.
4 1 step DFP.
5 1 step BFS.
Default: goptD ("modcg", 1);
maxr Maximum number of rows for the Hessian approximation
for the BFGS algorithm. GRG2 allocates memory to
accommodate maxr = nvars. Smaller settings of maxr
can be used to reduce memory requirements and to force
a CG method to be used (see modcg description above).
In particular, maxr = 0 forces a CG method to be used at
every iteration. Other values of maxr depend on the number
of superbasic variables to determine when to use a CG
method.
Default: goptD ("maxr", nvars);
doscale 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: goptD ("doscale", 0);
minimize Minimize the objective function.
Whether the objective function is minimized
or maximized is usually determined by the sign
of argument nobj.
Use goptD ("minimize", 1); to force
GRG2 to minimize.
Default: goptD ("minimize", 0)
maximize Maximize the objective function.
Use goptD ("maximize", 1); to force
GRG2 to maximize.
Default: goptD ("minimize", 0);
limeval Limit on the number of function evaluations. limeval=0
is treated as no limit on the number of evaluations.
Default: goptD ("limeval", 0);
report Normally, a report file is produced. The report file
can be turned off by, goptD ("report", 0.0);
Default: goptD ("report", 1.0);
flush This option causes the buffer for the report file to be
flushed after each iteration. Useful for debugging.
To envoke this option, goptD ("flush", 1.0);
Default: goptD ("flush", 0.0);
default This returns all goptD options to their
default setting.
Use goptD ("default", -1.0);
The GRG2 program contains two functions that compute derivatives numerically. One of these functions, parshf, computes derivatives by forward differences requiring nvars calls to GCOMP for each derivative calculation, and the other one, parshc, computes derivatives by central differences using 2*nvars calls to GCOMP for each derivative calculation. 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 goptD("kderiv", 1). If GRG2 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 GRG2. When the number of variables is in the hundreds, it takes hundreds 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.
A user prepared function to compute derivatives can be passed to GRG2 through goptP, goptpf, or goptpm. The user provides a C function, Fortran subroutine, or MATLAB function that computes the derivatives and returns them in an array stored in packed column format (two dimensional array format for Fortran and MATLAB). This format for providing derivatives is described in the next three sections for C, Fortran, and MATLAB.
The most common problem with user provided derivatives is that the code that performs the computations contains errors and returns incorrect derivative values. The behavior of GRG2 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".
The prototype for goptP is:
TYPE_export(void) goptP(
void (STDCALL *Pname)(
double g[],
double x[],
int nfuns,
int nvars,
float *grad));
Note: The macros TYPE_export(void) and STDCALL are defined through the header file grg2.h for the target system.
The user provides a C function to GRG2 by inserting a statement of the form,
goptP(pname);
before the call to grg2. The C function pname must compute derivatives of the constraint functions and objective function (nfuns functions) with respect to variables x[j] and return them in the array grad so that grad[i*(nvars+1)+j] = Dg[i]/Dx[j] for i = 1,..., nfuns and j = 1,..., nvars.
The C function that computes derivatives for the example problem is named h11P and is listed in Section 2.2 C Program for Example Problem. It is passed to GRG2 by calling goptP(h11P) prior to calling grg2.
The Fortran interface for subroutine goptpf is:
subroutine goptpf (pname)
external pname
subroutine pname (g, x, nfuns, nvars, grad)
interger nfuns, nvars
double precision g(nfuns), x(nvars), grad(nfuns, nvars)
The user provides a Fortran subroutine to GRG2 by inserting a statement of the form,
call goptpf(pname);
before the call to grg2f. The Fortran subroutine pname must compute derivatives of the constraint functions and objective function (nfuns functions) with respect to variables x(j), j = 1,..., nvars and return them in the array grad so that grad(i,j) = Dg(i)/Dx(j) for i = 1,..., nfuns and j = 1,..., nvars.
The Fortran subroutine that computes derivatives for the example problem is named h11Pf and is listed in Section 2.3 Fortran Program for Example Problem. It IS passed to GRG2 by calling goptpf(h11Pf) prior to calling grg2f.
The MATLAB interface for function goptpm is:
goptpm ('pname')
[grad] = pname (g, x, nfuns, nvars)
The user provides a MATLAB function to GRG2 by inserting a statement of the form,
goptpm('pname');
before the call to grg2m. The MATLAB function pname must compute derivatives of the constraint functions and objective function (nfuns functions) with respect to variables x(j), j = 1,..., nvars and return them in the array grad so that grad(i,j) = Dg(i)/Dx(j) for i = 1,..., nfuns and j = 1,..., nvars.
The Fortran subroutine that computes derivatives for the example problem is named hmbl11 and is listed in Section 2.4 MATLAB Program for Example Problem. It IS passed to GRG2 by calling goptpm(hmbl11p) prior to calling grg2m.
The function goptL is used to set labels for variables and for constraint functions. These labels are printed in the GRG2 report and help make it easier to associate GRG2 variables and functions with problem variables and functions. Labels are limited to 10 characters. If goptL is not used, variables are labeled X(j) and functions are labeled G(i) in the report.The C prototype for goptL is given below with a C program example.
TYPE_export(void) goptL (
int nvnames,
char varname[][11],
int nfnames,
char funname[][11]);
Example:
char vnames[NVARS][11], fnames[NFUNS][11];
int i, j;
for (j = 1; j <= nvars; j++) strcpy (vnames[j], "VarName");
for (i = 1; i <= nfuns; i++) strcpy (fnames[i], "FunName");
goptL(nvars, vnames, nfuns, fnames);
grg2(nvars, ...);
The ability to set labels is not available in the Fortran and MATLAB GRG2 interfaces.
GRG2 produces a report that summarizes the problem starting values, the optimization process, and the final results. The report is written to the file specified by the report argument. Most of the information in the report is self-explanatory, but some of the terms require an understanding of the Generalized Reduced Gradient (GRG) algorithm. The full report for the example problem is given below and serves to illustrate the terms defined here. The reader is referred to the book Linear and Nonlinear Programming by David Luenberger, Reference 10., for a complete explanation of the GRG algorithm.
The report is divided into five sections: Problem Description, Starting Values, Solution Process, Final Results, and Summary. These five sections are discussed below.
Problem Description:
The Problem Description section displays the version of GRG2 used along with the date and time of the run. The number of variables and the number of functions (this number is the number of constraints plus 1 for the objective function) are displayed along with an indication of whether this is a minimization or maximization problem. If algorithmic parameters have been changed by the user, all parameter values are listed here. When all default values are used, the parameter values are not listed. The information in this section should confirm that the problem has been setup correctly.
GRG232 97.06 Report: Date Mon Jun 02 12:00:00 1997 Problem Title: GRG2 C Example: Himmelblau Problem 11 Problem Parameters: Number of variables is 5 Number of functions is 4 Objective function will be MINimized
Starting Values:
This section lists the starting values for the problem variables and the computed values of the constraint functions and the objective function. The functions (constraints) table shows the values computed using the starting values of the variables. The notation used in the Status column, for variables and constraints, and in the Type column for constraints, is described in the following table.
|
Status |
Description |
Type |
Description |
||
|
UL |
Constraint or variable is at its upper bound |
EQ |
Equality constraint |
||
|
LL |
Constraint or variable is at its lower bound |
LE |
Upper bound constraint |
||
|
EQ |
Equality constraint |
GE |
Lower bound constraint |
||
|
**** |
Constraint is violated, value is less than lower bound or greater than upper bound |
RNGE |
Lower and upper bound constraint |
||
|
FREE |
Variable has no lower bound and no upper bound |
OBJ |
Objective function |
||
|
FX |
Fixed variable, lower bound equals upper bound |
NA |
Constraint function ignored because it has no lower bound and no upper bound |
||
Starting Values
Functions:
Function Initial Lower Upper
No. Name Status Type Value Bound Bound
___ __________ ______ ____ __________ __________ __________
1 G RNGE 90.1116 0 92
2 G RNGE 96.1674 90 110
3 G **** RNGE 16.7629 20 25
4 G OBJ -32217.4
Variables:
Variable Initial Lower Upper
No. Name Status Value Bound Bound
___ __________ ______ __________ __________ __________
1 X LL 78 78 102
2 X LL 33 33 45
3 X LL 27 27 45
4 X LL 27 27 45
5 X LL 27 27 45
Solution Process:
The solution process is iterative. It begins at the starting values and ends at the final values. GRG2 attempts to improve the objective function value at each iteration. If the problem is not feasible, phase I may result in an objective function that is worse than the starting value, but this is necessary in order to achieve feasibility. The iteration log shows the iteration number in column 1. Information in columns 2 to 9 is described in the following table.
|
Column |
Heading |
Description |
|
1 |
Itn |
Iteration number |
|
2 |
Objective Function |
Value of the objective function. If some constraints are not feasible (see column 5) the objective function value is the sum of the constraint violations. |
|
3 |
Binding |
Number of constraints that are at their lower bound or at their upper bound. The binding constraints provide a set of nonlinear equations that can be used to solve for some of the variables in terms of other variables. See Super Basics. |
|
4 |
Super |
Variables are classified by GRG2 as basic, nonbasic, and superbasic. See the Final Results section for a description. The superbasic variables are the independent variables of the reduced problem. |
|
5 |
Infeas |
Number of infeasible constraints, i.e., constraints that violate their bounds for the current values of the variables. |
|
6 |
Norm of |
Norm of the reduced gradient also referred to as the Kuhn-Tucker factor or the K-T factor. This number is the maximum absolute value of gradient of the reduced objective function with respect to the set of superbasic variables scaled by the value of the variable and the inverse of the value of the objective function. When the value of the K_T factor is less than epstop, the stopping criteria, GRG2 terminates the iteration and returns the current variable values. |
|
7 |
Hessian |
Condition number of the Hessian. The Hessian is the matrix of second derivatives of the objective function with respect to the variables. A Large condition number indicates an ill-conditioned problem that may be difficult to solve accurately. |
|
8 |
Step |
Step size for the current iteration. |
|
9 |
Degen |
T in this column indicates a degenerate step(iteration) otherwise this column is blank. When GRG2 reclassifies variables from basic to superbasic and from superbasic to basic, the objective function does not change and the iteration is flagged as degenerate to avoid the appearance of convergence. |
Itn Objective Binding Super Infeas Norm of Hessian Step Degen No. Function Constrs Basics Constr Red.Grad Cond.No. Size Step ___ _________ _______ ______ ______ ________ ________ ______ _____ 0 3.2371 0 4 1 2.3 1 0 1 -28906 0 4 0 0.48 1 0.25 2 -28971 1 3 0 0.16 1 0.0044 3 -30146 1 3 0 0.096 2.3 0.2 4 -30211 1 2 0 0.096 1 0.013 5 -30400 2 1 0 0.031 1 0.042 6 -30666 2 0 0 0 1 0.22
Final Results:
This section lists the final values for the problem variables and the computed values of the constraint functions and the objective function. The functions table shows the values of the constraint functions computed using the final values of the variables.
The GRG algorithm uses the binding constraints to solve for some variables, classified as basic, in terms of the remaining variables thus reducing the number of independent variables. The variables that are not basic are classified as nonbasic if they are at one of their bounds and superbasic otherwise. The classification is given in the Status column of the Final Results table. The reduced objective function is the objective function treated as a function of the nonbasic and superbasic variables. The gradient of the reduced objective function is called the reduced gradient.
The Reduced Gradient column gives the value of the reduced gradient for each nonbasic and each superbasic variable. Reduced gradient values for basic variables are zero, by definition. Since nonbasic variables are at one of their bounds, the reduced gradient with respect to a nonbasic variable should have the "right" sign. For example, when the problem is a minimization and the nonbasic variable is at its lower bound the reduced gradient should be greater than or equal to zero. This is an indication that a small increase in the value of the variable will result in an increase in the value of the objective (i.e., a move away from the optimum). The reduced gradient of a nonbasic variable that is at its upper bound should be less than or equal to zero. The scale of the variable and the scale of the objective function effect reduced gradient values. This is why GRG2 scales these reduced gradient values by the value of the variable and by the reciprocal of the value of the objective function and then computes the K-T factor as the maximum over all the superbasics. The K-T factor is printed in the solution process section of the report in the Norm of Red.Grad column. A small K-T factor is a good indication that a local optimum has been located. If this value is less than or equal to epstop, the stopping criteria has been satisfied and GRG2 terminates with inform=0. GRG2 also terminates, with inform=1, when the objective function converges to a relative error of epstop for nstop consecutive iterations.
Final Results
Functions:
Distance
Initial Final from Lagrange
No. Name Value Value Status Nearest Multiplier
Bound
___ __________ _________ ___________ __________ __________ ____________
1 G 90.112 92 UpperBnd 3.11e-006 :U -403.27
2 G 96.167 98.84 Free 8.84 :L
3 G 16.763 20 LowerBnd -8.08e-008 :L 809.43
4 G -32217 -30666 Objective
Variables:
Distance
Initial Final from Reduced
No. Name Value Value Status Nearest Gradient
Bound
___ __________ _________ ___________ ________ __________ __________
1 X 78 78 NonBasic LowerBnd 0.124
2 X 33 33 NonBasic LowerBnd 0.0907
3 X 27 29.995 Basic 2.995 :L
4 X 27 45 NonBasic UpperBnd -0.0391
5 X 27 36.776 Basic 8.224 :U
Summary:
This section of the GRG2 report provides a brief summary of the run. The minimized or maximized objective function value is given along with the GRG2 termination message. The messages are defined in Section 3.1.5 "Termination Messages (inform)". The number of calls to GCOMP, for evaluating constraint functions and the objective function is also given. If a user provided routine computes derivatives, parshP, the summary gives the number of times this routine was called. Finally, the summary section displays the time used by the run.
The termination message provides valuable information on the final results returned be GRG2. The message that follows inform=0 gives the K-T factor that has been described above. This is the best indication that GRG2 has been successful in locating a local optimum of the problem to the desired accuracy. When the termination message is inform=1, GRG2 also gives the K-T factor, which may be small enough to accept the final results. If the K-T factor is not viewed as being small enough, the user may restart GRG2 from the final results, but with a smaller value of epstop. This process usually results in a smaller K-T factor. If it does not, it may indicate that a better result is not possible (due to inaccurate function values or due to inaccurate derivative values).
Summary
MINimized objective function value is -30665.5
Termination: INFORM = 0
Kuhn-Tucker conditions are satisfied to
within 0 for the current variable values.
Relative change in the objective function value
is 0.0087 for the last iteration.
Number of function evaluations 54
Number of calls to user provided PARSH 7
Time used is 0.165 seconds.
It is possible for GRG2 to terminate at a point that is not optimal. If X is the vector of n variables when GRG2 terminates (output value of grg2 argument xx), then one of the following conditions must hold:
1. X closely approximates a global optimum,
2. X closely approximates a local optimum, but not a global optimum, or
3. X does not closely approximate even a local optimum.
These statements apply in both Phases I and II (see Section 1. "Introduction" for an explanation). This section deals with how to tell which of the above conditions hold and what to do if X is not a global optimum.
To aid in spotting errors in the input data and to see if the functions in GCOMP are coded correctly, all data values are written to the report file at the start of a run. If a problem is being solved for the first time, it may be useful to make a run that prints the input data and then stops. To do this, use goptD("limser", 0);
If derivatives are computed by user-supplied functions (see Section 3.3.1 "Providing Derivatives Through goptP") it may be that these functions contain errors. If goptD("ckgrad", 2) is used, the user-supplied derivatives, at the initial point, are compared with numerical estimates. Any differences are written to the report file. Large differences indicate errors in the user-supplied functions and GRG2 terminates.
GRG2 uses finite difference approximations to compute derivatives of the problem functions. These are computed by GRG2 routine parshf using forward differences or parshc using central differences. Use of parshf is the default option, while parshc is invoked by goptD("kderiv" , 1); See Section 3.2.4 "Descriptions of GRG2 Algorithmic Parameters".
Finite difference derivatives have roughly half the precision of the functions values g[i] computed in GCOMP. Hence, if each function is accurate to from 13 to 16 significant digits, then the derivatives have about 7 or 8 significant digits, which is adequate in most instances. However, if each function has 8 or fewer significant digits, then the derivatives have 4 or less, which may seriously hinder the optimizer. Low precision in GCOMP often occurs when numerical routines are used to evaluate one or more of the functions. This can occur in (a) chemical process models where recycle loops require iterative solution of a system of implicit nonlinear equations, or (b) when differential equations are solved numerically within GCOMP. Then, the accuracy of these numerical calculations determines the accuracy of the functions.
Inaccurate derivatives can cause GRG2 to terminate at a nonoptimal point, often with the message inform=2, "All remedies have failed to find a better point. Optimization process terminated" (see Section 3.1.5"Termination Messages" for an explanation). Possible remedies include (a) increasing the accuracy of the functions computed in GCOMP, (b) trying central differences, goptD("kderiv", 1), (c) trying different stepsizes, goptD("pstep", step_size), or (d) using analytic derivatives computed by a use-supplied function (see Section 3.3.1 "Providing Derivatives Through goptP"). Any of these options should be tested. See Section 4.1 "Errors in the Input Data or in GCOMP" for how to do this.
Proper scaling of both variables and problem functions is very important for successful operation of GRG2. It is recommended that all problem functions be scaled to have absolute values less than 100. In addition, their values should be large enough so that a 1.0e-4 error in computing them is not significant (the default value of epnewt is 1.0e-6, see Section 3.2.4 "Descriptions of GRG2 Algorithmic Parameters".
Variables should be scaled so that a unit change (changes of 1.0) represents a small but significant change in that variable. In addition, it is advisable to avoid having the constraint or objective functions much more sensitive to some variables than others. A symptom of bad scaling is the presence of very large derivative values. The "doscale" option provides one type of scaling that is based on gradient values. Use goptD("doscale", 1) to invoke this option.
It is possible for GRG2 to call GCOMP to be called with some x[j] slightly outside of its bounds. If this is likely to cause underflow or overflow or an error termination (such as attempting to take the square root or log of a negative number), preventive action must be taken in GCOMP. If, for example, x[j] all have positive lower bounds and there is a division by some of the x[j], a loop should be set up at the start of GCOMP comparing x[j] with a small positive number less than the lower bounds. If any x[j] is less than this value, set all the effected functions g[i] to an arbitrarily large number such as 1.e30 and return. This will cause the stepsize to be cut back.
Neither GRG2 nor any other nonlinear optimization package can guarantee finding a global optimum in cases where there are distinct local optima. If you know that your problem is convex (minimizing a convex objective over a convex feasible region), then any local optimum is global, so this problem cannot occur. Otherwise, one must use a combination of (a) knowledge about the problem and (b) a variety of starting points to help decide the local/global issue. If all starting points yield approximately the same final point, and that point is satisfactory to persons knowledgeable about the problem, one has a high level of confidence that the point found is (globally) optimal.
There are several ways to assess the accuracy of a solution.
In (2) above, some of the nonbasic variables can be fixed at their current values (see Section 3.1.3 "Arguments of grg2" for how to do this) and let GRG2 vary the others. If no variation produces a significantly improved feasible point, confidence in the current solution is increased. Confidence is also increased if different starting points lead to nearly the same final point. In (4), the best indication of accuracy occurs if the Kuhn-Tucker conditions are satisfied (inform=0). However, even if they are not, but the reduced gradient components of superbasic variables have been reduced from their initial values, by say 3 orders of magnitude or more, while reduced gradients of nonbasic variables at bound have the correct sign, this is symptomatic of reasonably high accuracy.
If perceived accuracy is not sufficient, reducing epstop and/or increasing nstop (see Section 3.2.4 "Descriptions of GRG2 Algorithmic Parameters") often helps. It is often helpful to reduce epnewt if epstop is reduced, since this increases the accuracy of the reduced gradient computation. Other parameter changes and/or option choices may also be helpful. See Section 5. "GRG2 Tolerances and Algorithmic Options".
All GRG2 tolerances and algorithmic options have default values and need not be set by the user. However, manipulation of these tolerances and options can sometimes improve the performance of GRG2. All parameters can be set with goptD. See see Section 3.2.4 "Descriptions of GRG2 Algorithmic Parameters".
The default parameter settings are appropriate when function values are accurate to a precision of 1.0e-12 or more. This is the case when function involve only algebraic operations (+, -, *, /) and other functions (sin, cos, exp) that can be computed accurately. When functions are complex and cannot be evaluated accurately it is necessary to set GRG2 parameters, epstop, epnewt, and epinit to values larger than their default values. It may be necessary to try a range of values to determine which ones work best for a given problem.
The most critical tolerance is epnewt. Increasing it can sometimes speed convergence (by requiring fewer Newton iterations) while decreasing it occasionally yields a more accurate solution or gets the iterations moving if the algorithm gets "stuck". Any values larger than 1.0e-2 should be treated cautiously, as should any value smaller than 1.0e-6. Choosing a value for epinit different from epnewt has helped solve a few problems that were not solved otherwise. We have used epinit = 1.0e-4, epnewt = 1.0e-6. Choosing a smaller value for epstop (or a larger value for nstop, Section 3.1 "GRG2 C Function Interface" describes these tolerances) usually improves the accuracy of the final solution. If the problem is degenerate and this is slowing computations, choosing a larger value for epspiv may help by allowing pivots on elements that were previously rejected. If convergence of the iterations is a problem, reducing epspiv and/or increasing epnewt may help. A nonzero value of eph1ep will cause Phase I to consider the true objective along with the sum of infeasibilities and may yield a better point at the end of Phase I.
Central differences, selected with goptD("kderiv", 1), are more accurate than forward differences. Central differences are exact for quadratic functions, while forward differences are exact only for linear functions. However, central differences require two function evaluations per derivative, while forward differences require only one. For further discussion, see Section 4.3 "Inaccurate Numerical Derivatives".
Quadratic extrapolation can often speed computations by providing better initial values for the iterations. This option is selected by goptD("iquad", 1). It is unnecessary if all constraints are linear.
The conjugate gradient (CG) method is useful in problems with many variables. It generally reduces memory requirements at the expense of more iterations and more computer time than the default quasi-Newton option. Use goptD("maxr", 0) and goptD("modcg", n) with 1£ n£ 5 to select one of the CG methods.
GRG2 can be used for unconstrained minimization or minimization subject only to upper and lower bounds. Simply compute the objective function in GCOMP. The program will then minimize the objective using the Broyden-Fletcher-Goldfarb-Shanno (BFGS) variable metric method [8] modified to deal with upper and lower bounds on the variables.
The program can also solve discrete-time optimal control problems. To do this, let nvars be the total number of control variables (over all time periods) and let nfuns be the total number of problem constraints other than bounds on the control, e.g., state variable inequality or terminal constraints. If there are none of these use nfuns = 1. Function GCOMP must be written so as to accept the control vector X, solve the system equations for the state variables, and evaluate the objective and constraint functions. If finite difference derivatives are to be used, the state variables need not appear anywhere except in GCOMP; GRG2 itself will be unaware that they exist. If analytic partial derivatives are used, then the state variables will be required by the user-supplied function that computes the derivatives.
1. Lasdon, L.S., Waren, A.D., Jain, A., and Ratner, M., "Design and Testing of a Generalized Reduced Gradient Code for Nonlinear Programming", ACM Transactions on Mathematical Software, Vol. 4, No. 1, March 1978, pp. 34-50.
2. Lasdon, L.S. and Waren, A.D., "Generalized Reduced Gradient Software for Linearly and Nonlinearly Constrained Problems", in "Design and Implementation of Optimization Software", H. Greenberg, ed., Sijthoff and Noordhoff, pubs, 1979.
3. Abadie, J. and Carpentier, J. "Generalization of the Wolfe Reduced Gradient Method to the Case of Nonlinear Constraints", in Optimization, R. Fletcher (ed.), Academic Press, London; 1969, pp. 37-47.
4. Murtagh, B.A. and Saunders, M.A. "Large-scale Linearly Constrained Optimization", Mathematical Programming, Vol. 14, No. 1, January 1978, pp. 41-72.
5. Powell, M.J.D., "Restart Procedures for the Conjugate Gradient Method," Mathematical Programming, Vol. 12, No. 2, April 1977, pp. 241-255.
6. Colville, A.R., "A Comparative Study of Nonlinear Programming Codes," I.B.M. T.R. no. 320-2949 (1968).
7. Himmelblau, D.M., Applied Nonlinear Programming, McGraw-Hill Book Co., New York, 1972.
8. Fletcher, R., "A New Approach to Variable Metric Algorithms", Computer Journal, Vol. 13, 1970, pp. 317-322.
9. Smith, S. and Lasdon, L.S., Solving Large Sparse Nonlinear Programs Using GRG, ORSA Journal on Computing, Vol. 4, No. 1,Winter 1992, pp. 1-15.
10. Luenbuerger, David G., Linear and Nonlinear Programming, Second Edition, Addison-Wesley, Reading Massachusetts, 1984.
Windward Technologies, Inc.
12039 Mulholland
Meadows Place, TX 77477-1620
Phone: 281-564-6523
Fax: 281-754-4022
E-mail: wti@aol.com
Web site: http://web.wt.net/~wti
Leon S. Lasdon
Department of Management Science and Information Systems
University of Texas
Austin, TX 78712
Allan D. Waren
Department of Computer and Information Science
College of Business Administration
Cleveland, OH 44115
Windward Technologies, Inc. and Optimal Methods, Inc. make no warranty of any kind with regard to this material, including, but not limited to, the implied warranties of merchantability and fitness for a particular purpose. Neither Windward Technologies nor Optimal Methods shall be liable for errors contained herein or for incidental, consequential, or other indirect damages in connection with the furnishing, performance, or use of this material.
Copyright ã 1979-1997 by Windward Technologies, Inc. and Optimal Methods, Inc.
All rights reserved
. No part of this document may be photocopied or reproduced without the prior written consent of Windward Technologies, Inc.Trademarks:
Microsoft and MS are trademarks of Microsoft Corporation;Use, duplication, or disclosure by the US Government is subject to restrictions as set forth in FAR 52.227-19, subparagraph (c)(i)(ii) of DOD FAR SUPP 252.227-7013, or equivalent government clause for other agencies.
Restricted Rights Notice: The version of the WTI product described in this document is sold under a per user license agreement. Its use, duplication, and disclosure are subject to the restrictions in the license agreement.
Revised: June 24, 1997