DBOCLS(3)

SYNOPSIS
SUBROUTINE DBOCLS (W, MDW, MCON, MROWS, NCOLS, BL, BU, IND, IOPT, + X, RNORMC, RNORM, MODE, RW, IW)

NAME
dbocls - Solve the bounded and constrained least squares problem consisting of solving the equation E*X = F (in the least squares sense) subject to the linear constraints C*X = Y.

LIBRARY
SLATEC

CATEGORY
K1A2A, G2E, G2H1, G2H2

TYPE
DOUBLE PRECISION (SBOCLS-S, DBOCLS-D)

KEYWORDS
BOUNDS, CONSTRAINTS, INEQUALITY, LEAST SQUARES, LINEAR

AUTHOR
Hanson, R. J., (SNLA)

DESCRIPTION
**** All INPUT and OUTPUT real variables are DOUBLE PRECISION **** This subprogram solves the bounded and constrained least squares problem. The problem statement is: Solve E*X = F (least squares sense), subject to constraints C*X=Y. In this formulation both X and Y are unknowns, and both may have bounds on any of their components. This formulation of the problem allows the user to have equality and inequality constraints as well as simple bounds on the solution components. This constrained linear least squares subprogram solves E*X=F subject to C*X=Y, where E is MROWS by NCOLS, C is MCON by NCOLS. The user must have dimension statements of the form DIMENSION W(MDW,NCOLS+MCON+1), BL(NCOLS+MCON), BU(NCOLS+MCON), * X(2*(NCOLS+MCON)+2+NX), RW(6*NCOLS+5*MCON) INTEGER IND(NCOLS+MCON), IOPT(17+NI), IW(2*(NCOLS+MCON)) (here NX=number of extra locations required for the options; NX=0 if no options are in use. Also NI=number of extra locations for options 1-9.) INPUT ----- ------------------------- W(MDW,*),MCON,MROWS,NCOLS ------------------------- The array W contains the (possibly null) matrix [C:*] followed by [E:F]. This must be placed in W as follows: [C : *] W = [ ] [E : F] The (*) after C indicates that this data can be undefined. The matrix [E:F] has MROWS rows and NCOLS+1 columns. The matrix C is placed in the first MCON rows of W(*,*) while [E:F] follows in rows MCON+1 through MCON+MROWS of W(*,*). The vector F is placed in rows MCON+1 through MCON+MROWS, column NCOLS+1. The values of MDW and NCOLS must be positive; the value of MCON must be nonnegative. An exception to this occurs when using option 1 for accumulation of blocks of equations. In that case MROWS is an OUTPUT variable only, and the matrix data for [E:F] is placed in W(*,*), one block of rows at a time. See IOPT(*) contents, option number 1, for further details. The row dimension, MDW, of the array W(*,*) must satisfy the inequality: If using option 1, MDW .ge. MCON + max(max. number of rows accumulated, NCOLS) + 1. If using option 8, MDW .ge. MCON + MROWS. Else MDW .ge. MCON + max(MROWS, NCOLS). Other values are errors, but this is checked only when using option=2. The value of MROWS is an output parameter when using option number 1 for accumulating large blocks of least squares equations before solving the problem. See IOPT(*) contents for details about option 1. ------------------ BL(*),BU(*),IND(*) ------------------ These arrays contain the information about the bounds that the solution values are to satisfy. The value of IND(J) tells the type of bound and BL(J) and BU(J) give the explicit values for the respective upper and lower bounds on the unknowns X and Y. The first NVARS entries of IND(*), BL(*) and BU(*) specify bounds on X; the next MCON entries specify bounds on Y. 1. For IND(J)=1, require X(J) .ge. BL(J); IF J.gt.NCOLS, Y(J-NCOLS) .ge. BL(J). (the value of BU(J) is not used.) 2. For IND(J)=2, require X(J) .le. BU(J); IF J.gt.NCOLS, Y(J-NCOLS) .le. BU(J). (the value of BL(J) is not used.) 3. For IND(J)=3, require X(J) .ge. BL(J) and X(J) .le. BU(J); IF J.gt.NCOLS, Y(J-NCOLS) .ge. BL(J) and Y(J-NCOLS) .le. BU(J). (to impose equality constraints have BL(J)=BU(J)= constraining value.) 4. For IND(J)=4, no bounds on X(J) or Y(J-NCOLS) are required. (the values of BL(J) and BU(J) are not used.) Values other than 1,2,3 or 4 for IND(J) are errors. In the case IND(J)=3 (upper and lower bounds) the condition BL(J) .gt. BU(J) is an error. The values BL(J), BU(J), J .gt. NCOLS, will be changed. Significant changes mean that the constraints are infeasible. (Users must make this decision themselves.) The new values for BL(J), BU(J), J .gt. NCOLS, define a region such that the perturbed problem is feasible. If users know that their problem is feasible, this step can be skipped by using option number 8 described below. See IOPT(*) description. ------- IOPT(*) ------- This is the array where the user can specify nonstandard options for DBOCLS( ). Most of the time this feature can be ignored by setting the input value IOPT(1)=99. Occasionally users may have needs that require use of the following subprogram options. For details about how to use the options see below: IOPT(*) CONTENTS. Option Number Brief Statement of Purpose ------ ------ ----- --------- -- ------- 1 Return to user for accumulation of blocks of least squares equations. The values of IOPT(*) are changed with this option. The changes are updates to pointers for placing the rows of equations into position for processing. 2 Check lengths of all arrays used in the subprogram. 3 Column scaling of the data matrix, [C]. [E] 4 User provides column scaling for matrix [C]. [E] 5 Provide option array to the low-level subprogram SBOLS( ). 6 Provide option array to the low-level subprogram SBOLSM( ). 7 Move the IOPT(*) processing pointer. 8 Do not preprocess the constraints to resolve infeasibilities. 9 Do not pretriangularize the least squares matrix. 99 No more options to change. ---- X(*) ---- This array is used to pass data associated with options 4,5 and 6. Ignore this parameter (on input) if no options are used. Otherwise see below: IOPT(*) CONTENTS. OUTPUT ------ ----------------- X(*),RNORMC,RNORM ----------------- The array X(*) contains a solution (if MODE .ge.0 or .eq.-22) for the constrained least squares problem. The value RNORMC is the minimum residual vector length for the constraints C*X - Y = 0. The value RNORM is the minimum residual vector length for the least squares equations. Normally RNORMC=0, but in the case of inconsistent constraints this value will be nonzero. The values of X are returned in the first NVARS entries of X(*). The values of Y are returned in the last MCON entries of X(*). ---- MODE ---- The sign of MODE determines whether the subprogram has completed normally, or encountered an error condition or abnormal status. A value of MODE .ge. 0 signifies that the subprogram has completed normally. The value of mode (.ge. 0) is the number of variables in an active status: not at a bound nor at the value zero, for the case of free variables. A negative value of MODE will be one of the cases (-57)-(-41), (-37)-(-22), (-19)-(-2). Values .lt. -1 correspond to an abnormal completion of the subprogram. These error messages are in groups for the subprograms DBOCLS(), SBOLSM(), and SBOLS(). An approximate solution will be returned to the user only when max. iterations is reached, MODE=-22. ----------- RW(*),IW(*) ----------- These are working arrays. (normally the user can ignore the contents of these arrays.) IOPT(*) CONTENTS ------- -------- The option array allows a user to modify some internal variables in the subprogram without recompiling the source code. A central goal of the initial software design was to do a good job for most people. Thus the use of options will be restricted to a select group of users. The processing of the option array proceeds as follows: a pointer, here called LP, is initially set to the value 1. At the pointer position the option number is extracted and used for locating other information that allows for options to be changed. The portion of the array IOPT(*) that is used for each option is fixed; the user and the subprogram both know how many locations are needed for each option. The value of LP is updated for each option based on the amount of storage in IOPT(*) that is required. A great deal of error checking is done by the subprogram on the contents of the option array. Nevertheless it is still possible to give the subprogram optional input that is meaningless. For example option 4 uses the locations X(NCOLS+IOFF),...,X(NCOLS+IOFF+NCOLS-1) for passing scaling data. The user must manage the allocation of these locations. 1 - This option allows the user to solve problems with a large number of rows compared to the number of variables. The idea is that the subprogram returns to the user (perhaps many times) and receives new least squares equations from the calling program unit. Eventually the user signals "that's all" and a solution is then computed. The value of MROWS is an output variable when this option is used. Its value is always in the range 0 .le. MROWS triangularization of the entire set of equations. If LP is the processing pointer for IOPT(*), the usage for the sequential processing of blocks of equations is IOPT(LP)=1 Move block of equations to W(*,*) starting at the first row of W(*,*). IOPT(LP+3)=# of rows in the block; user defined The user now calls DBOCLS( ) in a loop. The value of IOPT(LP+1) directs the user's action. The value of IOPT(LP+2) points to where the subsequent rows are to be placed in W(*,*). Both of these values are first defined in the subprogram. The user changes the value of IOPT(LP+1) (to 2) as a signal that all of the rows have been processed. <-------EXIT LOOP SOLUTION COMPUTED IF MODE .GE. 0 Use of this option adds 4 to the required length of IOPT(*). 2 - This option is useful for checking the lengths of all arrays used by DBOCLS( ) against their actual requirements for this problem. The idea is simple: the user's program unit passes the declared dimension information of the arrays. These values are compared against the problem-dependent needs within the subprogram. If any of the dimensions are too small an error message is printed and a negative value of MODE is returned, -41 to -47. The printed error message tells how long the dimension should be. If LP is the processing pointer for IOPT(*), IOPT(LP)=2 IOPT(LP+1)=Row dimension of W(*,*) IOPT(LP+2)=Col. dimension of W(*,*) IOPT(LP+3)=Dimensions of BL(*),BU(*),IND(*) IOPT(LP+4)=Dimension of X(*) IOPT(LP+5)=Dimension of RW(*) IOPT(LP+6)=Dimension of IW(*) IOPT(LP+7)=Dimension of IOPT(*) CALL DBOCLS( ) Use of this option adds 8 to the required length of IOPT(*). 3 - This option can change the type of scaling for the data matrix. Nominally each nonzero column of the matrix is scaled so that the magnitude of its largest entry is equal to the value ONE. If LP is the processing pointer for IOPT(*), IOPT(LP)=3 IOPT(LP+1)=1,2 or 3 1= Nominal scaling as noted; 2= Each nonzero column scaled to have length ONE; 3= Identity scaling; scaling effectively suppressed. CALL DBOCLS( ) Use of this option adds 2 to the required length of IOPT(*). 4 - This options allows the user to provide arbitrary (positive) column scaling for the matrix. If LP is the processing pointer for IOPT(*), IOPT(LP)=4 IOPT(LP+1)=IOFF X(NCOLS+IOFF),...,X(NCOLS+IOFF+NCOLS-1) = Positive scale factors for cols. of E. CALL DBOCLS( ) Use of this option adds 2 to the required length of IOPT(*) and NCOLS to the required length of X(*). 5 - This option allows the user to provide an option array to the low-level subprogram SBOLS( ). If LP is the processing pointer for IOPT(*), IOPT(LP)=5 IOPT(LP+1)= Position in IOPT(*) where option array data for SBOLS( ) begins. CALL DBOCLS( ) Use of this option adds 2 to the required length of IOPT(*). 6 - This option allows the user to provide an option array to the low-level subprogram SBOLSM( ). If LP is the processing pointer for IOPT(*), IOPT(LP)=6 IOPT(LP+1)= Position in IOPT(*) where option array data for SBOLSM( ) begins. CALL DBOCLS( ) Use of this option adds 2 to the required length of IOPT(*). 7 - Move the processing pointer (either forward or backward) to the location IOPT(LP+1). The processing pointer moves to locations LP+2 if option number 7 is used with the value -7. For example to skip over locations 3,...,NCOLS+2, IOPT(1)=7 IOPT(2)=NCOLS+3 (IOPT(I), I=3,...,NCOLS+2 are not defined here.) IOPT(NCOLS+3)=99 CALL DBOCLS( ) CAUTION: Misuse of this option can yield some very hard-to-find bugs. Use it with care. It is intended to be used for passing option arrays to other subprograms. 8 - This option allows the user to suppress the algorithmic feature of DBOCLS( ) that processes the constraint equations C*X = Y and resolves infeasibilities. The steps normally done are to solve C*X - Y = 0 in a least squares sense using the stated bounds on both X and Y. Then the "reachable" vector Y = C*X is computed using the solution X obtained. Finally the stated bounds for Y are enlarged to include C*X. To suppress the feature: IOPT(LP)=8 CALL DBOCLS( ) Use of this option adds 1 to the required length of IOPT(*). 9 - This option allows the user to suppress the pretriangularizing step of the least squares matrix that is done within DBOCLS( ). This is primarily a means of enhancing the subprogram efficiency and has little effect on accuracy. To suppress the step, set: IOPT(LP)=9 CALL DBOCLS( ) Use of this option adds 1 to the required length of IOPT(*). 99 -- There are no more options to change. Only option numbers -99, -9,-8,...,-1, 1,2,...,9, and 99 are permitted. Other values are errors. Options -99,-1,...,-9 mean that the respective options 99,1,...,9 are left at their default values. An example is the option to suppress the preprocessing of constraints: IOPT(1)=-8 Option is recognized but not changed IOPT(2)=99 CALL DBOCLS( ) Error Messages for DBOCLS() ----- -------- --- -------- WARNING in... DBOCLS(). THE ROW DIMENSION OF W(,)=(I1) MUST BE .GE. THE NUMBER OF EFFECTIVE ROWS=(I2). IN ABOVE MESSAGE, I1= 1 IN ABOVE MESSAGE, I2= 2 ERROR NUMBER = 41 WARNING IN... DBOCLS(). THE COLUMN DIMENSION OF W(,)=(I1) MUST BE .GE. NCOLS+ MCON+1=(I2). IN ABOVE MESSAGE, I1= 2 IN ABOVE MESSAGE, I2= 3 ERROR NUMBER = 42 WARNING IN... DBOCLS(). THE DIMENSIONS OF THE ARRAYS BL(),BU(), AND IND()=(I1) MUST BE .GE. NCOLS+MCON=(I2). IN ABOVE MESSAGE, I1= 1 IN ABOVE MESSAGE, I2= 2 ERROR NUMBER = 43 WARNING IN... DBOCLS(). THE DIMENSION OF X()=(I1) MUST BE IN ABOVE MESSAGE, I1= 1 IN ABOVE MESSAGE, I2= 2 ERROR NUMBER = 44 WARNING IN... DBOCLS(). THE . DBOCLS() THE DIMENSION OF IW()=(I1) MUST BE .GE. 2*NCOLS+2*MCON=(I2). IN ABOVE MESSAGE, I1= 1 IN ABOVE MESSAGE, I2= 4 ERROR NUMBER = 46 WARNING IN... DBOCLS(). THE DIMENSION OF IOPT()=(I1) MUST BE .GE. THE REQD. LEN.=(I2). IN ABOVE MESSAGE, I1= 16 IN ABOVE MESSAGE, I2= 18 ERROR NUMBER = 47 WARNING IN... DBOCLS(). ISCALE OPTION=(I1) MUST BE 1-3. IN ABOVE MESSAGE, I1= 0 ERROR NUMBER = 48 WARNING IN... DBOCLS(). OFFSET PAST X(NCOLS) (I1) FOR USER-PROVIDED COLUMN SCALING MUST BE POSITIVE. IN ABOVE MESSAGE, I1= 0 ERROR NUMBER = 49 WARNING IN... DBOCLS(). EACH PROVIDED COL. SCALE FACTOR MUST BE POSITIVE. COMPONENT (I1) NOW = (R1). IN ABOVE MESSAGE, I1= 1 IN ABOVE MESSAGE, R1= 0. ERROR NUMBER = 50 WARNING IN... DBOCLS(). THE OPTION NUMBER=(I1) IS NOT DEFINED. IN ABOVE MESSAGE, I1= 1001 ERROR NUMBER = 51 WARNING IN... DBOCLS(). NO. OF ROWS=(I1) MUST BE .GE. 0 .AND. .LE. MDW-MCON=(I2). IN ABOVE MESSAGE, I1= 2 IN ABOVE MESSAGE, I2= 1 ERROR NUMBER = 52 WARNING IN... DBOCLS(). MDW=(I1) MUST BE POSITIVE. IN ABOVE MESSAGE, I1= 0 ERROR NUMBER = 53 WARNING IN... DBOCLS(). MCON=(I1) MUST BE NONNEGATIVE. IN ABOVE MESSAGE, I1= -1 ERROR NUMBER = 54 WARNING IN... DBOCLS(). NCOLS=(I1) THE NO. OF VARIABLES MUST BE POSITIVE. IN ABOVE MESSAGE, I1= 0 ERROR NUMBER = 55 WARNING IN... DBOCLS(). FOR J=(I1), IND(J)=(I2) MUST BE 1-4. IN ABOVE MESSAGE, I1= 1 IN ABOVE MESSAGE, I2= 0 ERROR NUMBER = 56 WARNING IN... DBOCLS(). FOR J=(I1), BOUND BL(J)=(R1) IS .GT. BU(J)=(R2). IN ABOVE MESSAGE, I1= 1

IN ABOVE MESSAGE, R1= .1000000000E+01 IN ABOVE MESSAGE, R2= 0. ERROR NUMBER = 57 LINEAR CONSTRAINTS, SNLA REPT. SAND82-1517, AUG., (1982).

REFERENCES
R. J. Hanson, Linear least squares with bounds and linear constraints, Report SAND82-1517, Sandia Laboratories, August 1982.

ROUTINES CALLED
D1MACH, DASUM, DBOLS, DCOPY, DDOT, DNRM2, DSCAL, XERMSG

REVISION HISTORY (YYMMDD)
821220 DATE WRITTEN
891006 Cosmetic changes to prologue. (WRB)
891006 REVISION DATE from Version 3.2
891214 Prologue converted to Version 4.0 format. (BAB)
900510 Convert XERRWV calls to XERMSG calls. (RWC)
910819 Added variable M for MOUT+MCON in reference to DBOLS. (WRB)
920501 Reformatted the REFERENCES section. (WRB)