Documentation for module cp_lbfgs

LBFGS-B routine (version 3.0, April 25, 2011)

source: cp_lbfgs.F
Loading...

public Subroutines/Functions:

This subroutine partitions the working arrays wa and iwa, and then uses the limited memory BFGS method to solve the bound constrained optimization problem by calling mainlb. (The direct method will be used in the subspace minimization.)

SUBROUTINEsetulb(n, m, x, lower_bound, upper_bound, nbd, f, g, factr, pgtol, wa, iwa, task, iprint, csave, lsave, isave, dsave, trust_radius)

This subroutine partitions the working arrays wa and iwa, and then uses the limited memory BFGS method to solve the bound constrained optimization problem by calling mainlb. (The direct method will be used in the subspace minimization.)

Arguments:
INTEGER,
INTENT(in)
:: n n is the dimension of the problem.
INTEGER,
INTENT(in)
:: m m is the maximum number of variable metric corrections used to define the limited memory matrix.
REAL(dp),
INTENT(inout)
:: x(n) On entry x is an approximation to the solution. On exit x is the current approximation.
REAL(dp)
:: lower_bound(n) the lower bound on x.
REAL(dp)
:: upper_bound(n) the upper bound on x.
INTEGER
:: nbd(n) nbd represents the type of bounds imposed on the variables, and must be specified as follows: nbd(i)=0 if x(i) is unbounded, 1 if x(i) has only a lower bound, 2 if x(i) has both lower and upper bounds, and 3 if x(i) has only an upper bound.
REAL(dp)
:: f On first entry f is unspecified. On final exit f is the value of the function at x.
REAL(dp)
:: g(n) On first entry g is unspecified. On final exit g is the value of the gradient at x.
REAL(dp),
INTENT(in)
:: factr factr >= 0 is specified by the user. The iteration will stop when
REAL(dp),
INTENT(in)
:: pgtol pgtol >= 0 is specified by the user. The iteration will stop when
REAL(dp)
:: wa(2*m*n+5*n+11*m*m+8*m) working array
INTEGER
:: iwa(3*n) integer working array
CHARACTER(60)
:: task is a working string of characters of length 60 indicating the current job when entering and quitting this subroutine.
INTEGER
:: iprint iprint is a variable that must be set by the user. It controls the frequency and type of output generated: iprint<0 no output is generated; iprint=0 print only one line at the last iteration; 0100 print details of every iteration including x and g; When iprint > 0, the file iterate.dat will be created to summarize the iteration.
CHARACTER(60)
:: csave is a working string of characters
LOGICAL
:: lsave(4) lsave is a working array On exit with 'task' = NEW_X, the following information is available: If lsave(1) = .true. then the initial X has been replaced by its projection in the feasible set If lsave(2) = .true. then the problem is constrained; If lsave(3) = .true. then each variable has upper and lower bounds;
INTEGER
:: isave(44) isave is a working array On exit with 'task' = NEW_X, the following information is available: isave(22) = the total number of intervals explored in the search of Cauchy points; isave(26) = the total number of skipped BFGS updates before the current iteration; isave(30) = the number of current iteration; isave(31) = the total number of BFGS updates prior the current iteration; isave(33) = the number of intervals explored in the search of Cauchy point in the current iteration; isave(34) = the total number of function and gradient evaluations; isave(36) = the number of function value or gradient evaluations in the current iteration; if isave(37) = 0 then the subspace argmin is within the box; if isave(37) = 1 then the subspace argmin is beyond the box; isave(38) = the number of free variables in the current iteration; isave(39) = the number of active constraints in the current iteration; n + 1 - isave(40) = the number of variables leaving the set of active constraints in the current iteration; isave(41) = the number of variables entering the set of active constraints in the current iteration.
REAL(dp)
:: dsave(29) dsave is a working array of dimension 29. On exit with 'task' = NEW_X, the following information is available: dsave(1) = current 'theta' in the BFGS matrix; dsave(2) = f(x) in the previous iteration; dsave(3) = factr*epsmch; dsave(4) = 2-norm of the line search direction vector; dsave(5) = the machine precision epsmch generated by the code; dsave(7) = the accumulated time spent on searching for Cauchy points; dsave(8) = the accumulated time spent on subspace minimization; dsave(9) = the accumulated time spent on line search; dsave(11) = the slope of the line search function at the current point of line search; dsave(12) = the maximum relative step length imposed in line search; dsave(13) = the infinity norm of the projected gradient; dsave(14) = the relative step length in the line search; dsave(15) = the slope of the line search function at the starting point of the line search; dsave(16) = the square of the 2-norm of the line search direction vector.
REAL(dp),
INTENT(in)
:: trust_radius ...