Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Handles all functions related to the CELL
10 : !> \par History
11 : !> 11.2008 Teodoro Laino [tlaino] - deeply cleaning cell_type from units
12 : !> 10.2014 Moved many routines from cell_types.F here.
13 : !> \author Matthias KracK (16.01.2002, based on a earlier version of CJM, JGH)
14 : ! **************************************************************************************************
15 : MODULE cell_types
16 : USE cp_units, ONLY: cp_unit_to_cp2k
17 : USE kinds, ONLY: dp
18 : USE mathconstants, ONLY: degree
19 : USE mathlib, ONLY: angle
20 : #include "../base/base_uses.f90"
21 :
22 : IMPLICIT NONE
23 :
24 : PRIVATE
25 :
26 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cell_types'
27 :
28 : ! Impose cell symmetry
29 : INTEGER, PARAMETER, PUBLIC :: cell_sym_none = 0, &
30 : cell_sym_triclinic = 1, &
31 : cell_sym_monoclinic = 2, &
32 : cell_sym_monoclinic_gamma_ab = 3, &
33 : cell_sym_orthorhombic = 4, &
34 : cell_sym_tetragonal_ab = 5, &
35 : cell_sym_tetragonal_ac = 6, &
36 : cell_sym_tetragonal_bc = 7, &
37 : cell_sym_rhombohedral = 8, &
38 : cell_sym_hexagonal_gamma_60 = 9, &
39 : cell_sym_hexagonal_gamma_120 = 10, &
40 : cell_sym_cubic = 11
41 :
42 : INTEGER, PARAMETER, PUBLIC :: use_perd_none = 0, &
43 : use_perd_x = 1, &
44 : use_perd_y = 2, &
45 : use_perd_z = 3, &
46 : use_perd_xy = 4, &
47 : use_perd_xz = 5, &
48 : use_perd_yz = 6, &
49 : use_perd_xyz = 7
50 :
51 : CHARACTER(LEN=3), DIMENSION(7), &
52 : PARAMETER, PUBLIC :: periodicity_string = [" X", " Y", " Z", &
53 : " XY", " XZ", " YZ", &
54 : "XYZ"]
55 :
56 : ! **************************************************************************************************
57 : !> \brief Type defining parameters related to the simulation cell
58 : !> \version 1.0
59 : ! **************************************************************************************************
60 : TYPE cell_type
61 : CHARACTER(LEN=12) :: tag = "CELL"
62 : INTEGER :: ref_count = -1, &
63 : symmetry_id = use_perd_none
64 : LOGICAL :: orthorhombic = .FALSE. ! actually means a diagonal hmat
65 : REAL(KIND=dp) :: deth = 0.0_dp
66 : INTEGER, DIMENSION(3) :: perd = -1
67 : REAL(KIND=dp), DIMENSION(3, 3) :: hmat = 0.0_dp, &
68 : h_inv = 0.0_dp
69 : END TYPE cell_type
70 :
71 : TYPE cell_p_type
72 : TYPE(cell_type), POINTER :: cell => NULL()
73 : END TYPE cell_p_type
74 :
75 : ! Public data types
76 : PUBLIC :: cell_type, &
77 : cell_p_type
78 :
79 : ! Public subroutines
80 : PUBLIC :: cell_clone, &
81 : cell_copy, &
82 : cell_release, &
83 : cell_retain, &
84 : get_cell, &
85 : parse_cell_line
86 :
87 : #if defined (__PLUMED2)
88 : PUBLIC :: pbc_cp2k_plumed_getset_cell
89 : #endif
90 :
91 : ! Public functions
92 : PUBLIC :: plane_distance, &
93 : pbc, &
94 : real_to_scaled, &
95 : scaled_to_real
96 :
97 : INTERFACE pbc
98 : MODULE PROCEDURE pbc1, pbc2, pbc3, pbc4
99 : END INTERFACE
100 :
101 : CONTAINS
102 :
103 : ! **************************************************************************************************
104 : !> \brief Clone cell variable
105 : !> \param cell_in Cell variable to be clone
106 : !> \param cell_out Cloned cell variable
107 : !> \param tag Optional new tag for cloned cell variable
108 : !> \par History
109 : !> - Optional tag added (17.05.2023, MK)
110 : ! **************************************************************************************************
111 31513 : SUBROUTINE cell_clone(cell_in, cell_out, tag)
112 :
113 : TYPE(cell_type), POINTER :: cell_in, cell_out
114 : CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: tag
115 :
116 31513 : cell_out = cell_in
117 31513 : cell_out%ref_count = 1
118 20331 : IF (PRESENT(tag)) cell_out%tag = tag
119 :
120 31513 : END SUBROUTINE cell_clone
121 :
122 : ! **************************************************************************************************
123 : !> \brief Copy cell variable
124 : !> \param cell_in Cell variable to be copied
125 : !> \param cell_out Copy of cell variable
126 : !> \param tag Optional new tag
127 : !> \par History
128 : !> - Optional tag added (17.05.2023, MK)
129 : ! **************************************************************************************************
130 228972 : SUBROUTINE cell_copy(cell_in, cell_out, tag)
131 :
132 : TYPE(cell_type), POINTER :: cell_in, cell_out
133 : CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: tag
134 :
135 228972 : cell_out%deth = cell_in%deth
136 1831776 : cell_out%perd = cell_in%perd
137 5953272 : cell_out%hmat = cell_in%hmat
138 5953272 : cell_out%h_inv = cell_in%h_inv
139 228972 : cell_out%orthorhombic = cell_in%orthorhombic
140 228972 : cell_out%symmetry_id = cell_in%symmetry_id
141 228972 : IF (PRESENT(tag)) THEN
142 9312 : cell_out%tag = tag
143 : ELSE
144 219660 : cell_out%tag = cell_in%tag
145 : END IF
146 :
147 228972 : END SUBROUTINE cell_copy
148 :
149 : ! **************************************************************************************************
150 : !> \brief Read cell info from a line (parsed from a file)
151 : !> \param input_line ...
152 : !> \param cell_itimes ...
153 : !> \param cell_time ...
154 : !> \param h ...
155 : !> \param vol ...
156 : !> \date 19.02.2008
157 : !> \author Teodoro Laino [tlaino] - University of Zurich
158 : !> \version 1.0
159 : ! **************************************************************************************************
160 368 : SUBROUTINE parse_cell_line(input_line, cell_itimes, cell_time, h, vol)
161 :
162 : CHARACTER(LEN=*), INTENT(IN) :: input_line
163 : INTEGER, INTENT(OUT) :: cell_itimes
164 : REAL(KIND=dp), INTENT(OUT) :: cell_time
165 : REAL(KIND=dp), DIMENSION(3, 3), INTENT(OUT) :: h
166 : REAL(KIND=dp), INTENT(OUT) :: vol
167 :
168 : INTEGER :: i, j
169 :
170 368 : READ (input_line, *) cell_itimes, cell_time, &
171 736 : h(1, 1), h(2, 1), h(3, 1), h(1, 2), h(2, 2), h(3, 2), h(1, 3), h(2, 3), h(3, 3), vol
172 1472 : DO i = 1, 3
173 4784 : DO j = 1, 3
174 4416 : h(j, i) = cp_unit_to_cp2k(h(j, i), "angstrom")
175 : END DO
176 : END DO
177 :
178 368 : END SUBROUTINE parse_cell_line
179 :
180 : ! **************************************************************************************************
181 : !> \brief Get informations about a simulation cell.
182 : !> \param cell ...
183 : !> \param alpha ...
184 : !> \param beta ...
185 : !> \param gamma ...
186 : !> \param deth ...
187 : !> \param orthorhombic ...
188 : !> \param abc ...
189 : !> \param periodic ...
190 : !> \param h ...
191 : !> \param h_inv ...
192 : !> \param symmetry_id ...
193 : !> \param tag ...
194 : !> \date 16.01.2002
195 : !> \author Matthias Krack
196 : !> \version 1.0
197 : ! **************************************************************************************************
198 133264069 : SUBROUTINE get_cell(cell, alpha, beta, gamma, deth, orthorhombic, abc, periodic, &
199 : h, h_inv, symmetry_id, tag)
200 :
201 : TYPE(cell_type), POINTER :: cell
202 : REAL(KIND=dp), INTENT(OUT), OPTIONAL :: alpha, beta, gamma, deth
203 : LOGICAL, INTENT(OUT), OPTIONAL :: orthorhombic
204 : REAL(KIND=dp), DIMENSION(3), INTENT(OUT), OPTIONAL :: abc
205 : INTEGER, DIMENSION(3), INTENT(OUT), OPTIONAL :: periodic
206 : REAL(KIND=dp), DIMENSION(3, 3), INTENT(OUT), &
207 : OPTIONAL :: h, h_inv
208 : INTEGER, INTENT(OUT), OPTIONAL :: symmetry_id
209 : CHARACTER(LEN=*), INTENT(OUT), OPTIONAL :: tag
210 :
211 0 : CPASSERT(ASSOCIATED(cell))
212 :
213 133264069 : IF (PRESENT(deth)) deth = cell%deth ! the volume
214 133264069 : IF (PRESENT(orthorhombic)) orthorhombic = cell%orthorhombic
215 511231690 : IF (PRESENT(periodic)) periodic(:) = cell%perd(:)
216 133431865 : IF (PRESENT(h)) h(:, :) = cell%hmat(:, :)
217 133264189 : IF (PRESENT(h_inv)) h_inv(:, :) = cell%h_inv(:, :)
218 :
219 : ! Calculate the lengths of the cell vectors a, b, and c
220 133264069 : IF (PRESENT(abc)) THEN
221 : abc(1) = SQRT(cell%hmat(1, 1)*cell%hmat(1, 1) + &
222 : cell%hmat(2, 1)*cell%hmat(2, 1) + &
223 7125996 : cell%hmat(3, 1)*cell%hmat(3, 1))
224 : abc(2) = SQRT(cell%hmat(1, 2)*cell%hmat(1, 2) + &
225 : cell%hmat(2, 2)*cell%hmat(2, 2) + &
226 7125996 : cell%hmat(3, 2)*cell%hmat(3, 2))
227 : abc(3) = SQRT(cell%hmat(1, 3)*cell%hmat(1, 3) + &
228 : cell%hmat(2, 3)*cell%hmat(2, 3) + &
229 7125996 : cell%hmat(3, 3)*cell%hmat(3, 3))
230 : END IF
231 :
232 : ! Angles between the cell vectors a, b, and c
233 : ! alpha = <(b,c)
234 133264069 : IF (PRESENT(alpha)) alpha = angle(cell%hmat(:, 2), cell%hmat(:, 3))*degree
235 : ! beta = <(a,c)
236 133264069 : IF (PRESENT(beta)) beta = angle(cell%hmat(:, 1), cell%hmat(:, 3))*degree
237 : ! gamma = <(a,b)
238 133264069 : IF (PRESENT(gamma)) gamma = angle(cell%hmat(:, 1), cell%hmat(:, 2))*degree
239 133264069 : IF (PRESENT(symmetry_id)) symmetry_id = cell%symmetry_id
240 133264069 : IF (PRESENT(tag)) tag = cell%tag
241 :
242 133264069 : END SUBROUTINE get_cell
243 :
244 : ! **************************************************************************************************
245 : !> \brief Calculate the distance between two lattice planes as defined by
246 : !> a triple of Miller indices (hkl).
247 : !> \param h ...
248 : !> \param k ...
249 : !> \param l ...
250 : !> \param cell ...
251 : !> \return ...
252 : !> \date 18.11.2004
253 : !> \author Matthias Krack
254 : !> \version 1.0
255 : ! **************************************************************************************************
256 7026477 : FUNCTION plane_distance(h, k, l, cell) RESULT(distance)
257 :
258 : INTEGER, INTENT(IN) :: h, k, l
259 : TYPE(cell_type), POINTER :: cell
260 : REAL(KIND=dp) :: distance
261 :
262 : REAL(KIND=dp) :: a, alpha, b, beta, c, cosa, cosb, cosg, &
263 : d, gamma, x, y, z
264 : REAL(KIND=dp), DIMENSION(3) :: abc
265 :
266 7026477 : x = REAL(h, KIND=dp)
267 7026477 : y = REAL(k, KIND=dp)
268 7026477 : z = REAL(l, KIND=dp)
269 :
270 7026477 : CALL get_cell(cell=cell, abc=abc)
271 :
272 7026477 : a = abc(1)
273 7026477 : b = abc(2)
274 7026477 : c = abc(3)
275 :
276 7026477 : IF (cell%orthorhombic) THEN
277 :
278 6873912 : d = (x/a)**2 + (y/b)**2 + (z/c)**2
279 :
280 : ELSE
281 :
282 : CALL get_cell(cell=cell, &
283 : alpha=alpha, &
284 : beta=beta, &
285 152565 : gamma=gamma)
286 :
287 152565 : alpha = alpha/degree
288 152565 : beta = beta/degree
289 152565 : gamma = gamma/degree
290 :
291 152565 : cosa = COS(alpha)
292 152565 : cosb = COS(beta)
293 152565 : cosg = COS(gamma)
294 :
295 : d = ((x*b*c*SIN(alpha))**2 + &
296 : (y*c*a*SIN(beta))**2 + &
297 : (z*a*b*SIN(gamma))**2 + &
298 : 2.0_dp*a*b*c*(x*y*c*(cosa*cosb - cosg) + &
299 : z*x*b*(cosg*cosa - cosb) + &
300 : y*z*a*(cosb*cosg - cosa)))/ &
301 : ((a*b*c)**2*(1.0_dp - cosa**2 - cosb**2 - cosg**2 + &
302 152565 : 2.0_dp*cosa*cosb*cosg))
303 :
304 : END IF
305 :
306 7026477 : distance = 1.0_dp/SQRT(d)
307 :
308 7026477 : END FUNCTION plane_distance
309 :
310 : ! **************************************************************************************************
311 : !> \brief Apply the periodic boundary conditions defined by a simulation
312 : !> cell to a position vector r.
313 : !> \param r ...
314 : !> \param cell ...
315 : !> \return ...
316 : !> \date 16.01.2002
317 : !> \author Matthias Krack
318 : !> \version 1.0
319 : ! **************************************************************************************************
320 384433150 : FUNCTION pbc1(r, cell) RESULT(r_pbc)
321 :
322 : REAL(KIND=dp), DIMENSION(3), INTENT(IN) :: r
323 : TYPE(cell_type), POINTER :: cell
324 : REAL(KIND=dp), DIMENSION(3) :: r_pbc
325 :
326 : REAL(KIND=dp), DIMENSION(3) :: s
327 :
328 384433150 : CPASSERT(ASSOCIATED(cell))
329 :
330 384433150 : IF (cell%orthorhombic) THEN
331 358245433 : r_pbc(1) = r(1) - cell%hmat(1, 1)*cell%perd(1)*ANINT(cell%h_inv(1, 1)*r(1))
332 358245433 : r_pbc(2) = r(2) - cell%hmat(2, 2)*cell%perd(2)*ANINT(cell%h_inv(2, 2)*r(2))
333 358245433 : r_pbc(3) = r(3) - cell%hmat(3, 3)*cell%perd(3)*ANINT(cell%h_inv(3, 3)*r(3))
334 : ELSE
335 26187717 : s(1) = cell%h_inv(1, 1)*r(1) + cell%h_inv(1, 2)*r(2) + cell%h_inv(1, 3)*r(3)
336 26187717 : s(2) = cell%h_inv(2, 1)*r(1) + cell%h_inv(2, 2)*r(2) + cell%h_inv(2, 3)*r(3)
337 26187717 : s(3) = cell%h_inv(3, 1)*r(1) + cell%h_inv(3, 2)*r(2) + cell%h_inv(3, 3)*r(3)
338 26187717 : s(1) = s(1) - cell%perd(1)*ANINT(s(1))
339 26187717 : s(2) = s(2) - cell%perd(2)*ANINT(s(2))
340 26187717 : s(3) = s(3) - cell%perd(3)*ANINT(s(3))
341 26187717 : r_pbc(1) = cell%hmat(1, 1)*s(1) + cell%hmat(1, 2)*s(2) + cell%hmat(1, 3)*s(3)
342 26187717 : r_pbc(2) = cell%hmat(2, 1)*s(1) + cell%hmat(2, 2)*s(2) + cell%hmat(2, 3)*s(3)
343 26187717 : r_pbc(3) = cell%hmat(3, 1)*s(1) + cell%hmat(3, 2)*s(2) + cell%hmat(3, 3)*s(3)
344 : END IF
345 :
346 384433150 : END FUNCTION pbc1
347 :
348 : ! **************************************************************************************************
349 : !> \brief Apply the periodic boundary conditions defined by a simulation
350 : !> cell to a position vector r subtracting nl from the periodic images
351 : !> \param r ...
352 : !> \param cell ...
353 : !> \param nl ...
354 : !> \return ...
355 : !> \date 16.01.2002
356 : !> \author Matthias Krack
357 : !> \version 1.0
358 : ! **************************************************************************************************
359 142966578 : FUNCTION pbc2(r, cell, nl) RESULT(r_pbc)
360 :
361 : REAL(KIND=dp), DIMENSION(3), INTENT(IN) :: r
362 : TYPE(cell_type), POINTER :: cell
363 : INTEGER, DIMENSION(3), INTENT(IN) :: nl
364 : REAL(KIND=dp), DIMENSION(3) :: r_pbc
365 :
366 : REAL(KIND=dp), DIMENSION(3) :: s
367 :
368 142966578 : CPASSERT(ASSOCIATED(cell))
369 :
370 142966578 : IF (cell%orthorhombic) THEN
371 : r_pbc(1) = r(1) - cell%hmat(1, 1)*cell%perd(1)* &
372 132353154 : REAL(NINT(cell%h_inv(1, 1)*r(1)) - nl(1), dp)
373 : r_pbc(2) = r(2) - cell%hmat(2, 2)*cell%perd(2)* &
374 132353154 : REAL(NINT(cell%h_inv(2, 2)*r(2)) - nl(2), dp)
375 : r_pbc(3) = r(3) - cell%hmat(3, 3)*cell%perd(3)* &
376 132353154 : REAL(NINT(cell%h_inv(3, 3)*r(3)) - nl(3), dp)
377 : ELSE
378 10613424 : s(1) = cell%h_inv(1, 1)*r(1) + cell%h_inv(1, 2)*r(2) + cell%h_inv(1, 3)*r(3)
379 10613424 : s(2) = cell%h_inv(2, 1)*r(1) + cell%h_inv(2, 2)*r(2) + cell%h_inv(2, 3)*r(3)
380 10613424 : s(3) = cell%h_inv(3, 1)*r(1) + cell%h_inv(3, 2)*r(2) + cell%h_inv(3, 3)*r(3)
381 10613424 : s(1) = s(1) - cell%perd(1)*REAL(NINT(s(1)) - nl(1), dp)
382 10613424 : s(2) = s(2) - cell%perd(2)*REAL(NINT(s(2)) - nl(2), dp)
383 10613424 : s(3) = s(3) - cell%perd(3)*REAL(NINT(s(3)) - nl(3), dp)
384 10613424 : r_pbc(1) = cell%hmat(1, 1)*s(1) + cell%hmat(1, 2)*s(2) + cell%hmat(1, 3)*s(3)
385 10613424 : r_pbc(2) = cell%hmat(2, 1)*s(1) + cell%hmat(2, 2)*s(2) + cell%hmat(2, 3)*s(3)
386 10613424 : r_pbc(3) = cell%hmat(3, 1)*s(1) + cell%hmat(3, 2)*s(2) + cell%hmat(3, 3)*s(3)
387 : END IF
388 :
389 142966578 : END FUNCTION pbc2
390 :
391 : ! **************************************************************************************************
392 : !> \brief Apply the periodic boundary conditions defined by the simulation
393 : !> cell cell to the vector pointing from atom a to atom b.
394 : !> \param ra ...
395 : !> \param rb ...
396 : !> \param cell ...
397 : !> \return ...
398 : !> \date 11.03.2004
399 : !> \author Matthias Krack
400 : !> \version 1.0
401 : ! **************************************************************************************************
402 125807792 : FUNCTION pbc3(ra, rb, cell) RESULT(rab_pbc)
403 :
404 : REAL(KIND=dp), DIMENSION(3), INTENT(IN) :: ra, rb
405 : TYPE(cell_type), POINTER :: cell
406 : REAL(KIND=dp), DIMENSION(3) :: rab_pbc
407 :
408 : INTEGER :: icell, jcell, kcell
409 : INTEGER, DIMENSION(3) :: periodic
410 : REAL(KIND=dp) :: rab2, rab2_pbc
411 : REAL(KIND=dp), DIMENSION(3) :: r, ra_pbc, rab, rb_image, rb_pbc, s2r
412 :
413 125807792 : CALL get_cell(cell=cell, periodic=periodic)
414 :
415 125807792 : ra_pbc(:) = pbc(ra(:), cell)
416 125807792 : rb_pbc(:) = pbc(rb(:), cell)
417 :
418 125807792 : rab2_pbc = HUGE(1.0_dp)
419 :
420 498855862 : DO icell = -periodic(1), periodic(1)
421 1613626922 : DO jcell = -periodic(2), periodic(2)
422 4827753028 : DO kcell = -periodic(3), periodic(3)
423 13359735592 : r = REAL([icell, jcell, kcell], dp)
424 3339933898 : CALL scaled_to_real(s2r, r, cell)
425 13359735592 : rb_image(:) = rb_pbc(:) + s2r
426 13359735592 : rab(:) = rb_image(:) - ra_pbc(:)
427 3339933898 : rab2 = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
428 4454704958 : IF (rab2 < rab2_pbc) THEN
429 2613591020 : rab2_pbc = rab2
430 2613591020 : rab_pbc(:) = rab(:)
431 : END IF
432 : END DO
433 : END DO
434 : END DO
435 :
436 125807792 : END FUNCTION pbc3
437 :
438 : !if positive_range == true, r(i) (or s(i)) in range [0, hmat(i,i)],
439 : !else, r(i) (s(i)) in range [-hmat(i,i)/2, hmat(i,i)/2]
440 : ! **************************************************************************************************
441 : !> \brief ...
442 : !> \param r ...
443 : !> \param cell ...
444 : !> \param positive_range ...
445 : !> \return ...
446 : ! **************************************************************************************************
447 238 : FUNCTION pbc4(r, cell, positive_range) RESULT(r_pbc)
448 :
449 : REAL(KIND=dp), DIMENSION(3), INTENT(IN) :: r
450 : TYPE(cell_type), POINTER :: cell
451 : LOGICAL :: positive_range
452 : REAL(KIND=dp), DIMENSION(3) :: r_pbc
453 :
454 : REAL(KIND=dp), DIMENSION(3) :: s
455 :
456 238 : CPASSERT(ASSOCIATED(cell))
457 :
458 238 : IF (positive_range) THEN
459 238 : IF (cell%orthorhombic) THEN
460 238 : r_pbc(1) = r(1) - cell%hmat(1, 1)*cell%perd(1)*FLOOR(cell%h_inv(1, 1)*r(1))
461 238 : r_pbc(2) = r(2) - cell%hmat(2, 2)*cell%perd(2)*FLOOR(cell%h_inv(2, 2)*r(2))
462 238 : r_pbc(3) = r(3) - cell%hmat(3, 3)*cell%perd(3)*FLOOR(cell%h_inv(3, 3)*r(3))
463 : ELSE
464 0 : s(1) = cell%h_inv(1, 1)*r(1) + cell%h_inv(1, 2)*r(2) + cell%h_inv(1, 3)*r(3)
465 0 : s(2) = cell%h_inv(2, 1)*r(1) + cell%h_inv(2, 2)*r(2) + cell%h_inv(2, 3)*r(3)
466 0 : s(3) = cell%h_inv(3, 1)*r(1) + cell%h_inv(3, 2)*r(2) + cell%h_inv(3, 3)*r(3)
467 0 : s(1) = s(1) - cell%perd(1)*FLOOR(s(1))
468 0 : s(2) = s(2) - cell%perd(2)*FLOOR(s(2))
469 0 : s(3) = s(3) - cell%perd(3)*FLOOR(s(3))
470 0 : r_pbc(1) = cell%hmat(1, 1)*s(1) + cell%hmat(1, 2)*s(2) + cell%hmat(1, 3)*s(3)
471 0 : r_pbc(2) = cell%hmat(2, 1)*s(1) + cell%hmat(2, 2)*s(2) + cell%hmat(2, 3)*s(3)
472 0 : r_pbc(3) = cell%hmat(3, 1)*s(1) + cell%hmat(3, 2)*s(2) + cell%hmat(3, 3)*s(3)
473 : END IF
474 : ELSE
475 0 : r_pbc = pbc1(r, cell)
476 : END IF
477 :
478 238 : END FUNCTION pbc4
479 :
480 : ! **************************************************************************************************
481 : !> \brief Transform real to scaled cell coordinates.
482 : !> s=h_inv*r
483 : !> \param s ...
484 : !> \param r ...
485 : !> \param cell ...
486 : !> \date 16.01.2002
487 : !> \author Matthias Krack
488 : !> \version 1.0
489 : ! **************************************************************************************************
490 93226106 : SUBROUTINE real_to_scaled(s, r, cell)
491 :
492 : REAL(KIND=dp), DIMENSION(3), INTENT(OUT) :: s
493 : REAL(KIND=dp), DIMENSION(3), INTENT(IN) :: r
494 : TYPE(cell_type), POINTER :: cell
495 :
496 93226106 : CPASSERT(ASSOCIATED(cell))
497 :
498 93226106 : IF (cell%orthorhombic) THEN
499 86429500 : s(1) = cell%h_inv(1, 1)*r(1)
500 86429500 : s(2) = cell%h_inv(2, 2)*r(2)
501 86429500 : s(3) = cell%h_inv(3, 3)*r(3)
502 : ELSE
503 6796606 : s(1) = cell%h_inv(1, 1)*r(1) + cell%h_inv(1, 2)*r(2) + cell%h_inv(1, 3)*r(3)
504 6796606 : s(2) = cell%h_inv(2, 1)*r(1) + cell%h_inv(2, 2)*r(2) + cell%h_inv(2, 3)*r(3)
505 6796606 : s(3) = cell%h_inv(3, 1)*r(1) + cell%h_inv(3, 2)*r(2) + cell%h_inv(3, 3)*r(3)
506 : END IF
507 :
508 93226106 : END SUBROUTINE real_to_scaled
509 :
510 : ! **************************************************************************************************
511 : !> \brief Transform scaled cell coordinates real coordinates.
512 : !> r=h*s
513 : !> \param r ...
514 : !> \param s ...
515 : !> \param cell ...
516 : !> \date 16.01.2002
517 : !> \author Matthias Krack
518 : !> \version 1.0
519 : ! **************************************************************************************************
520 3478055948 : SUBROUTINE scaled_to_real(r, s, cell)
521 :
522 : REAL(KIND=dp), DIMENSION(3), INTENT(OUT) :: r
523 : REAL(KIND=dp), DIMENSION(3), INTENT(IN) :: s
524 : TYPE(cell_type), POINTER :: cell
525 :
526 3478055948 : CPASSERT(ASSOCIATED(cell))
527 :
528 3478055948 : IF (cell%orthorhombic) THEN
529 3201654673 : r(1) = cell%hmat(1, 1)*s(1)
530 3201654673 : r(2) = cell%hmat(2, 2)*s(2)
531 3201654673 : r(3) = cell%hmat(3, 3)*s(3)
532 : ELSE
533 276401275 : r(1) = cell%hmat(1, 1)*s(1) + cell%hmat(1, 2)*s(2) + cell%hmat(1, 3)*s(3)
534 276401275 : r(2) = cell%hmat(2, 1)*s(1) + cell%hmat(2, 2)*s(2) + cell%hmat(2, 3)*s(3)
535 276401275 : r(3) = cell%hmat(3, 1)*s(1) + cell%hmat(3, 2)*s(2) + cell%hmat(3, 3)*s(3)
536 : END IF
537 :
538 3478055948 : END SUBROUTINE scaled_to_real
539 : ! **************************************************************************************************
540 : !> \brief retains the given cell (see doc/ReferenceCounting.html)
541 : !> \param cell the cell to retain
542 : !> \par History
543 : !> 09.2003 created [fawzi]
544 : !> \author Fawzi Mohamed
545 : ! **************************************************************************************************
546 57185 : SUBROUTINE cell_retain(cell)
547 :
548 : TYPE(cell_type), POINTER :: cell
549 :
550 57185 : CPASSERT(ASSOCIATED(cell))
551 57185 : CPASSERT(cell%ref_count > 0)
552 57185 : cell%ref_count = cell%ref_count + 1
553 :
554 57185 : END SUBROUTINE cell_retain
555 :
556 : ! **************************************************************************************************
557 : !> \brief releases the given cell (see doc/ReferenceCounting.html)
558 : !> \param cell the cell to release
559 : !> \par History
560 : !> 09.2003 created [fawzi]
561 : !> \author Fawzi Mohamed
562 : ! **************************************************************************************************
563 152491 : SUBROUTINE cell_release(cell)
564 :
565 : TYPE(cell_type), POINTER :: cell
566 :
567 152491 : IF (ASSOCIATED(cell)) THEN
568 118464 : CPASSERT(cell%ref_count > 0)
569 118464 : cell%ref_count = cell%ref_count - 1
570 118464 : IF (cell%ref_count == 0) THEN
571 61279 : DEALLOCATE (cell)
572 : END IF
573 118464 : NULLIFY (cell)
574 : END IF
575 :
576 152491 : END SUBROUTINE cell_release
577 :
578 : #if defined (__PLUMED2)
579 : ! **************************************************************************************************
580 : !> \brief For the interface with plumed, pass a cell pointer and retrieve it
581 : !> later. It's a hack, but avoids passing the cell back and forth
582 : !> across the Fortran/C++ interface
583 : !> \param cell ...
584 : !> \param set ...
585 : !> \date 28.02.2013
586 : !> \author RK
587 : !> \version 1.0
588 : ! **************************************************************************************************
589 2 : SUBROUTINE pbc_cp2k_plumed_getset_cell(cell, set)
590 :
591 : TYPE(cell_type), POINTER :: cell
592 : LOGICAL :: set
593 :
594 : TYPE(cell_type), POINTER, SAVE :: stored_cell
595 :
596 2 : IF (set) THEN
597 2 : stored_cell => cell
598 : ELSE
599 0 : cell => stored_cell
600 : END IF
601 :
602 2 : END SUBROUTINE pbc_cp2k_plumed_getset_cell
603 : #endif
604 :
605 0 : END MODULE cell_types
|
