The main program and the model defining procedures can now be created. They are for our small example tutorial.f90 (tutorial.c).

The Main Program

The main program can in most cases be used directly or with few changes. Before the optimization routine of CONOPT (COI_Solve) is called you must go through the steps mentioned in section 5 of the Reference Manual: Create and initialize a Control vector, define the size of the model, define the Base, i.e. select between Fortran and C calling conventions, define the objective, and register the callback routines. The declarations and statements can almost be used without change in most models. Only the actual sizes and the names of the callback routines in the COIDEF_Xx calls should be changed. Note that all arguments to CONOPT calls in general are passed by reference. You must therefore define local (or global) variables with the sizes, initialize these variables, and pass their address to CONOPT.

Some of the callback routines are messaging routines, where it often is possible to use some standard routines. The file std.f90 contains standard versions of the Message, ErrMsg, Status, and Solution routines with the names Std_Message, Std_ErrMsg, Std_Status, and Std_Solution. We will use these routines by including std.f90, and by declaring the four routines as External. The Std_Message routine feeds three files (see section 11.1 in the Reference Manual): The Screen file is represented by unit *, the documentation file is associated with unit 10, and the status file is associated with unit 11. The last two units are therefore opened in the main program before calling COI_Solve.

Similarly, for the C interface, the file std.c contains standard versions of the Message, ErrMsg, Status, and Solution routines with the names Std_Message, Std_ErrMsg, Std_Status, and Std_Solution. These routines can be used by including std.c, and by registering these callback routines. As above, the Std_Message routine feeds three files; however, in the C interface the Screen file is represented by Stdout, the documentation file is associated with file *fd, and the status file is associated with file *fs. These file handles are global variables, and the files are opened in the main program before calling COI_Solve.

Once CONOPT has been started, it will call a number of callback routines that must be supplied and registered by the user. In addition to the messaging routines these include a ReadMatrix and FDEval routine. The implementations of these callback routines are for this model called Tut_ReadMatrix and Tut_FDEval. Tut_ReadMatrix provides CONOPT with detailed information about the model and Tut_FDEval evaluates the nonlinear terms of the model and their derivatives. Tut_ReadMatrix is called once and Tut_FDEval is called repeatedly during the optimization. Tut_ReadMatrix and Tut_FDEval are described next.

The problem data and matrix reading callback – Tut_ReadMatrix

Based on the information provided in the control vector, CONOPT will allocate a number of vectors into which the user must define the details of the model. This must be done in the user written function called Tut_ReadMatrix. The argument list can be seen in Appendix A. The usage, type, and default content of the individual vectors are described below. For the variables we define:

lower (LOWER): lower bounds (real*8 (double), default: -Inf),
curr (CURR): current or initial values (real*8 (double), default: 0),
upper (UPPER): upper bounds (real*8 (double), default: +Inf),
vsta (VSTA): variable status (integer, (int), default; not used) for each variable.

In our model, the default value of lower of -Inf (minus infinity or unbounded from below) cannot be used so we must define lower for all four variables. To show the use of default initial values we provide non-default initial values for X(1) and X(2) only. In practice, it is highly recommended to provide as many and as good initial values as possible. It is particularly important when variables appear in products (such as here). We use the default upper bounds of +Inf so upper is not defined. By default, vsta is not used and is therefore not defined in our example (See COIDEF_IniStat for details on status values).

For the constraints we define:

type (TYPE): constraint type (integer (int), see the codes below),
rhs (RHS): the right hand sides (real*8 (double), default: 0),
esta (ESTA): equation status (integer (int), default: not used) for each equation.

The type vector must contain integer codes as follows:

0: equal, \(g(x) = RHS\),
1: greater than or equal, \(g(x) \geq RHS\),
2: less than or equal, \(g(x) \leq RHS\),
3: non-constraining, \(g(x) =\!\!N\!\!= RHS\).

In our model, we have one non-constraining equation, namely the objective at index 1, and two equalities. Only constraint 1 and 3 have nonzero right hand sides so only these two elements of rhs are defined. By default esta is not used and is therefore not defined in our example.

rowno (ROWNO): row numbers (integer (int)),
value (VALUE): derivative values (real*8 (double)),
nlflag (NLFLAG): nonlinearity flags (integer (int)) for each derivative or Jacobian elements.

