BCNLS
Solves a nonlinear least-squares problem subject to bounds on the variables and general linear constraints.
Required Arguments
FCN — User-supplied subroutine to evaluate the function to be minimized. The usage is CALL FCN (M, N, X, F), where
M ‑ Number of functions. (Input)
N ‑ Number of variables. (Input)
X ‑ Array of length N containing the point at which the function will be evaluated. (Input)
F ‑ Array of length M containing the computed function at the point X. (Output)
The routine FCN must be declared EXTERNAL in the calling program.
M — Number of functions. (Input)
C — MCON × N matrix containing the coefficients of the MCON general linear constraints. (Input)
BL — Vector of length MCON containing the lower limit of the general constraints. (Input).
BU — Vector of length MCON containing the upper limit of the general constraints. (Input).
IRTYPE — Vector of length MCON indicating the types of general constraints in the matrix C. (Input)
Let R(I) = C(I, 1)*X(1) + … + C(I, N)*X(N). Then the value of IRTYPE(I) signifies the following:
Let R(I) = C(I, 1)*X(1) + … + C(I, N)*X(N). Then the value of IRTYPE(I) signifies the following:
IRTYPE(I) | I-th CONSTRAINT |
|---|---|
0 | BL(I).EQ.R(I).EQ.BU(I) |
1 | R(I).LE.BU(I) |
2 | R(I).GE.BL(I) |
3 | BL(I).LE.R(I).LE.BU(I) |
XLB — Vector of length N containing the lower bounds on variables; if there is no lower bound on a variable, then 1.0E30 should be set as the lower bound. (Input)
XUB — Vector of length N containing the upper bounds on variables; if there is no upper bound on a variable, then ‑1.0E30 should be set as the upper bound. (Input)
X — Vector of length N containing the approximate solution. (Output)
Optional Arguments
N — Number of variables. (Input)
Default: N = SIZE (C,2).
Default: N = SIZE (C,2).
MCON — The number of general linear constraints for the system, not including simple bounds. (Input)
Default: MCON = SIZE (C,1).
Default: MCON = SIZE (C,1).
LDC — Leading dimension of C exactly as specified in the dimension statement of the calling program. (Input)
LDC must be at least MCON.
LDC must be at least MCON.
Default: LDC = SIZE (C,1).
XGUESS — Vector of length N containing the initial guess. (Input)
Default: XGUESS = 0.0.
Default: XGUESS = 0.0.
RNORM — The Euclidean length of components of the function f (x) after the approximate solution has been found. (Output).
ISTAT — Scalar indicating further information about the approximate solution X. (Output)
See the Comments section for a description of the tolerances and the vectors IPARAM and RPARAM.
See the Comments section for a description of the tolerances and the vectors IPARAM and RPARAM.
ISTAT | Meaning |
|---|---|
1 | The function f(x) has a length less than TOLF = RPARAM(1). This is the expected value for ISTAT when an actual zero value of f (x) is anticipated. |
2 | The function f(x) has reached a local minimum. This is the expected value for ISTAT when a nonzero value of f(x) is anticipated. |
3 | A small change (absolute) was noted for the vector x. A full model problem step was taken. The condition for ISTAT = 2 may also be satisfied, so that a minimum has been found. However, this test is made before the test for ISTAT = 2. |
4 | A small change (relative) was noted for the vector x. A full model problem step was taken. The condition for ISTAT = 2 may also be satisfied, so that a minimum has been found. However, this test is made before the test for ISTAT = 2. |
5 | The number of terms in the quadratic model is being restricted by the amount of storage allowed for that purpose. It is suggested, but not required, that additional storage be given for the quadratic model parameters. This is accessed through the vector IPARAM, documented below. |
6 | Return for evaluation of function and Jacobian if reverse communication is desired. See the Comments below. |
FORTRAN 90 Interface
Generic: CALL BCNLS (FCN, M, C, BL, BU, IRTYPE, XLB, XUB, X [, …])
Specific: The specific interface names are S_BCNLS and D_BCNLS.
FORTRAN 77 Interface
Single: CALL BCNLS (FCN, M, N, MCON, C, LDC, BL, BU, IRTYPE, XLB, XUB, XGUESS, X, RNORM, ISTAT)
Double: The double precision name is DBCNLS.
Description
The routine BCNLS solves the nonlinear least squares problem
subject to
BCNLS is based on the routine DQED by R.J. Hanson and F.T. Krogh. The section of BCNLS that approximates, using finite differences, the Jacobian of f(x) is a modification of JACBF by D.E. Salane.
Comments
1. Workspace may be explicitly provided, if desired, by use of B2NLS/DB2NLS. The reference is:
CALL B2NLS (FCN, M, N, MCON, C, LDC, BL, BU, IRTYPE, XLB, XUB, XGUESS, X, RNORM,ISTAT, IPARAM, RPARAM, JAC, F, FJ, LDFJ, IWORK, LIWORK, WORK, LWORK)
The additional arguments are as follows:
IPARAM — Integer vector of length six used to change certain default attributes of BCNLS. (Input).
If the default parameters are desired for B2NLS, set IPARAM(1) to zero. Otherwise, if any nondefault parameters are desired for IPARAM or RPARAM, the following steps should be taken before calling B2NLS:
If the default parameters are desired for B2NLS, set IPARAM(1) to zero. Otherwise, if any nondefault parameters are desired for IPARAM or RPARAM, the following steps should be taken before calling B2NLS:
CALL B7NLS (IPARAM, RPARAM)
Set nondefault values for IPARAM and RPARAM.
If double precision is being used, DB7NLS should be called instead. Following is a list of parameters and the default values.
If double precision is being used, DB7NLS should be called instead. Following is a list of parameters and the default values.
IPARAM(1) = Initialization flag.
IPARAM(2) = ITMAX, the maximum number of iterations allowed.
Default: 75
Default: 75
IPARAM(3) = a flag that suppresses the use of the quadratic model in the inner loop. If set to one, then the quadratic model is never used. Otherwise use the quadratic model where appropriate. This option decreases the amount of workspace as well as the computing overhead required. A user may wish to determine if the application really requires the use of the quadratic model.
Default: 0
Default: 0
IPARAM(4) = NTERMS, one more than the maximum number of terms used in the quadratic model.
Default: 5
Default: 5
IPARAM(5) = RCSTAT, a flag that determines whether forward or reverse communication is used. If set to zero, forward communication through functions FCN and JAC is used. If set to one, reverse communication is used, and the dummy routines B10LS/DB10LS and B11LS/DB11LS may be used in place of FCN and JAC, respectively. When BCNLS returns with ISTAT = 6, arrays F and FJ are filled with f(x) and the Jacobian of f(x), respectively. BCNLS is then called again.
Default: 0
IPARAM(6) = a flag that determines whether the analytic Jacobian, as supplied in JAC, is used, or if a finite difference approximation is computed. If set to zero, JAC is not accessed and finite differences are used. If set to one, JAC is used to compute the Jacobian.
Default: 0
Default: 0
RPARAM — Real vector of length 7 used to change certain default attributes of
BCNLS. (Input)
For the description of RPARAM, we make the following definitions:
FC current value of the length of f (x)
FB best value of length of f (x)
FL value of length of f (x) at the previous step
PV predicted value of length of f (x), after the step is taken, using the approximating model ɛ machine epsilon = amach(4).
The conditions ∣FB ‑ PV∣ ≤ TOLSNR*FB and ∣FC ‑ PV∣ ≤ TOLP*FB and ∣FC ‑ FL∣ ≤ TOLSNR*FB together with taking a full model step, must be satisfied before the condition ISTAT = 2 is returned. (Decreasing any of the values for TOLF, TOLD, TOLX, TOLSNR, or TOLP will likely increase the number of iterations required for convergence.)
BCNLS. (Input)
For the description of RPARAM, we make the following definitions:
FC current value of the length of f (x)
FB best value of length of f (x)
FL value of length of f (x) at the previous step
PV predicted value of length of f (x), after the step is taken, using the approximating model ɛ machine epsilon = amach(4).
The conditions ∣FB ‑ PV∣ ≤ TOLSNR*FB and ∣FC ‑ PV∣ ≤ TOLP*FB and ∣FC ‑ FL∣ ≤ TOLSNR*FB together with taking a full model step, must be satisfied before the condition ISTAT = 2 is returned. (Decreasing any of the values for TOLF, TOLD, TOLX, TOLSNR, or TOLP will likely increase the number of iterations required for convergence.)
RPARAM(1) = TOLF, tolerance used for stopping when FC ≤ TOLF.
RPARAM(2) = TOLX, tolerance for stopping when change to x values has length less than or equal to TOLX*length of x values.
RPARAM(3) = TOLD, tolerance for stopping when change to x values has length less than or equal to TOLD.
RPARAM(4) = TOLSNR, tolerance used in stopping condition ISTAT = 2.
Default: 1.E5
Default: 1.E5
RPARAM(5) = TOLP, tolerance used in stopping condition ISTAT = 2.
Default: 1.E5
Default: 1.E5
RPARAM(6) = TOLUSE, tolerance used to avoid values of x in the quadratic model's interpolation of previous points. Decreasing this value may result in more terms being included in the quadratic model.
RPARAM(7) = COND, largest condition number to allow when solving for the quadratic model coefficients. Increasing this value may result in more terms being included in the quadratic model.
Default: 30
Default: 30
JAC — User-supplied subroutine to evaluate the Jacobian. The usage is CALL JAC(M, N, X, FJAC, LDFJAC), where
M ‑ Number of functions. (Input)
N ‑ Number of variables. (Input)
X ‑ Array of length N containing the point at which the Jacobian will be evaluated. (Input)
FJAC ‑ The computed M × N Jacobian at the point X. (Output)
LDFJAC ‑ Leading dimension of the array FJAC. (Input)
The routine JAC must be declared EXTERNAL in the calling program.
M ‑ Number of functions. (Input)
N ‑ Number of variables. (Input)
X ‑ Array of length N containing the point at which the Jacobian will be evaluated. (Input)
FJAC ‑ The computed M × N Jacobian at the point X. (Output)
LDFJAC ‑ Leading dimension of the array FJAC. (Input)
The routine JAC must be declared EXTERNAL in the calling program.
F — Real vector of length N used to pass f(x) if reverse communication (IPARAM(4)) is enabled. This array must be allocated regardless of the setting of (IPARAM(4)). (Input)
FJ — Real array of size M × N. It is used to store the Jacobian matrix of f(x) if reverse communication (IPARAM(4)) is enabled. This array must be allocated regardless of the setting of (IPARAM(4)). (Input)
Specifically,
Specifically,
LDFJ — Leading dimension of FJ exactly as specified in the dimension statement of the calling program. (Input)
IWORK — Integer work vector of length LIWORK.
LIWORK — Length of work vector IWORK. LIWORK must be at least
5MCON + 12N + 47 + MAX(M, N)
5MCON + 12N + 47 + MAX(M, N)
WORK — Real work vector of length LWORK
LWORK — Length of work vector WORK. LWORK must be at least 41N + 6M + 11MCON + (M + MCON)(N + 1) + NA(NA + 7) + 8 MAX(M, N) + 99. Where NA = MCON + 2N + 6.
2. Informational errors
Type | Code | Description |
|---|---|---|
3 | 1 | The function f(x) has reached a value that may be a local minimum. However, the bounds on the trust region defining the size of the step are being hit at each step. Thus, the situation is suspect. (Situations of this type can occur when the solution is at infinity at some of the components of the unknowns, x). |
3 | 2 | The model problem solver has noted a value for the linear or quadratic model problem residual vector length that is greater than or equal to the current value of the function, i.e. the Euclidean length of f(x). This situation probably means that the evaluation of f (x) has more uncertainty or noise than is possible to account for in the tolerances used to not a local minimum. The value of x is suspect, but a minimum has probably been found. |
3 | 3 | More than ITMAX iterations were taken to obtain the solution. The value obtained for x is suspect, although it is the best set of x values that occurred in the entire computation. The value of ITMAX can be increased though the IPARAM vector. |
Examples
Example 1
This example finds the four variables x1, x2, x3, x4 that are in the model function
There are values of h(t) at five values of t.
h(0.05) = 2.206
h(0.1) = 1.994
h(0.4) = 1.35
h(0.5) = 1.216
h(1.0) = 0.7358
There are also the constraints that x2, x4 ≤ 0, x1, x3 ≥ 0, and x2 and x4 must be separated by at least 0.05. Nothing more about the values of the parameters is known so the initial guess is 0.
USE BCNLS_INT
USE UMACH_INT
USE WRRRN_INT
IMPLICIT NONE
INTEGER MCON, N
PARAMETER (MCON=1, N=4)
! SPECIFICATIONS FOR PARAMETERS
INTEGER LDC, M
PARAMETER (M=5, LDC=MCON)
! SPECIFICATIONS FOR LOCAL VARIABLES
INTEGER IRTYPE(MCON), NOUT
REAL BL(MCON), C(MCON,N), RNORM, X(N), XLB(N), &
XUB(N)
! SPECIFICATIONS FOR SUBROUTINES
! SPECIFICATIONS FOR FUNCTIONS
EXTERNAL FCN
!
CALL UMACH (2, NOUT)
! Define the separation between x(2)
! and x(4)
C(1,1) = 0.0
C(1,2) = 1.0
C(1,3) = 0.0
C(1,4) = -1.0
BL(1) = 0.05
IRTYPE(1) = 2
! Set lower bounds on variables
XLB(1) = 0.0
XLB(2) = 1.0E30
XLB(3) = 0.0
XLB(4) = 1.0E30
! Set upper bounds on variables
XUB(1) = -1.0E30
XUB(2) = 0.0
XUB(3) = -1.0E30
XUB(4) = 0.0
!
CALL BCNLS (FCN, M, C, BL, BL, IRTYPE, XLB, XUB, X, RNORM=RNORM)
CALL WRRRN ('X', X, 1, N, 1)
WRITE (NOUT,99999) RNORM
99999 FORMAT (/, 'rnorm = ', E10.5)
END
!
SUBROUTINE FCN (M, N, X, F)
! SPECIFICATIONS FOR ARGUMENTS
INTEGER M, N
REAL X(*), F(*)
! SPECIFICATIONS FOR LOCAL VARIABLES
INTEGER I
! SPECIFICATIONS FOR SAVE VARIABLES
REAL H(5), T(5)
SAVE H, T
! SPECIFICATIONS FOR INTRINSICS
INTRINSIC EXP
REAL EXP
!
DATA T/0.05, 0.1, 0.4, 0.5, 1.0/
DATA H/2.206, 1.994, 1.35, 1.216, 0.7358/
!
DO 10 I=1, M
F(I) = X(1)*EXP(X(2)*T(I)) + X(3)*EXP(X(4)*T(I)) - H(I)
10 CONTINUE
RETURN
END
Output
X
1 2 3 4
1.999 -1.000 0.500 -9.954
rnorm = .42425E-03
Example 2
This example solves the same problem as the last example, but reverse communication is used to evaluate f(x) and the Jacobian of f(x). The use of the quadratic model is turned off.
USE B2NLS_INT
USE UMACH_INT
USE WRRRN_INT
IMPLICIT NONE
INTEGER LDC, LDFJ, M, MCON, N
PARAMETER (M=5, MCON=1, N=4, LDC=MCON, LDFJ=M)
! Specifications for local variables
INTEGER I, IPARAM(6), IRTYPE(MCON), ISTAT, IWORK(1000), &
LIWORK, LWORK, NOUT
REAL BL(MCON), C(MCON,N), F(M), FJ(M,N), RNORM, RPARAM(7), &
WORK(1000), X(N), XGUESS(N), XLB(N), XUB(N)
REAL H(5), T(5)
SAVE H, T
INTRINSIC EXP
REAL EXP
! Specifications for subroutines
EXTERNAL B7NLS
! Specifications for functions
EXTERNAL B10LS, B11LS
!
DATA T/0.05, 0.1, 0.4, 0.5, 1.0/
DATA H/2.206, 1.994, 1.35, 1.216, 0.7358/
!
CALL UMACH (2, NOUT)
! Define the separation between x(2)
! and x(4)
C(1,1) = 0.0
C(1,2) = 1.0
C(1,3) = 0.0
C(1,4) = -1.0
BL(1) = 0.05
IRTYPE(1) = 2
! Set lower bounds on variables
XLB(1) = 0.0
XLB(2) = 1.0E30
XLB(3) = 0.0
XLB(4) = 1.0E30
! Set upper bounds on variables
XUB(1) = -1.0E30
XUB(2) = 0.0
XUB(3) = -1.0E30
XUB(4) = 0.0
! Set initial guess to 0.0
XGUESS = 0.0E0
! Call B7NLS to set default parameters
CALL B7NLS (IPARAM, RPARAM)
! Suppress the use of the quadratic
! model, evaluate functions and
! Jacobian by reverse communication
IPARAM(3) = 1
IPARAM(5) = 1
IPARAM(6) = 1
LWORK = 1000
LIWORK = 1000
! Specify dummy routines for FCN
! and JAC since we are using reverse
! communication
10 CONTINUE
CALL B2NLS (B10LS, M, N, MCON, C, LDC, BL, BL, IRTYPE, XLB, &
XUB, XGUESS, X, RNORM, ISTAT, IPARAM, RPARAM, &
B11LS, F, FJ, LDFJ, IWORK, LIWORK, WORK, LWORK)
!
! Evaluate functions if the routine
! returns with ISTAT = 6
IF (ISTAT .EQ. 6) THEN
DO 20 I=1, M
FJ(I,1) = EXP(X(2)*T(I))
FJ(I,2) = T(I)*X(1)*FJ(I,1)
FJ(I,3) = EXP(X(4)*T(I))
FJ(I,4) = T(I)*X(3)*FJ(I,3)
F(I) = X(1)*FJ(I,1) + X(3)*FJ(I,3) - H(I)
20 CONTINUE
GO TO 10
END IF
!
CALL WRRRN ('X', X, 1, N, 1)
WRITE (NOUT,99999) RNORM
99999 FORMAT (/, 'rnorm = ', E10.5)
END
Output
X
1 2 3 4
1.999 -1.000 0.500 -9.954
rnorm = .42413E-03
Published date: 03/19/2020
Last modified date: 03/19/2020
