Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief module handles definition of the tree nodes for the global and
10 : !> the subtrees binary tree
11 : !> --------------------------------------------------
12 : !> | parent element |
13 : !> | / \ |
14 : !> | accepted (acc) / \ not accepted (nacc) |
15 : !> | / \ |
16 : !> | child child |
17 : !> | / \ / \ |
18 : !> | ... |
19 : !> --------------------------------------------------
20 : !> tree creation assuming acceptance (acc) AND rejectance (nacc)
21 : !> of configuration
22 : !> if configuration is accepted: new configuration (child on acc) on basis
23 : !> of last configuration (one level up)
24 : !> if configuration is rejected: child on nacc on basis of last accepted
25 : !> element (last element which is on acc brach of its parent element)
26 : !> The global tree handles all configurations of different subtrees.
27 : !> The structure element "conf" is an array related to the temperature
28 : !> (sorted) and points to the subtree elements.
29 : !> \par History
30 : !> 11.2012 created [Mandes Schoenherr]
31 : !> \author Mandes
32 : ! **************************************************************************************************
33 :
34 : MODULE tmc_tree_types
35 : USE kinds, ONLY: dp
36 : #include "../base/base_uses.f90"
37 :
38 : IMPLICIT NONE
39 :
40 : PRIVATE
41 :
42 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'tmc_tree_types'
43 :
44 : PUBLIC :: tree_type, global_tree_type
45 : PUBLIC :: elem_list_type, elem_array_type, gt_elem_list_type
46 : PUBLIC :: add_to_list, clean_list
47 : PUBLIC :: read_subtree_elem_unformated, write_subtree_elem_unformated
48 :
49 : !-- tree element status
50 : INTEGER, PARAMETER, PUBLIC :: status_created = 100
51 : INTEGER, PARAMETER, PUBLIC :: status_calculate_energy = 101
52 : INTEGER, PARAMETER, PUBLIC :: status_calc_approx_ener = 102
53 :
54 : INTEGER, PARAMETER, PUBLIC :: status_calculate_NMC_steps = 111
55 : INTEGER, PARAMETER, PUBLIC :: status_calculate_MD = 112
56 : INTEGER, PARAMETER, PUBLIC :: status_calculated = 113
57 :
58 : INTEGER, PARAMETER, PUBLIC :: status_accepted_result = 123
59 : INTEGER, PARAMETER, PUBLIC :: status_accepted = 122
60 : INTEGER, PARAMETER, PUBLIC :: status_rejected = 121
61 : INTEGER, PARAMETER, PUBLIC :: status_rejected_result = 120
62 :
63 : INTEGER, PARAMETER, PUBLIC :: status_cancel_nmc = 133
64 : INTEGER, PARAMETER, PUBLIC :: status_cancel_ener = 132
65 : INTEGER, PARAMETER, PUBLIC :: status_canceled_nmc = 131
66 : INTEGER, PARAMETER, PUBLIC :: status_canceled_ener = 130
67 :
68 : INTEGER, PARAMETER, PUBLIC :: status_deleted = 140
69 : INTEGER, PARAMETER, PUBLIC :: status_deleted_result = 141
70 :
71 : !-- dimension status (for e.g. dividing atoms in sub box)
72 : INTEGER, PARAMETER, PUBLIC :: status_ok = 42
73 : INTEGER, PARAMETER, PUBLIC :: status_frozen = -1
74 : INTEGER, PARAMETER, PUBLIC :: status_proton_disorder = 1
75 :
76 : !-- subtree element
77 : TYPE tree_type
78 : TYPE(tree_type), POINTER :: parent => NULL() ! points to element one level up
79 : !-- acc..accepted goes to next level (next step),
80 : ! nacc..not accepted takes an alternative configutation
81 : TYPE(tree_type), POINTER :: acc => NULL(), nacc => NULL()
82 : !-- type of MC move (swap is handled only in global tree)
83 : INTEGER :: move_type = -1
84 : !-- status (e.g. calculated, MD calculation, accepted...)
85 : INTEGER :: stat = status_created
86 : REAL(KIND=dp), DIMENSION(:), POINTER :: subbox_center => NULL()
87 : REAL(KIND=dp), DIMENSION(:), POINTER :: pos => NULL() ! position array
88 : INTEGER, DIMENSION(:), POINTER :: mol => NULL() ! specifies the molecules the atoms participate
89 : REAL(KIND=dp), DIMENSION(:), POINTER :: vel => NULL() ! velocity array
90 : REAL(KIND=dp), DIMENSION(:), POINTER :: frc => NULL() ! force array
91 : REAL(KIND=dp), DIMENSION(:), POINTER :: dipole => NULL() ! dipole moments array
92 : INTEGER, DIMENSION(:), POINTER :: elem_stat => NULL() ! status for every dimension
93 : INTEGER :: nr = -1 ! tree node number
94 : REAL(KIND=dp), DIMENSION(3, 2, 3) :: rng_seed = 0 ! random seed for childs
95 : !-- remembers which subtree number element is from
96 : INTEGER :: sub_tree_nr = -1
97 : !-- remembers the temperature the configurational change (NMC) is done with
98 : INTEGER :: temp_created = 0
99 : !-- pointer to counter of next subtree element number
100 : INTEGER, POINTER :: next_elem_nr => NULL()
101 : !-- for calculating the NPT ensamble, variable box sizes are necessary.
102 : REAL(KIND=dp), DIMENSION(:), POINTER :: box_scale => NULL()
103 : REAL(KIND=dp) :: potential = 0.0_dp ! potential energy
104 : !-- potential energy calculated using (MD potential) cp2k input file
105 : REAL(KIND=dp) :: e_pot_approx = 0.0_dp
106 : !-- kinetic energy (espacially for HMC, where the velocities are respected)
107 : REAL(KIND=dp) :: ekin = 0.0_dp
108 : !-- kinetic energy before md steps (after gaussian velocity change)
109 : REAL(KIND=dp) :: ekin_before_md = 0.0_dp
110 : !-- estimated energies are stored in loop order in this array
111 : REAL(KIND=dp), DIMENSION(4) :: scf_energies = 0.0_dp
112 : !-- counter to get last position in the array loop
113 : INTEGER :: scf_energies_count = 0
114 : !-- list of global tree elements referint to that node (reference back to global tree)
115 : ! if no reference exist anymore, global tree element can be deleted
116 : TYPE(gt_elem_list_type), POINTER :: gt_nodes_references => NULL()
117 : END TYPE tree_type
118 :
119 : ! type for global tree element list in tree elements
120 : TYPE gt_elem_list_type
121 : TYPE(global_tree_type), POINTER :: gt_elem => NULL()
122 : TYPE(gt_elem_list_type), POINTER :: next => NULL()
123 : END TYPE gt_elem_list_type
124 :
125 : TYPE elem_list_type
126 : TYPE(tree_type), POINTER :: elem => NULL()
127 : TYPE(elem_list_type), POINTER :: next => NULL()
128 : INTEGER :: temp_ind = 0
129 : INTEGER :: nr = -1
130 : END TYPE elem_list_type
131 :
132 : !-- array with subtree elements
133 : TYPE elem_array_type
134 : TYPE(tree_type), POINTER :: elem => NULL()
135 : LOGICAL :: busy = .FALSE.
136 : LOGICAL :: canceled = .FALSE.
137 : REAL(KIND=dp) :: start_time = 0.0_dp
138 : END TYPE elem_array_type
139 :
140 : !-- global tree element
141 : TYPE global_tree_type
142 : TYPE(global_tree_type), POINTER :: parent => NULL() ! points to element one level up
143 : !-- acc..accepted goes to next level (next step),
144 : ! nacc..not accepted takes an alternative configutation
145 : TYPE(global_tree_type), POINTER :: acc => NULL(), nacc => NULL()
146 : !-- status (e.g. calculated, MD calculation, accepted...)
147 : INTEGER :: stat = -99
148 : !-- remember if configuration in node are swaped
149 : LOGICAL :: swaped = .FALSE.
150 : !-- stores the index of the configuration (temperature)
151 : ! which is changed
152 : INTEGER :: mv_conf = -54321
153 : !-- stores the index of the configuration (temp.) which should change next
154 : INTEGER :: mv_next_conf = -2345
155 : !-- list of pointes to subtree elements (Temp sorting)
156 : TYPE(elem_array_type), DIMENSION(:), ALLOCATABLE :: conf
157 : !-- remembers if last configuration is assumed to be accepted or rejected (next branc in tree);
158 : ! In case of swaping, it shows if the configuration of a certain temperature is assumed
159 : ! to be acc/rej (which branch is followed at the last modification of the conf of this temp.
160 : !TODO store conf_n_acc in a bitshifted array to decrease the size (1Logical = 1Byte)
161 : LOGICAL, DIMENSION(:), ALLOCATABLE :: conf_n_acc
162 : INTEGER :: nr = 0 ! tree node number
163 : REAL(KIND=dp), DIMENSION(3, 2, 3) :: rng_seed = 0.0_dp ! random seed for childs
164 : !-- random number for acceptance check
165 : REAL(KIND=dp) :: rnd_nr = 0.0_dp
166 : !-- approximate probability of acceptance will be adapted while calculating the exact energy
167 : REAL(KIND=dp) :: prob_acc = 0.0_dp ! estimated acceptance probability
168 : REAL(KIND=dp) :: Temp = 0.0_dp ! temperature for simulated annealing
169 : END TYPE global_tree_type
170 :
171 : CONTAINS
172 :
173 : ! **************************************************************************************************
174 : !> \brief add a certain element to the specified element list at the beginning
175 : !> \param elem the sub tree element, to be added
176 : !> \param list ...
177 : !> \param temp_ind ...
178 : !> \param nr ...
179 : !> \author Mandes 11.2012
180 : ! **************************************************************************************************
181 1 : SUBROUTINE add_to_list(elem, list, temp_ind, nr)
182 : TYPE(tree_type), POINTER :: elem
183 : TYPE(elem_list_type), POINTER :: list
184 : INTEGER, OPTIONAL :: temp_ind, nr
185 :
186 : TYPE(elem_list_type), POINTER :: last, list_elem_tmp
187 :
188 1 : NULLIFY (list_elem_tmp, last)
189 :
190 1 : CPASSERT(ASSOCIATED(elem))
191 :
192 1 : ALLOCATE (list_elem_tmp)
193 1 : list_elem_tmp%elem => elem
194 : list_elem_tmp%next => NULL()
195 1 : IF (PRESENT(temp_ind)) THEN
196 0 : list_elem_tmp%temp_ind = temp_ind
197 : ELSE
198 1 : list_elem_tmp%temp_ind = -1
199 : END IF
200 :
201 1 : IF (PRESENT(nr)) THEN
202 0 : list_elem_tmp%nr = nr
203 : ELSE
204 : list_elem_tmp%nr = -1
205 : END IF
206 :
207 1 : IF (ASSOCIATED(list) .EQV. .FALSE.) THEN
208 1 : list => list_elem_tmp
209 : ELSE
210 : last => list
211 0 : DO WHILE (ASSOCIATED(last%next))
212 0 : last => last%next
213 : END DO
214 0 : last%next => list_elem_tmp
215 : END IF
216 :
217 1 : END SUBROUTINE add_to_list
218 :
219 : ! **************************************************************************************************
220 : !> \brief clean a certain element element list
221 : !> \param list ...
222 : !> \author Mandes 11.2012
223 : ! **************************************************************************************************
224 28 : SUBROUTINE clean_list(list)
225 : TYPE(elem_list_type), POINTER :: list
226 :
227 : TYPE(elem_list_type), POINTER :: list_elem_tmp
228 :
229 28 : NULLIFY (list_elem_tmp)
230 :
231 28 : DO WHILE (ASSOCIATED(list))
232 0 : list_elem_tmp => list%next
233 0 : DEALLOCATE (list)
234 0 : list => list_elem_tmp
235 : END DO
236 28 : END SUBROUTINE clean_list
237 :
238 : ! **************************************************************************************************
239 : !> \brief prints out the TMC sub tree structure element unformated in file
240 : !> \param elem ...
241 : !> \param io_unit ...
242 : !> \param
243 : !> \author Mandes 11.2012
244 : ! **************************************************************************************************
245 6 : SUBROUTINE write_subtree_elem_unformated(elem, io_unit)
246 : TYPE(tree_type), POINTER :: elem
247 : INTEGER :: io_unit
248 :
249 6 : CPASSERT(ASSOCIATED(elem))
250 6 : CPASSERT(io_unit > 0)
251 6 : WRITE (io_unit) elem%nr, &
252 6 : elem%sub_tree_nr, &
253 6 : elem%stat, &
254 168 : elem%rng_seed, &
255 6 : elem%move_type, &
256 6 : elem%temp_created, &
257 6 : elem%potential, &
258 6 : elem%e_pot_approx, &
259 6 : elem%ekin, &
260 12 : elem%ekin_before_md
261 6 : CALL write_subtree_elem_darray(elem%pos, io_unit)
262 6 : CALL write_subtree_elem_darray(elem%vel, io_unit)
263 6 : CALL write_subtree_elem_darray(elem%frc, io_unit)
264 6 : CALL write_subtree_elem_darray(elem%box_scale, io_unit)
265 6 : CALL write_subtree_elem_darray(elem%dipole, io_unit)
266 6 : END SUBROUTINE write_subtree_elem_unformated
267 :
268 : ! **************************************************************************************************
269 : !> \brief reads the TMC sub tree structure element unformated in file
270 : !> \param elem ...
271 : !> \param io_unit ...
272 : !> \param
273 : !> \author Mandes 11.2012
274 : ! **************************************************************************************************
275 3 : SUBROUTINE read_subtree_elem_unformated(elem, io_unit)
276 : TYPE(tree_type), POINTER :: elem
277 : INTEGER :: io_unit
278 :
279 3 : CPASSERT(ASSOCIATED(elem))
280 3 : CPASSERT(io_unit > 0)
281 :
282 3 : READ (io_unit) elem%nr, &
283 3 : elem%sub_tree_nr, &
284 3 : elem%stat, &
285 84 : elem%rng_seed, &
286 3 : elem%move_type, &
287 3 : elem%temp_created, &
288 3 : elem%potential, &
289 3 : elem%e_pot_approx, &
290 3 : elem%ekin, &
291 6 : elem%ekin_before_md
292 3 : CALL read_subtree_elem_darray(elem%pos, io_unit)
293 3 : CALL read_subtree_elem_darray(elem%vel, io_unit)
294 3 : CALL read_subtree_elem_darray(elem%frc, io_unit)
295 3 : CALL read_subtree_elem_darray(elem%box_scale, io_unit)
296 3 : CALL read_subtree_elem_darray(elem%dipole, io_unit)
297 3 : END SUBROUTINE read_subtree_elem_unformated
298 :
299 : ! **************************************************************************************************
300 : !> \brief ...
301 : !> \param array ...
302 : !> \param io_unit ...
303 : ! **************************************************************************************************
304 30 : SUBROUTINE write_subtree_elem_darray(array, io_unit)
305 : REAL(KIND=dp), DIMENSION(:), POINTER :: array
306 : INTEGER :: io_unit
307 :
308 30 : WRITE (io_unit) ASSOCIATED(array)
309 30 : IF (ASSOCIATED(array)) THEN
310 18 : WRITE (io_unit) SIZE(array)
311 792 : WRITE (io_unit) array
312 : END IF
313 30 : END SUBROUTINE write_subtree_elem_darray
314 :
315 : ! **************************************************************************************************
316 : !> \brief ...
317 : !> \param array ...
318 : !> \param io_unit ...
319 : ! **************************************************************************************************
320 15 : SUBROUTINE read_subtree_elem_darray(array, io_unit)
321 : REAL(KIND=dp), DIMENSION(:), POINTER :: array
322 : INTEGER :: io_unit
323 :
324 : INTEGER :: i_tmp
325 : LOGICAL :: l_tmp
326 :
327 15 : READ (io_unit) l_tmp
328 15 : IF (l_tmp) THEN
329 9 : READ (io_unit) i_tmp
330 9 : IF (ASSOCIATED(array)) THEN
331 9 : CPASSERT(SIZE(array) == i_tmp)
332 : ELSE
333 0 : ALLOCATE (array(i_tmp))
334 : END IF
335 396 : READ (io_unit) array
336 : END IF
337 15 : END SUBROUTINE read_subtree_elem_darray
338 :
339 0 : END MODULE tmc_tree_types
|