(git:58e3e09)
atom_optimization.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief Optimizer for the atomic code
10 ! **************************************************************************************************
12  USE atom_types, ONLY: atom_optimization_type
13  USE kinds, ONLY: dp
14  USE lapack, ONLY: lapack_sgelss
15 #include "./base/base_uses.f90"
16 
17  IMPLICIT NONE
18 
19  PRIVATE
20 
21  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'atom_optimization'
22 
23  TYPE hmat_type
24  REAL(dp) :: energy = 0.0_dp
25  REAL(dp) :: error = 0.0_dp
26  REAL(dp), DIMENSION(:, :, :), POINTER :: emat => null()
27  REAL(dp), DIMENSION(:, :, :), POINTER :: fmat => null()
28  REAL(dp), DIMENSION(:, :, :), POINTER :: pmat => null()
29  END TYPE hmat_type
30 
31  TYPE atom_history_type
32  INTEGER :: max_history = 0
33  INTEGER :: hlen = 0
34  INTEGER :: hpos = 0
35  REAL(dp) :: damping = 0.0_dp
36  REAL(dp) :: eps_diis = 0.0_dp
37  REAL(dp), DIMENSION(:, :), POINTER :: dmat => null()
38  TYPE(hmat_type), DIMENSION(:), POINTER :: hmat => null()
39  END TYPE atom_history_type
40 
41  PUBLIC :: atom_opt_fmat, &
43 
44 CONTAINS
45 
46 ! **************************************************************************************************
47 !> \brief Initialise a circular buffer to keep Kohn-Sham and density matrices from previous iteration.
48 !> \param history object to initialise
49 !> \param optimization optimisation parameters
50 !> \param matrix reference matrix. Historic matrices will have the same size as
51 !> this reference matrix
52 !> \par History
53 !> * 08.2016 new structure element: density matrix [Juerg Hutter]
54 !> * 08.2008 created [Juerg Hutter]
55 ! **************************************************************************************************
56  PURE SUBROUTINE atom_history_init(history, optimization, matrix)
57  TYPE(atom_history_type), INTENT(INOUT) :: history
58  TYPE(atom_optimization_type), INTENT(IN) :: optimization
59  REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: matrix
60 
61  INTEGER :: i, n1, n2, n3, ndiis
62  REAL(dp) :: damp, eps
63 
64  ndiis = optimization%n_diis
65  eps = optimization%eps_diis
66  damp = optimization%damping
67 
68  CALL atom_history_release(history)
69 
70  history%max_history = ndiis
71  history%hlen = 0
72  history%hpos = 0
73  history%damping = damp
74  history%eps_diis = eps
75  ALLOCATE (history%dmat(ndiis + 1, ndiis + 1))
76 
77  ALLOCATE (history%hmat(ndiis))
78  n1 = SIZE(matrix, 1)
79  n2 = SIZE(matrix, 2)
80  n3 = SIZE(matrix, 3)
81  DO i = 1, ndiis
82  history%hmat(i)%energy = 0.0_dp
83  history%hmat(i)%error = 0.0_dp
84  ALLOCATE (history%hmat(i)%emat(n1, n2, n3))
85  ALLOCATE (history%hmat(i)%fmat(n1, n2, n3))
86  ALLOCATE (history%hmat(i)%pmat(n1, n2, n3))
87  END DO
88 
89  END SUBROUTINE atom_history_init
90 
91 ! **************************************************************************************************
92 !> \brief Add matrices from the current iteration into the circular buffer.
93 !> \param history object to keep historic matrices
94 !> \param pmat density matrix
95 !> \param fmat Kohn-Sham matrix
96 !> \param emat error matrix
97 !> \param energy total energy
98 !> \param error convergence
99 !> \par History
100 !> * 08.2016 new formal argument: density matrix [Juerg Hutter]
101 !> * 08.2008 created [Juerg Hutter]
102 ! **************************************************************************************************
103  PURE SUBROUTINE atom_history_update(history, pmat, fmat, emat, energy, error)
104  TYPE(atom_history_type), INTENT(INOUT) :: history
105  REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: pmat, fmat, emat
106  REAL(dp), INTENT(IN) :: energy, error
107 
108  INTEGER :: nlen, nmax, nnow
109 
110  nmax = history%max_history
111  nlen = min(history%hlen + 1, nmax)
112  nnow = history%hpos + 1
113  IF (nnow > nmax) nnow = 1
114 
115  history%hmat(nnow)%energy = energy
116  history%hmat(nnow)%error = error
117  history%hmat(nnow)%pmat = pmat
118  history%hmat(nnow)%fmat = fmat
119  history%hmat(nnow)%emat = emat
120 
121  history%hlen = nlen
122  history%hpos = nnow
123 
124  END SUBROUTINE atom_history_update
125 
126 ! **************************************************************************************************
127 !> \brief Release circular buffer to keep historic matrices.
128 !> \param history object to release
129 !> \par History
130 !> * 08.2008 created [Juerg Hutter]
131 ! **************************************************************************************************
132  PURE SUBROUTINE atom_history_release(history)
133  TYPE(atom_history_type), INTENT(INOUT) :: history
134 
135  INTEGER :: i
136 
137  history%max_history = 0
138  history%hlen = 0
139  history%hpos = 0
140  history%damping = 0._dp
141  history%eps_diis = 0._dp
142  IF (ASSOCIATED(history%dmat)) THEN
143  DEALLOCATE (history%dmat)
144  END IF
145  IF (ASSOCIATED(history%hmat)) THEN
146  DO i = 1, SIZE(history%hmat)
147  IF (ASSOCIATED(history%hmat(i)%emat)) THEN
148  DEALLOCATE (history%hmat(i)%emat)
149  END IF
150  IF (ASSOCIATED(history%hmat(i)%fmat)) THEN
151  DEALLOCATE (history%hmat(i)%fmat)
152  END IF
153  IF (ASSOCIATED(history%hmat(i)%pmat)) THEN
154  DEALLOCATE (history%hmat(i)%pmat)
155  END IF
156  END DO
157  DEALLOCATE (history%hmat)
158  END IF
159 
160  END SUBROUTINE atom_history_release
161 
162 ! **************************************************************************************************
163 !> \brief Construct a Kohn-Sham matrix for the next iteration based on the historic data.
164 !> \param fmat new Kohn-Sham matrix
165 !> \param history historic matrices
166 !> \param err convergence
167 !> \par History
168 !> * 08.2016 renamed to atom_opt_fmat() [Juerg Hutter]
169 !> * 08.2008 created as atom_opt() [Juerg Hutter]
170 ! **************************************************************************************************
171  SUBROUTINE atom_opt_fmat(fmat, history, err)
172  REAL(dp), DIMENSION(:, :, :), INTENT(INOUT) :: fmat
173  TYPE(atom_history_type), INTENT(INOUT) :: history
174  REAL(dp), INTENT(IN) :: err
175 
176  INTEGER :: i, info, j, lwork, na, nb, nlen, nm, &
177  nmax, nnow, rank
178  REAL(dp) :: a, rcond, t
179  REAL(dp), ALLOCATABLE, DIMENSION(:) :: s, work
180  REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: vec
181 
182  nmax = history%max_history
183  nnow = history%hpos
184  a = history%damping
185  IF (history%hlen > 1) THEN
186  IF (err < history%eps_diis) THEN
187  ! DIIS
188  rcond = 1.e-10_dp
189  lwork = 25*nmax
190  ALLOCATE (vec(nmax + 1, 2), s(nmax + 1), work(lwork))
191  nlen = history%hlen
192  vec = 0._dp
193  vec(nlen + 1, 1) = 1._dp
194  history%dmat(1:nlen, nlen + 1) = 1._dp
195  history%dmat(nlen + 1, 1:nlen) = 1._dp
196  history%dmat(nlen + 1, nlen + 1) = 0._dp
197  DO i = 1, nlen
198  na = nnow + 1 - i
199  IF (na < 1) na = nmax + na
200  DO j = i, nlen
201  nb = nnow + 1 - j
202  IF (nb < 1) nb = nmax + nb
203  t = sum(history%hmat(na)%emat*history%hmat(nb)%emat)
204  history%dmat(i, j) = t
205  history%dmat(j, i) = t
206  END DO
207  END DO
208  CALL lapack_sgelss(nlen + 1, nlen + 1, 1, history%dmat, nmax + 1, vec, nmax + 1, s, &
209  rcond, rank, work, lwork, info)
210  cpassert(info == 0)
211  fmat = 0._dp
212  DO i = 1, nlen
213  na = nnow + 1 - i
214  IF (na < 1) na = nmax + na
215  fmat = fmat + vec(i, 1)*history%hmat(na)%fmat
216  END DO
217 
218  DEALLOCATE (vec, s, work)
219  ELSE
220  ! damping
221  nm = nnow - 1
222  IF (nm < 1) nm = history%max_history
223  fmat = a*history%hmat(nnow)%fmat + (1._dp - a)*history%hmat(nm)%fmat
224  END IF
225  ELSEIF (history%hlen == 1) THEN
226  fmat = history%hmat(nnow)%fmat
227  ELSE
228  cpabort("")
229  END IF
230 
231  END SUBROUTINE atom_opt_fmat
232 
233 END MODULE atom_optimization
Optimizer for the atomic code.
subroutine, public atom_opt_fmat(fmat, history, err)
Construct a Kohn-Sham matrix for the next iteration based on the historic data.
pure subroutine, public atom_history_init(history, optimization, matrix)
Initialise a circular buffer to keep Kohn-Sham and density matrices from previous iteration.
pure subroutine, public atom_history_release(history)
Release circular buffer to keep historic matrices.
pure subroutine, public atom_history_update(history, pmat, fmat, emat, energy, error)
Add matrices from the current iteration into the circular buffer.
Define the atom type and its sub types.
Definition: atom_types.F:15
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
Interface to the LAPACK F77 library.
Definition: lapack.F:17