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 Contains type used for a Simulation Cell Optimization
10 : !> \par History
11 : !> none
12 : !> \author Teodoro Laino - created [tlaino] - 03.2008 - Zurich University
13 : ! **************************************************************************************************
14 : MODULE cell_opt_types
15 :
16 : USE cell_methods, ONLY: canonicalize_cell_matrix,&
17 : cell_create
18 : USE cell_opt_utils, ONLY: cell_opt_constraint_label,&
19 : keep_angles_conflict,&
20 : keep_symmetry_conflict,&
21 : keep_symmetry_coupled_set,&
22 : mask_shear,&
23 : parse_cell_constraint,&
24 : pres_crit_shear_pinned,&
25 : read_external_press_tensor
26 : USE cell_types, ONLY: cell_clone,&
27 : cell_release,&
28 : cell_type
29 : USE cp_log_handling, ONLY: cp_get_default_logger,&
30 : cp_logger_type
31 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
32 : cp_print_key_unit_nr
33 : USE cp_units, ONLY: cp_unit_from_cp2k
34 : USE force_env_types, ONLY: force_env_get,&
35 : force_env_type
36 : USE input_constants, ONLY: fix_none
37 : USE input_section_types, ONLY: section_vals_type,&
38 : section_vals_val_get,&
39 : section_vals_val_set
40 : USE kinds, ONLY: default_string_length,&
41 : dp
42 : USE string_utilities, ONLY: uppercase
43 : #include "../base/base_uses.f90"
44 :
45 : IMPLICIT NONE
46 : PRIVATE
47 :
48 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'cell_opt_types'
49 :
50 : PUBLIC :: cell_opt_env_type, &
51 : cell_opt_env_create, &
52 : cell_opt_env_release
53 :
54 : ! **************************************************************************************************
55 : !> \brief Type containing all informations abour the simulation cell optimization
56 : !> \par History
57 : !> none
58 : !> \author Teodoro Laino - created [tlaino] - 03.2008 - Zurich University
59 : ! **************************************************************************************************
60 : TYPE cell_opt_env_type
61 : ! Simulation cell optimization parameters
62 : ! constraint_id is a mask of fixed cell components (comp_* / mask_*):
63 : ! a set bit means fixed, 0 (= fix_none) leaves all components relaxable
64 : INTEGER :: constraint_id = fix_none
65 : LOGICAL :: keep_angles = .FALSE., &
66 : keep_symmetry = .FALSE., &
67 : keep_volume = .FALSE.
68 : REAL(KIND=dp) :: pres_ext = 0.0_dp, pres_int = 0.0_dp, pres_tol = 0.0_dp, pres_constr = 0.0_dp
69 : REAL(KIND=dp), DIMENSION(3, 3) :: mtrx = 0.0_dp
70 : ! Fixed rotations between the user's Cartesian frame and the canonical
71 : ! frame used for the six independent cell optimization variables.
72 : REAL(KIND=dp), DIMENSION(3, 3) :: input_to_opt = 0.0_dp, &
73 : opt_to_input = 0.0_dp
74 : TYPE(cell_type), POINTER :: opt_cell => NULL()
75 : END TYPE cell_opt_env_type
76 :
77 : CONTAINS
78 :
79 : ! **************************************************************************************************
80 : !> \brief ...
81 : !> \param cell_env ...
82 : !> \param force_env ...
83 : !> \param geo_section ...
84 : !> \par History
85 : !> none
86 : !> \author Teodoro Laino - created [tlaino] - 03.2008 - Zurich University
87 : ! **************************************************************************************************
88 8658 : SUBROUTINE cell_opt_env_create(cell_env, force_env, geo_section)
89 : TYPE(cell_opt_env_type), INTENT(OUT) :: cell_env
90 : TYPE(force_env_type), POINTER :: force_env
91 : TYPE(section_vals_type), POINTER :: geo_section
92 :
93 : CHARACTER(LEN=default_string_length) :: label
94 : CHARACTER(LEN=default_string_length), &
95 234 : DIMENSION(:), POINTER :: canonical_tokens, constraint_tokens
96 : INTEGER :: output_unit
97 : LOGICAL :: constraint_explicit
98 : TYPE(cell_type), POINTER :: cell
99 : TYPE(cp_logger_type), POINTER :: logger
100 :
101 234 : NULLIFY (cell_env%opt_cell, cell, canonical_tokens, constraint_tokens)
102 234 : CALL force_env_get(force_env, cell=cell)
103 234 : CALL cell_create(cell_env%opt_cell)
104 234 : CALL cell_clone(cell, cell_env%opt_cell, tag="CELL_OPT_CANONICAL")
105 234 : CALL canonicalize_cell_matrix(cell_env%opt_cell)
106 : ! Keep the physical cell in the input orientation while the optimizer
107 : ! works with the canonical cell and its six non-redundant components.
108 12168 : cell_env%input_to_opt = MATMUL(cell_env%opt_cell%hmat, cell%h_inv)
109 12168 : cell_env%opt_to_input = MATMUL(cell%hmat, cell_env%opt_cell%h_inv)
110 234 : CALL section_vals_val_get(geo_section, "KEEP_VOLUME", l_val=cell_env%keep_volume)
111 234 : CALL section_vals_val_get(geo_section, "KEEP_ANGLES", l_val=cell_env%keep_angles)
112 234 : CALL section_vals_val_get(geo_section, "KEEP_SYMMETRY", l_val=cell_env%keep_symmetry)
113 234 : CALL section_vals_val_get(geo_section, "PRESSURE_TOLERANCE", r_val=cell_env%pres_tol)
114 : CALL section_vals_val_get(geo_section, "CONSTRAINT", &
115 234 : c_vals=constraint_tokens, explicit=constraint_explicit)
116 : CALL parse_cell_constraint(constraint_tokens, cell, constraint_explicit, &
117 234 : cell_env%constraint_id)
118 : ! Only for explicitly set values: an unset keyword must stay unparsed so
119 : ! the restart echo keeps omitting the line, as it always has.
120 234 : IF (constraint_explicit) THEN
121 102 : ALLOCATE (canonical_tokens(SIZE(constraint_tokens)))
122 114 : canonical_tokens = constraint_tokens
123 74 : CALL uppercase(canonical_tokens)
124 34 : CALL section_vals_val_set(geo_section, "CONSTRAINT", c_vals_ptr=canonical_tokens)
125 : END IF
126 :
127 234 : CALL cell_opt_validate_constraint(cell_env)
128 :
129 : ! The masked pressure criterion is structurally off target for this
130 : ! combination on a non-orthorhombic cell. The run stays legal: the
131 : ! geometry criteria still work, and the pressure criterion may still
132 : ! be met when the pinned shear stresses vanish by symmetry. The
133 : ! decision lives in pres_crit_shear_pinned, which the unit test
134 : ! checks over all 64 masks.
135 234 : IF (pres_crit_shear_pinned(cell_env%constraint_id, &
136 : cell_env%opt_cell%orthorhombic, &
137 : cell_env%keep_angles)) THEN
138 4 : IF (IAND(cell_env%constraint_id, mask_shear) /= 0) THEN
139 : CALL cp_warn(__LOCATION__, &
140 : "CELL_OPT: CONSTRAINT ("// &
141 : TRIM(cell_opt_constraint_label(cell_env%constraint_id))// &
142 : ") pins a shear component of a non-orthorhombic cell. "// &
143 : "The masked pressure criterion is met at the constrained minimum "// &
144 : "only if the stresses conjugate to the pinned components vanish "// &
145 : "there, for example when the pinned values are the stress-free "// &
146 : "ones; otherwise it can stay off the external pressure. Check "// &
147 : "the per-step ""Pressure is converged"" line: with L-BFGS the "// &
148 : """run CONVERGED"" banner can appear while it is NO. The "// &
149 2 : "projected-stress criterion planned in #6031 removes this limitation.")
150 : ELSE
151 : CALL cp_warn(__LOCATION__, &
152 : "CELL_OPT: CONSTRAINT ("// &
153 : TRIM(cell_opt_constraint_label(cell_env%constraint_id))// &
154 : ") fixes a diagonal component while KEEP_ANGLES pins the "// &
155 : "angles of a non-orthorhombic cell. The masked pressure criterion "// &
156 : "is met at the constrained minimum only if the stresses conjugate "// &
157 : "to the pinned components vanish there, for example when the "// &
158 : "pinned values are the stress-free ones; otherwise it can stay "// &
159 : "off the external pressure. Check the per-step ""Pressure is "// &
160 : "converged"" line: with L-BFGS the ""run CONVERGED"" banner can "// &
161 : "appear while it is NO. The projected-stress criterion planned "// &
162 2 : "in #6031 removes this limitation.")
163 : END IF
164 : END IF
165 :
166 : ! Get the external pressure
167 : CALL read_external_press_tensor(geo_section, cell_env%opt_cell, cell_env%pres_ext, &
168 234 : cell_env%mtrx, cell_env%input_to_opt)
169 :
170 : ! Print cell optimisation setup
171 234 : NULLIFY (logger)
172 234 : logger => cp_get_default_logger()
173 234 : output_unit = cp_print_key_unit_nr(logger, geo_section, "PRINT%CELL", extension=".Log")
174 234 : IF (output_unit > 0) THEN
175 : WRITE (UNIT=output_unit, FMT="(/,T2,A,T61,F20.1)") &
176 117 : "CELL_OPT| Pressure tolerance [bar]: ", cp_unit_from_cp2k(cell_env%pres_tol, "bar")
177 117 : IF (cell_env%keep_volume) THEN
178 : WRITE (UNIT=output_unit, FMT="(T2,A,T78,A3)") &
179 3 : "CELL_OPT| Keep volume of cell: ", "YES"
180 : ELSE
181 : WRITE (UNIT=output_unit, FMT="(T2,A,T78,A3)") &
182 114 : "CELL_OPT| Keep volume of cell: ", " NO"
183 : END IF
184 117 : IF (cell_env%keep_angles) THEN
185 : WRITE (UNIT=output_unit, FMT="(T2,A,T78,A3)") &
186 21 : "CELL_OPT| Keep angles between the cell vectors: ", "YES"
187 : ELSE
188 : WRITE (UNIT=output_unit, FMT="(T2,A,T78,A3)") &
189 96 : "CELL_OPT| Keep angles between the cell vectors: ", " NO"
190 : END IF
191 117 : IF (cell_env%keep_symmetry) THEN
192 : WRITE (UNIT=output_unit, FMT="(T2,A,T78,A3)") &
193 21 : "CELL_OPT| Keep cell symmetry: ", "YES"
194 : ELSE
195 : WRITE (UNIT=output_unit, FMT="(T2,A,T78,A3)") &
196 96 : "CELL_OPT| Keep cell symmetry: ", " NO"
197 : END IF
198 117 : label = cell_opt_constraint_label(cell_env%constraint_id)
199 468 : IF (COUNT(cell%perd /= 0) == 2 .AND. .NOT. constraint_explicit) THEN
200 1 : label = TRIM(label)//" (auto)"
201 : END IF
202 : WRITE (UNIT=output_unit, FMT="(T2,A,T77,A)") &
203 117 : "CELL_OPT| Constraint: ", TRIM(label)
204 : END IF
205 234 : CALL cp_print_key_finished_output(output_unit, logger, geo_section, "PRINT%CELL")
206 :
207 468 : END SUBROUTINE cell_opt_env_create
208 :
209 : ! **************************************************************************************************
210 : !> \brief Abort on constraint combinations that corrupt the KEEP_* projections
211 : !> \param cell_env ...
212 : ! **************************************************************************************************
213 234 : SUBROUTINE cell_opt_validate_constraint(cell_env)
214 : TYPE(cell_opt_env_type), INTENT(IN) :: cell_env
215 :
216 : LOGICAL :: bad_constraint
217 :
218 : ! Reject CONSTRAINT combinations that would silently corrupt the
219 : ! KEEP_ANGLES projection or the KEEP_SYMMETRY symmetrization.
220 : ! Slot order of apply_cell_constraints: AX, BX, BY, CX, CY, CZ.
221 : ! KEEP_ANGLES projects the gradient slots {BX,BY} onto b and {CX,CY,CZ} onto c.
222 : ! Gradient slots that a legacy value zeroes:
223 : ! x = {AX,BX,CX} B partial, C partial -> breaks parallelism
224 : ! y = {BX,BY,CY} B complete, C partial -> breaks parallelism
225 : ! z = {CX,CY,CZ} B untouched, C complete -> safe
226 : ! xy = {AX..CY} B complete, C partial -> breaks parallelism
227 : ! xz = {AX,BX,CX,CY,CZ} B partial, C complete -> breaks parallelism
228 : ! yz = {BX..CZ} B complete, C complete -> safe
229 : ! A partial group breaks the parallelism. A complete or untouched group keeps it.
230 234 : IF (cell_env%keep_angles) THEN
231 : ! The decision lives in keep_angles_conflict, which the unit test
232 : ! checks over all 64 masks.
233 42 : IF (keep_angles_conflict(cell_env%constraint_id, &
234 : cell_env%opt_cell%orthorhombic)) THEN
235 : CALL cp_abort(__LOCATION__, &
236 : "CELL_OPT: CONSTRAINT ("// &
237 : TRIM(cell_opt_constraint_label(cell_env%constraint_id))// &
238 : ") fixes only part of the b or c cell vector and would "// &
239 0 : "destroy the KEEP_ANGLES projection.")
240 : END IF
241 : END IF
242 :
243 234 : IF (cell_env%keep_symmetry .AND. cell_env%constraint_id /= fix_none) THEN
244 : ! The decision lives in keep_symmetry_conflict, which the unit test
245 : ! checks over all 64 masks for every crystal system.
246 : bad_constraint = keep_symmetry_conflict(cell_env%constraint_id, &
247 0 : cell_env%opt_cell%symmetry_id)
248 0 : IF (bad_constraint) THEN
249 : CALL cp_abort(__LOCATION__, &
250 : "CELL_OPT: CONSTRAINT ("// &
251 : TRIM(cell_opt_constraint_label(cell_env%constraint_id))// &
252 : ") fixes only part of the components coupled by "// &
253 : "KEEP_SYMMETRY ("// &
254 : TRIM(cell_opt_constraint_label( &
255 : keep_symmetry_coupled_set(cell_env%opt_cell%symmetry_id)))// &
256 0 : ") and would break the cell symmetry.")
257 : END IF
258 : END IF
259 :
260 234 : END SUBROUTINE cell_opt_validate_constraint
261 :
262 : ! **************************************************************************************************
263 : !> \brief ...
264 : !> \param cell_env ...
265 : !> \par History
266 : !> none
267 : !> \author Teodoro Laino - created [tlaino] - 03.2008 - Zurich University
268 : ! **************************************************************************************************
269 234 : SUBROUTINE cell_opt_env_release(cell_env)
270 : TYPE(cell_opt_env_type), INTENT(INOUT) :: cell_env
271 :
272 234 : CALL cell_release(cell_env%opt_cell)
273 :
274 234 : END SUBROUTINE cell_opt_env_release
275 :
276 0 : END MODULE cell_opt_types
|
