Statistiques
| Révision :

root / src / Hupdate_all.f90 @ 8

Historique | Voir | Annoter | Télécharger (4,29 ko)

1
 subroutine hupdate_all(n, ds, dgrad, Hess)
2

    
3
  use Path_module, only :  HUpdate,Hinv
4
  use Io_module, only : IoOut
5

    
6
  IMPLICIT NONE
7

    
8
  INTEGER, PARAMETER :: KINT=KIND(1)
9
  INTEGER, PARAMETER :: KREAL=KIND(1.0D0)
10

    
11
  INTEGER(KINT), INTENT(IN) :: n
12
  real(KREAL), INTENT(IN) :: ds(n), dgrad(n)
13
  REAL(KREAL),INTENT(INOUT) :: Hess(n:n)
14

    
15
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16
!
17
! This subroutine is a driver for the update of the Hessian 
18
! in the Quasi-Newton optimization step.
19
! We propose mostly update of the inverse Hessian.
20
!
21
! input:
22
!  - n: number of coordinates
23
!  - ds: GeomOld-Geom (expressed in the optimization coord.)
24
!  - dgrad: GradOld-Grad  (expressed in the optimization coord.)
25
!
26
!  input/output:
27
!  - Hess: old Hessian in input, replaced by updated Hessian on output
28
!
29
!!!!!!!!!!!!!!!!
30
!
31
! For now, we propose :
32
! - For update of inverse Hessian:
33
!      - BFGS (Broyden, Fletcher, Goldfarb, and Shanno)
34
!      - DFP (Davidson - Fletcher - Powell )
35
!      - MS (Murtagh-Sargent)
36
!
37
! - For update of hessian:
38
!      - Bofill (combination of Murtagh-Sargent and Powell-symmetric-Broyden
39
!      - MS (Murtagh-Sargent)
40
!      - BFGS
41
!      - DFP
42
!----------------------------------------------------
43
! We use the fact that some formula are the dual of other
44
! that is, one can convert the formula for update of the Hessian 
45
! into the update of the inverse Hessian by replacing the Hessian
46
! (denoted by B in the Fletcher book and in the Nocedal&Wrigth book)
47
! by the inverse hessian (denoted by H), and replacing ds by dgrad in the
48
! formula.
49
! See for example:
50
! 1) Numerical Optimization, Springer 2nd Ed., 2006
51
! by Jorge Nocedal & Stephen J. Wright
52
! 2) ractical Methods of Optimization, John Wiley & Sons, second ed., 1987.
53
! by R. Fletcher
54
!
55
!
56
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
57

    
58

    
59

    
60
  interface
61
     function valid(string) result (isValid)
62
       logical                  :: isValid
63
       character(*), intent(in) :: string
64
     end function valid
65

    
66

    
67
    SUBROUTINE die(routine, msg, file, line, unit)
68

    
69
      Use VarTypes
70
      Use io_module
71

    
72
      implicit none
73
!--------------------------------------------------------------- Input Variables
74
      character(len=*), intent(in)           :: routine, msg
75
      character(len=*), intent(in), optional :: file
76
      integer(KINT), intent(in), optional      :: line, unit
77

    
78
    END SUBROUTINE die
79

    
80
    SUBROUTINE Warning(routine, msg, file, line, unit)
81

    
82
      Use VarTypes
83
      Use io_module
84

    
85
      implicit none
86

    
87
      character(len=*), intent(in)           :: routine, msg
88
      character(len=*), intent(in), optional :: file
89
      integer(KINT), intent(in), optional      :: line, unit
90

    
91
    END SUBROUTINE Warning
92

    
93

    
94
  end interface
95

    
96
  ! ======================================================================
97
  LOGICAL :: Debug
98

    
99
  ! ======================================================================
100

    
101
  debug = valid ('HESS') .OR. Valid('HUPDATE')
102

    
103
  if (debug) WRITE(*,*) "================================= Entering Hupdate_all ==============="
104

    
105
  SELECT CASE (Hinv)
106
     CASE (.TRUE.)
107
        SELECT CASE (HUpdate)
108
           CASE ("BFGS")
109
              Call Hinvup_BFGS(n,ds,dgrad,hess)
110
           CASE ("MS")
111
! we use the fact that MS is self-dual 
112
              Call Hupdate_MS(n,dgrad,ds,hess)
113
           CASE ("DFP")
114
              Call Hinvup_DFP(n,ds,dgrad,hess)
115
           CASE DEFAULT
116
              Call Warning("Hupdate_all Hinv=T",Trim(HUpdate) // " not known: using BFGS" &
117
                   ,Unit=IoOut)
118
              Call Hinvup_BFGS(n,ds,dgrad,hess)
119
         END SELECT
120
     CASE (.FALSE.)
121
        SELECT CASE (HUpdate)
122
           CASE ("BOFILL")
123
              Call Hupdate_Bofill(n,ds,dgrad,hess)
124
           CASE ("MS")
125
              Call Hupdate_MS(n,ds,dgrad,hess)
126
           CASE ("DFP")
127
! we use the fact that DFP is the dual of BFGS
128
              Call Hinvup_BFGS(n,dgrad,ds,hess)
129
           CASE ("BFGS")
130
! we use the fact that DFP is the dual of BFGS
131
              Call Hinvup_DFP(n,dgrad,ds,hess)
132
           CASE DEFAULT
133
              Call Warning("Hupdate_all Hinv=F",Trim(HUpdate) // " not known: using BoFill" &
134
                   ,Unit=IoOut)
135
              Call Hupdate_Bofill(n,ds,dgrad,hess)
136
           END SELECT
137
    END SELECT
138

    
139
  if (debug) WRITE(*,*) "================================= Exiting Hupdate_all ==============="
140
end subroutine hupdate_all