Statistiques
| Révision :

root / src / Hupdate_all.f90 @ 12

Historique | Voir | Annoter | Télécharger (5,59 ko)

1 1 pfleura2
 subroutine hupdate_all(n, ds, dgrad, Hess)
2 1 pfleura2
3 12 pfleura2
!----------------------------------------------------------------------
4 12 pfleura2
!  Copyright 2003-2014 Ecole Normale Supérieure de Lyon,
5 12 pfleura2
!  Centre National de la Recherche Scientifique,
6 12 pfleura2
!  Université Claude Bernard Lyon 1. All rights reserved.
7 12 pfleura2
!
8 12 pfleura2
!  This work is registered with the Agency for the Protection of Programs
9 12 pfleura2
!  as IDDN.FR.001.100009.000.S.P.2014.000.30625
10 12 pfleura2
!
11 12 pfleura2
!  Authors: P. Fleurat-Lessard, P. Dayal
12 12 pfleura2
!  Contact: optnpath@gmail.com
13 12 pfleura2
!
14 12 pfleura2
! This file is part of "Opt'n Path".
15 12 pfleura2
!
16 12 pfleura2
!  "Opt'n Path" is free software: you can redistribute it and/or modify
17 12 pfleura2
!  it under the terms of the GNU Affero General Public License as
18 12 pfleura2
!  published by the Free Software Foundation, either version 3 of the License,
19 12 pfleura2
!  or (at your option) any later version.
20 12 pfleura2
!
21 12 pfleura2
!  "Opt'n Path" is distributed in the hope that it will be useful,
22 12 pfleura2
!  but WITHOUT ANY WARRANTY; without even the implied warranty of
23 12 pfleura2
!
24 12 pfleura2
!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
25 12 pfleura2
!  GNU Affero General Public License for more details.
26 12 pfleura2
!
27 12 pfleura2
!  You should have received a copy of the GNU Affero General Public License
28 12 pfleura2
!  along with "Opt'n Path". If not, see <http://www.gnu.org/licenses/>.
29 12 pfleura2
!
30 12 pfleura2
! Contact The Office of Technology Licensing, valorisation@ens-lyon.fr,
31 12 pfleura2
! for commercial licensing opportunities.
32 12 pfleura2
!----------------------------------------------------------------------
33 12 pfleura2
34 1 pfleura2
  use Path_module, only :  HUpdate,Hinv
35 1 pfleura2
  use Io_module, only : IoOut
36 1 pfleura2
37 1 pfleura2
  IMPLICIT NONE
38 1 pfleura2
39 1 pfleura2
  INTEGER, PARAMETER :: KINT=KIND(1)
40 1 pfleura2
  INTEGER, PARAMETER :: KREAL=KIND(1.0D0)
41 1 pfleura2
42 1 pfleura2
  INTEGER(KINT), INTENT(IN) :: n
43 1 pfleura2
  real(KREAL), INTENT(IN) :: ds(n), dgrad(n)
44 1 pfleura2
  REAL(KREAL),INTENT(INOUT) :: Hess(n:n)
45 1 pfleura2
46 1 pfleura2
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
47 1 pfleura2
!
48 1 pfleura2
! This subroutine is a driver for the update of the Hessian
49 1 pfleura2
! in the Quasi-Newton optimization step.
50 1 pfleura2
! We propose mostly update of the inverse Hessian.
51 1 pfleura2
!
52 1 pfleura2
! input:
53 1 pfleura2
!  - n: number of coordinates
54 1 pfleura2
!  - ds: GeomOld-Geom (expressed in the optimization coord.)
55 1 pfleura2
!  - dgrad: GradOld-Grad  (expressed in the optimization coord.)
56 1 pfleura2
!
57 1 pfleura2
!  input/output:
58 1 pfleura2
!  - Hess: old Hessian in input, replaced by updated Hessian on output
59 1 pfleura2
!
60 1 pfleura2
!!!!!!!!!!!!!!!!
61 1 pfleura2
!
62 1 pfleura2
! For now, we propose :
63 1 pfleura2
! - For update of inverse Hessian:
64 1 pfleura2
!      - BFGS (Broyden, Fletcher, Goldfarb, and Shanno)
65 1 pfleura2
!      - DFP (Davidson - Fletcher - Powell )
66 1 pfleura2
!      - MS (Murtagh-Sargent)
67 1 pfleura2
!
68 1 pfleura2
! - For update of hessian:
69 1 pfleura2
!      - Bofill (combination of Murtagh-Sargent and Powell-symmetric-Broyden
70 1 pfleura2
!      - MS (Murtagh-Sargent)
71 1 pfleura2
!      - BFGS
72 1 pfleura2
!      - DFP
73 1 pfleura2
!----------------------------------------------------
74 1 pfleura2
! We use the fact that some formula are the dual of other
75 1 pfleura2
! that is, one can convert the formula for update of the Hessian
76 1 pfleura2
! into the update of the inverse Hessian by replacing the Hessian
77 1 pfleura2
! (denoted by B in the Fletcher book and in the Nocedal&Wrigth book)
78 1 pfleura2
! by the inverse hessian (denoted by H), and replacing ds by dgrad in the
79 1 pfleura2
! formula.
80 1 pfleura2
! See for example:
81 1 pfleura2
! 1) Numerical Optimization, Springer 2nd Ed., 2006
82 1 pfleura2
! by Jorge Nocedal & Stephen J. Wright
83 1 pfleura2
! 2) ractical Methods of Optimization, John Wiley & Sons, second ed., 1987.
84 1 pfleura2
! by R. Fletcher
85 1 pfleura2
!
86 1 pfleura2
!
87 1 pfleura2
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
88 1 pfleura2
89 1 pfleura2
90 1 pfleura2
91 1 pfleura2
  interface