Only nonzero derivatives should be mentioned, and the derivatives must be sorted column-wise so all elements from column 0 appear first followed by all elements from column 1, then column 2, and finally column 3. Within each column the elements may appear in any order; the rowno vector tells CONOPT where each derivative belong. The nlflag vector tells CONOPT which derivatives are linear (nlflag = 0) and which are nonlinear (nlflag = 1). And the value vector provides CONOPT with the values of the derivatives. CONOPT will only use the linear derivatives, i.e. those with nlflag = 0.

colsta (COLSTA): Vector of start of column indices (integer (int)) into the rowno, value, and nlflag vectors.

colsta points to the start of each column within the sorted order. An extra element at position [Number of variables] must state the number of nonzeros.

The values of colsta, rowno, nlflag and the values of value for the linear derivatives are all derived from the Jacobian. The details can be read from the comments in Tut_ReadMatrix in tutorial.f90 or tutorial.c.

The remaining arguments in Tut_ReadMatrix define the sizes of the vectors:

n (*N): number of variables,
m (*M): number of constraints, and
nz (*NZ): number of nonzeros.

And usrmem (USRMEM) can be used to communicate between the users subroutines.

The return value of Tut_ReadMatrix can be used to communicate errors to CONOPT. Return 0 if all went well, and return a nonzero value if there was an error during model setup. If Tut_ReadMatrix returns a nonzero value then CONOPT, or rather COI_Solve, will return immediately.

The first derivative evaluation callback – Tut_FDEval

The last callback routine that must be provided by the user is a function evaluation routine, here called Tut_FDEval. This routine provides CONOPT with numerical values of the nonlinear terms and numerical values of their derivatives. Seen from CONOPT, Tut_FDEval is a black box: CONOPT provides a vector of values for the variables, x, plus some control information, and it expects to get function and/or derivative information back.

The complete argument list for Tut_FDEval can be seen in tutorial.f90 or tutorial.c. For most purposes, only x (X), g (G), jac (JAC), rowno (ROWNO), and mode (MODE) are needed so we will only discuss these five arguments here. Detailed information about the other arguments can be found in the Reference Manual:

Argument	Description
`x` (`X`)	A `real*8` (`double`) vector with the current values of the variables. This vector is provided by CONOPT and must not be changed by the user. The values reported in this vector are guaranteed by CONOPT to be within the upper and lower bounds
`mode` (`*MODE`)	An integer scalar control parameter. It is defined by CONOPT and it must not be changed by the user. `mode == 1` instructs the user to compute the value of the nonlinear term in constraint number `rowno` and to return it in `g`. `mode == 2` instructs the user to compute the nonzero derivatives of the nonlinear terms in constraint number `rowno` and to return them in `jac`. `mode == 3` is equivalent to both `mode == 1` and `mode == 2`.
`rowno` (`*ROWNO`)	`Tut_FDEval` is called for one row or constraint at a time. The integer `rowno` tells the user which row to work on in the current call. `rowno` is defined by CONOPT and must not be changed by the user. Note that `Tut_FDEval` only will be called for nonlinear constraints so `rowno` will never take values corresponding to linear constraints. In our model, `Tut_FDEval` will never be called with `rowno == 2`. (This behavior can be changed with various options defined using `COIDEF_FVforAll` and `COIDEF_FVincLin` to allow `Tut_FDEval` to compute constant and/or linear terms).
`g` (`*G`)	If `mode == 1` or `mode == 3`, then the user should return the value of the nonlinear terms in the `real*8` (`double`) scalar `g`. If `mode == 2` then CONOPT will not use `g`. Note that `g` is a scalar. There is no index corresponding to `rowno`.
`jac` (`*JAC`)	If `mode == 2` or `mode == 3`, then the user should return the derivatives of the nonlinear terms in constraint number `rowno` in the `real*8` (`double`) vectors `jac`. The derivative with respect to variable `x(i)` should be returned in `jac(i)`. Based on information from `Tut_ReadMatrix` CONOPT will only extract the relevant nonlinear derivatives from `jac` so it is not necessary to set the remaining positions to zero.

The constraints and their derivatives could be programmed using x(0), x(1), ..., directly to denote the variables. In the example, we have created local variables with the names L, Inp, Out, and P and we have moved the values from the X-vector into these four local variables to make Tut_FDEval more readable. All later statements use these local variables. We have also declared all parameters as local variables with the same names as those used in the mathematical model formulation and assigned them their proper values to make the model easier to read and easier to update.

