MODULE norm_min !*..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..* !| .**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**. | !| | | | !| | THE SOLVER FOR THE NORM MINIMIZATION PROBLEM | | !| | | | !| .**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**. | !*..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..* !| | !| Utilizes new version of PLQDF1 by Ladislav Luksan as a quadratic solver. | !| | !| Utilizes PVMM by Ladislav Luksan as a norm minimization solver. This subroutine | !| uses PQSUBS and MQSUBS by Ladislav Luksan. PVMM is a VARIABLE METRIC ALGORITHM | !| for UNCONSTRAINED and LINEARLY CONSTRAINED OPTIMIZATION. | !| | !| The subroutine PVMM together with PQSUBS and MQSUBS are licensed by | !| the GNU Lesser General Public License (LGPL). | !| | !| | !| | !| .**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**. | !*..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..* USE constants ! Double precision (i.e. accuracy) and varibles: NRES,NDEC,NREM,NADD,NIT,NFV,NFG,NFH USE bundle1 ! The BUNDLE of the DC component f_1 USE bundle2 ! The BUNDLE of the DC component f_2 USE functions, ONLY : user_size_b1, user_size_b2, user_n ! Contain INFORMATION from the USER: ! The biggest possible size of both bundles and also the number of variables. IMPLICIT NONE EXTERNAL PVMM ! Contains the Variable metric method PVVML by Ladislav Luksan and it is needed to solve ! the norm minimization problem at Step 3 of the 'main iteration' algorithm. EXTERNAL PQSUBS ! Basic modules for PVVM (by Ladislav Luksan) EXTERNAL MQSUBS ! Matrix modules for PVVM (by Ladislav Luksan) INTEGER, PARAMETER :: m_size = user_size_b1 + user_size_b2 ! The biggest possible size of the matrix (and the vector) used to define the subgradient ! (and the objective function) of the norm minimization problem. (NOTICE: aggregated elements of B_1 and B_2 are NOT used) INTEGER, PARAMETER :: lkm = user_n ! The parameter 'lkm' is the number of variables in the original minimization problem ! and it is needed/used in the EXTERNAL subroutine FUNDER. ! Matrices which are used to define the objective function of the norm minimization problem and its gradient. ! Both needed/used in the EXTERNAL subroutine FUNDER. REAL(KIND=dp), DIMENSION(lkm, m_size), SAVE :: f_matrix ! used to define the objective function in the norm minimization problem. REAL(KIND=dp), DIMENSION(m_size, m_size), SAVE :: g_matrix ! used to define gradient of the objective function in the norm minimization problem. CONTAINS !*..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..* !| .**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**. | !| | | | !| | CONTAINS SUBROUTINES: | | !| | | | !| | SOLVER FOR NORM PROBLEM : norm_solver(bxi, obj, NF, n1, n2, B1, B2) | | !| | | | !| | FORMS MATRICES USED IN NORM PROBLEM : norm_matrix(dimensio, n1, n2, B1, B2) | | !| | | | !| | | | !| .**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**. | !*..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..**..* !*************************************************************************************** ! ---------------------------------------------------------------------------------- | ! | | | ! | SOLVER FOR THE NORM MINIMIZATION PROBLEM | | ! | | | ! ---------------------------------------------------------------------------------- | !*************************************************************************************** SUBROUTINE norm_solver(bxi, obj, NF, n1, n2, B1, B2) ! ! Solves the norm minimization problem at Step 3 of the 'main iteration' algorithm. ! ! INPUT: * 'B1' and 'B2': The bundles B_1 and B_2 of the DC component f_1 and f_2 ! * 'n1' and 'n2': The size of the bundles B_1 and B_2 (also the 'current bundle element' is (i.e. has to be) taken into account in both bundles) ! * 'NF' : The number of variables in the objective function of the norm minimization problem (NOTICE: NF = n1 + n2) ! ! OUTPUT: * 'bxi': The solution vector bxi = bxi_1^* - bxi_2^* = f_matrix * lambda ('lambda' is approximation to the minimum) ! * 'obj': The value of the objective funtion at the solution 'lambda' (i.e. '|| bxi ||') ! !*********************************************************************************** IMPLICIT NONE !**************************** NEEDED FROM USER ************************************* TYPE(kimppu1), INTENT(IN) :: B1 ! the bundle B_1 for the DC component f_1 TYPE(kimppu2), INTENT(IN) :: B2 ! the bundle B_2 for the DC component f_2 REAL(KIND=dp), DIMENSION(give_n_b1(B1)), INTENT(OUT) :: bxi ! Output: bxi = bxi_1^* - bxi_2^* = f_matrix * lambda REAL(KIND=dp), INTENT(OUT) :: obj ! Output: the value of the objective function at the solution 'lambda' INTEGER, INTENT(IN) :: NF ! the number of variables in objective function of the norm min problem (NOTICE: NF = n1 + n2) INTEGER, INTENT(IN) :: n1, n2 ! the size of the bundle B_1 and B_2 (also the 'current element' is (i.e. has to be) taken into account in both bundles) !***************************** LOCAL VARIABLES WHEN CALLING PVMML ****************** ! .. Help variables .. INTEGER :: dimensio ! the length of subgradients INTEGER :: i,j, ind ! .. Other variables .. INTEGER, PARAMETER :: NB = 3 ! the simple bounds are accepted because NB > 0 INTEGER, PARAMETER :: NC = 2 ! the number of linear constraints (2 equality constraints) INTEGER :: IPRNT ! Specifies print in PVMML: IPRNT = 0 print is suppressed ! IPRNT = 1 basic print of final results ! IPRNT = -1 extended print of final results ! IPRNT = 2 basic print of intermediate and final results ! IPRNT = -2 extended print of intermediate and final results INTEGER :: ITERM ! OUTPUT from PVMML: indicates the cause of termination REAL(KIND=dp) :: GMAX, F ! OUTPUT from PVMML : F - the value of the objective at the solution ! GMAX - maximum absolute value of a partial derivative of the objective function INTEGER, DIMENSION(7) :: IPAR ! INTEGER parameters (If IPAR = 0, then default values are used in PVMML) REAL(KIND=dp), DIMENSION(NF) :: lambda ! ON INPUT: a vector with the initial estimate to the solution ! ON OUTPUT: the approximation to the minimum REAL(KIND=dp), DIMENSION(7) :: RPAR ! REAL parameters (IF RPAR = 0.0_dp, then default values are used in PVMML) ! .. Array Arguments in PVMML.. INTEGER, DIMENSION(NF) :: IX ! the vector containing types of simple bounds (IX=1 because XL(I) <= X(I) for each variable) REAL(KIND=dp), DIMENSION(NF) :: XL ! the vector containing lower bounds (XL(I) = 0.0_dp for each variable) REAL(KIND=dp), DIMENSION(NF) :: XU ! vector containing upper bounds (This is not needed) REAL(KIND=dp), DIMENSION(NC) :: CF ! optional/auxiliary: the vector containing values of constraint functions INTEGER, DIMENSION(NC) :: IC ! the vector containing constraint types: 5 - the equality constraint (NOTICE: CL(I) = C(I) = CU(I) ) REAL(KIND=dp), DIMENSION(NC) :: CL ! the vector containing lower bounds for constraints REAL(KIND=dp), DIMENSION(NC) :: CU ! the vector containing upper bounds for constraints REAL(KIND=dp), DIMENSION(NF*NC) :: CG ! the vector whose columns are normals of linear contsraints ! .. Initialization .. IPAR = 0 ! default values are used RPAR = 0.0_dp ! default values are used IPRNT = 0 ! print (because IPRNT=0, print is suppressed) IX = 1 ! the simple bounds are lower bounds XL = 0.0_dp ! values of the lower bounds XU = 10.0_dp ! values of the upper bounds (NOT needed here!) !CF = (/ /) ! optional IC = (/ 5, 5 /) ! constraints are equality constraints CL = (/ 1.0_dp, 1.0_dp /) ! lower bounds for constraints CU = (/ 1.0_dp, 1.0_dp /) ! upper bounds for constraints CG = 0.0_dp ! Initialization of CG DO i = 1, n1 ! the first constraint CG(i) = 1.0_dp END DO ind = 2*n1 + n2 + 1 DO i = ind, (ind + n2) ! the second constraint CG(i) = 1.0_dp END DO ! initial estimate to the solution lambda = 0.0_dp ! initialization of lambda lambda(1) = 1.0_dp lambda(n1+1) = 1.0_dp ! the length of subgradients (same as the number of variables in the original minimization problem). dimensio = give_n_b1(B1) CALL norm_matrix(dimensio, n1, n2, B1, B2) ! Calculates the matrices used to define the objective function and ! its gradient in the norm minimization problem ! Execution of the solver PVMML by Ladislav Luksan CALL PVMML(NF,NB,NC,lambda,IX,XL,XU,CF,IC,CL,CU,CG,IPAR,& ! objective F = || f_matrix * lambda ||^2 & RPAR, F, GMAX, IPRNT, ITERM ) IF ( F < 0.0_dp) THEN obj = 0.0_dp ELSE obj = SQRT(F) ! the original objective function 'obj' = || f_matrix * lambda || END IF bxi = 0.0_dp ! initialization DO j = 1, NF DO i = 1, dimensio bxi(i) = bxi(i) + f_matrix(i,j)*lambda(j) END DO END DO END SUBROUTINE norm_solver !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - !*************************************************************************************** ! ---------------------------------------------------------------------------------- | ! | | | ! | CALCULATE MATRICES f_matrix and g_matrix USED IN FUNDER | | ! | | | ! ---------------------------------------------------------------------------------- | !*************************************************************************************** SUBROUTINE norm_matrix(dimensio, n1, n2, B1, B2) ! Calculates the matrices used to define the objective of norm minimization problem ! and its gradient. These matrices are saved to the variables 'f_matrix' and 'g_matrix' ! ! NOTICE: * 'dimensio' : the length of the subgradient ! * 'n1' and 'n2' : the size of the bundles B_1 and B_2 (the 'current bundle element' is also taken into account in both bundles) (aggregation is NOT used) ! * 'B1' and 'B2* : the bundles B_1 and B_2 of the DC components f_1 and f_2 ! INTEGER, INTENT(IN) :: dimensio, n1, n2 ! dimensio - the length of the subgradients ! n1 and n2 - the size of the bundles B_1 and B_2 (the 'current bundle element' is also taken into account in both bundles) INTEGER :: i, j, k, ind ! .. Help variables .. REAL(KIND=dp), DIMENSION(dimensio) :: grad_b1, grad_b2 REAL(KIND=dp), DIMENSION(n1+n2,n1+n2) :: M ! M = grad^T * grad REAL(KIND=dp), DIMENSION(dimensio, n1+n2) :: grad ! the matrix whose columns are subgradients ! contains all the subgradients (both from B_1 and B_2) (aggregated element is NOT included) ! first there are the ones from B_1 and then the ones from B_2 ! .. Bundles .. TYPE(kimppu1), INTENT(IN) :: B1 ! the bundle of the DC component f_1 TYPE(kimppu2), INTENT(IN) :: B2 ! the bundle of the DC component f_2 ! .. Matrix 'grad' is formed .. DO i = 0, give_size_b1(B1) ! all subgradients from the bundle B_1 are looked through (except for the aggregated element) grad_b1 = give_subgrad_b1(B1,i) ! the subgradient of the bundle element i from the bundle B_1 DO j = 1, dimensio ! inserts the current grad_b1 into the matrix 'grad' (the starting position is (i+1):th column) grad(j, i+1) = grad_b1(j) END DO END DO ind = give_size_b1(B1) + 2 ! the starting position for subgradients from the bundle B_2 DO i = 0, give_size_b2(B2) ! all subgradients from the bundle B_2 are looked through (except for the aggregated element) grad_b2 = give_subgrad_b2(B2,i) ! the subgradient of the bundle element i from the bundle B_2 DO j = 1, dimensio ! inserts the current grad_b2 into the matrix 'grad' (NOTICE: the negative sign) grad(j, ind + i) = - grad_b2(j) ! (the position is (ind+i):th column in the matrix) END DO END DO ! .. Matrix 'grad' is saved to 'f_matrix' .. DO j = 1, n1+n2 DO i = 1, dimensio f_matrix(i,j) = grad(i,j) END DO END DO M = MATMUL( TRANSPOSE(grad), grad ) ! M = grad^T * grad ! .. Matrix 'M' is saved to 'g_matrix' .. DO j = 1, n1+n2 DO i = 1, n1+n2 g_matrix(i,j) = M(i,j) END DO END DO END SUBROUTINE norm_matrix !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - END MODULE norm_min