92 1 pfleura2
     function valid(string) result (isValid)
93 1 pfleura2
       logical                  :: isValid
94 1 pfleura2
       character(*), intent(in) :: string
95 1 pfleura2
     end function valid
96 1 pfleura2
97 8 pfleura2
98 8 pfleura2
    SUBROUTINE die(routine, msg, file, line, unit)
99 8 pfleura2
100 8 pfleura2
      Use VarTypes
101 8 pfleura2
      Use io_module
102 8 pfleura2
103 8 pfleura2
      implicit none
104 8 pfleura2
!--------------------------------------------------------------- Input Variables
105 8 pfleura2
      character(len=*), intent(in)           :: routine, msg
106 8 pfleura2
      character(len=*), intent(in), optional :: file
107 8 pfleura2
      integer(KINT), intent(in), optional      :: line, unit
108 8 pfleura2
109 8 pfleura2
    END SUBROUTINE die
110 8 pfleura2
111 8 pfleura2
    SUBROUTINE Warning(routine, msg, file, line, unit)
112 8 pfleura2
113 8 pfleura2
      Use VarTypes
114 8 pfleura2
      Use io_module
115 8 pfleura2
116 8 pfleura2
      implicit none
117 8 pfleura2
118 8 pfleura2
      character(len=*), intent(in)           :: routine, msg
119 8 pfleura2
      character(len=*), intent(in), optional :: file
120 8 pfleura2
      integer(KINT), intent(in), optional      :: line, unit
121 8 pfleura2
122 8 pfleura2
    END SUBROUTINE Warning
123 8 pfleura2
124 8 pfleura2
125 1 pfleura2
  end interface
126 1 pfleura2
127 1 pfleura2
  ! ======================================================================
128 1 pfleura2
  LOGICAL :: Debug
129 1 pfleura2
130 1 pfleura2
  ! ======================================================================
131 1 pfleura2
132 1 pfleura2
  debug = valid ('HESS') .OR. Valid('HUPDATE')
133 1 pfleura2
134 1 pfleura2
  if (debug) WRITE(*,*) "================================= Entering Hupdate_all ==============="
135 1 pfleura2
136 1 pfleura2
  SELECT CASE (Hinv)
137 1 pfleura2
     CASE (.TRUE.)
138 1 pfleura2
        SELECT CASE (HUpdate)
139 1 pfleura2
           CASE ("BFGS")
140 1 pfleura2
              Call Hinvup_BFGS(n,ds,dgrad,hess)
141 1 pfleura2
           CASE ("MS")
142 1 pfleura2
! we use the fact that MS is self-dual
143 1 pfleura2
              Call Hupdate_MS(n,dgrad,ds,hess)
144 1 pfleura2
           CASE ("DFP")
145 1 pfleura2
              Call Hinvup_DFP(n,ds,dgrad,hess)
146 1 pfleura2
           CASE DEFAULT
147 8 pfleura2
              Call Warning("Hupdate_all Hinv=T",Trim(HUpdate) // " not known: using BFGS" &
148 8 pfleura2
                   ,Unit=IoOut)
149 1 pfleura2
              Call Hinvup_BFGS(n,ds,dgrad,hess)
150 1 pfleura2
         END SELECT
151 1 pfleura2
     CASE (.FALSE.)
152 1 pfleura2
        SELECT CASE (HUpdate)
153 1 pfleura2
           CASE ("BOFILL")
154 1 pfleura2
              Call Hupdate_Bofill(n,ds,dgrad,hess)
155 1 pfleura2
           CASE ("MS")
156 1 pfleura2
              Call Hupdate_MS(n,ds,dgrad,hess)
157 1 pfleura2
           CASE ("DFP")
158 1 pfleura2
! we use the fact that DFP is the dual of BFGS
159 1 pfleura2
              Call Hinvup_BFGS(n,dgrad,ds,hess)
160 1 pfleura2
           CASE ("BFGS")
161 1 pfleura2
! we use the fact that DFP is the dual of BFGS
162 1 pfleura2
              Call Hinvup_DFP(n,dgrad,ds,hess)
163 1 pfleura2
           CASE DEFAULT
164 8 pfleura2
              Call Warning("Hupdate_all Hinv=F",Trim(HUpdate) // " not known: using BoFill" &
165 8 pfleura2
                   ,Unit=IoOut)
166 1 pfleura2
              Call Hupdate_Bofill(n,ds,dgrad,hess)
167 1 pfleura2
           END SELECT
168 1 pfleura2
    END SELECT
169 1 pfleura2
170 1 pfleura2
  if (debug) WRITE(*,*) "================================= Exiting Hupdate_all ==============="
171 1 pfleura2
end subroutine hupdate_all