The calculations are separated into one block for each constraint. Within each constraint, they are separated into two parts, one for constraint values and one for derivative values. rowno is tested to determine the constraint, and mode is tested to determine whether to compute constraint values or derivative values or both. Note, that it is illegal to change the derivative vector jac when mode == 1 and it is not advised to change the constraint g when mode == 2.

Only the nonlinear parts of each constraint as defined in the Constraint Reorganization section is included in Tut_FDEval. The nonlinear parts of the objective function was identified as

\[ Out P \]

so the only term that is included in Tut_FDEval is \(Out P\), and constraint 1 is therefore programmed as g = P * Out. The derivative with respect to \(Out\) is \(P\) and it is placed in jac(3) because Out is variable 3. Similarly, the derivative with respect to P is Out and it is placed in jac(4) because P is variable 4. The constraint and/or its derivatives are evaluated when rowno is 1.

The second constraint is programmed in the same way. In our example, we have calculated some terms that the function and its derivatives have in common. The function value is

\begin{equation} \begin{aligned} \texttt{g} &= (A_{l} L^{-\rho} + A_{k} K^{-\rho} + A_{Inp} Inp^{-\rho})^{-1/\rho}\\ &= \texttt{Hold1}(L, Inp)^{-1/\rho} \\ &= \texttt{Hold2} \end{aligned} \end{equation}

where the argument for the outer power function is defined as Hold1. The derivative with respect to \(L\) (variable 0) is

\begin{equation} \begin{aligned} \texttt{jac(1)} &= \frac{d\texttt{G}}{d\texttt{Hold1}} \frac{d\texttt{Hold1}}{dL}&& \\ &= \frac{-1}{\rho}\texttt{Hold1}^{-1/\rho}&&\times-\rho A_{l}L^{-\rho - 1} \\ &= \frac{-1}{\rho}\frac{\texttt{Hold2}}{\texttt{Hold1}}&&\times-\rho A_{l}L^{-\rho - 1} \\ &= \frac{\texttt{Hold2}}{\texttt{Hold1}}&&\times A_{l}L^{-\rho - 1}\\ &= \texttt{Hold3}&& \times A_{l}L^{-\rho - 1} \end{aligned} \end{equation}

and the derivative with respect to the \(Inp\) (variable 1) is similarly

\begin{equation} \begin{aligned} \texttt{jac(2)} &= \frac{d\texttt{G}}{d\texttt{Hold1}} \frac{d\texttt{Hold1}}{dInp}&& \\ &= \frac{-1}{\rho}\texttt{Hold1}^{-1/\rho}&&\times-\rho A_{l}Inp^{-\rho - 1} \\ &= \frac{-1}{\rho}\frac{\texttt{Hold2}}{\texttt{Hold1}}&&\times-\rho A_{l}Inp^{-\rho - 1} \\ &= \frac{\texttt{Hold2}}{\texttt{Hold1}}&&\times A_{l}Inp^{-\rho - 1}\\ &= \texttt{Hold3}&& \times A_{l}Inp^{-\rho - 1} \end{aligned} \end{equation}

Hold1 and Hold2 are computed initially and used when the function and derivatives themselves are calculated and stored in g and/or jac. This programming practice can make the calculations of the derivatives simpler and more efficient. In general, you can arrange the computations inside Tut_FDEval in any way you find convenient as long as the values returned in g and jac are correct. You could for example read values for the parameters from a file, or you could call other subroutines, both subroutine written specifically for the model at hand and general-purpose routines.

The derivation of the derivatives shown about was lengthy and therefore also error prone. It is a good idea to test that the function values and derivatives are consistent by turning the Function and Derivative Debugger on during the first runs. This is done by calling COIDEF_DebugFV. You may also consider writing the nonlinear expressions in a simpler form. In the current model, you can split the complicated equation into two by introducing an intermediate variable. This is discussed again in second derivatives.

Note, that there is no code for constraint number 3. This constraint is linear and is completely defined through subroutine Tut_ReadMatrix. Tut_FDEval will therefore never be called with rowno == 3.

The return value of FDEval should usually be zero. A nonzero return value will tell CONOPT that something is seriously wrong and CONOPT, or rather COI_Solve, will return immediately.