Line data Source code
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 : !> \par History
10 : !> Subroutine input_torsions changed (DG) 05-Dec-2000
11 : !> Output formats changed (DG) 05-Dec-2000
12 : !> JGH (26-01-2002) : force field parameters stored in tables, not in
13 : !> matrices. Input changed to have parameters labeled by the position
14 : !> and not atom pairs (triples etc)
15 : !> Teo (11.2005) : Moved all information on force field pair_potential to
16 : !> a much lighter memory structure
17 : !> Teo 09.2006 : Split all routines force_field I/O in a separate file
18 : !> \author CJM
19 : ! **************************************************************************************************
20 : MODULE force_fields_input
21 : USE bibliography, ONLY: Clabaut2020,&
22 : Clabaut2021,&
23 : Siepmann1995,&
24 : Tersoff1988,&
25 : Tosi1964a,&
26 : Tosi1964b,&
27 : Yamada2000,&
28 : cite_reference
29 : USE cp_files, ONLY: discover_file
30 : USE cp_linked_list_input, ONLY: cp_sll_val_next,&
31 : cp_sll_val_type
32 : USE cp_log_handling, ONLY: cp_get_default_logger,&
33 : cp_logger_type,&
34 : cp_to_string
35 : USE cp_output_handling, ONLY: cp_print_key_finished_output,&
36 : cp_print_key_unit_nr
37 : USE cp_parser_methods, ONLY: parser_get_next_line
38 : USE cp_parser_types, ONLY: cp_parser_type,&
39 : parser_create,&
40 : parser_release
41 : USE cp_units, ONLY: cp_unit_to_cp2k
42 : USE damping_dipole_types, ONLY: damping_info_type
43 : USE force_field_kind_types, ONLY: do_ff_amber,&
44 : do_ff_charmm,&
45 : do_ff_g87,&
46 : do_ff_g96,&
47 : do_ff_opls,&
48 : do_ff_undef,&
49 : legendre_data_type
50 : USE force_field_types, ONLY: force_field_type,&
51 : input_info_type
52 : USE force_fields_util, ONLY: get_generic_info
53 : USE input_section_types, ONLY: section_vals_get,&
54 : section_vals_get_subs_vals,&
55 : section_vals_list_get,&
56 : section_vals_type,&
57 : section_vals_val_get
58 : USE input_val_types, ONLY: val_get,&
59 : val_type
60 : USE kinds, ONLY: default_path_length,&
61 : default_string_length,&
62 : dp
63 : USE mathconstants, ONLY: pi
64 : USE mathlib, ONLY: invert_matrix
65 : USE memory_utilities, ONLY: reallocate
66 : USE message_passing, ONLY: mp_para_env_type
67 : USE pair_potential_types, ONLY: &
68 : allegro_pot_type, allegro_type, b4_type, bm_type, deepmd_type, &
69 : do_potential_single_allocation, ea_type, eam_pot_type, ft_pot_type, ft_type, ftd_type, &
70 : gal21_type, gal_type, gp_type, gw_type, ip_type, ipbv_pot_type, lj_charmm_type, &
71 : nequip_pot_type, nequip_type, no_potential_single_allocation, pair_potential_p_type, &
72 : pair_potential_reallocate, potential_single_allocation, quip_type, siepmann_type, &
73 : tab_pot_type, tab_type, tersoff_type, wl_type
74 : USE shell_potential_types, ONLY: shell_p_create,&
75 : shell_p_type
76 : USE string_utilities, ONLY: uppercase
77 : USE torch_api, ONLY: torch_allow_tf32,&
78 : torch_model_read_metadata
79 : #include "./base/base_uses.f90"
80 :
81 : IMPLICIT NONE
82 :
83 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'force_fields_input'
84 :
85 : PRIVATE
86 : PUBLIC :: read_force_field_section, &
87 : read_lj_section, &
88 : read_wl_section, &
89 : read_gd_section, &
90 : read_gp_section, &
91 : read_chrg_section
92 :
93 : CONTAINS
94 :
95 : ! **************************************************************************************************
96 : !> \brief Reads the force_field input section
97 : !> \param ff_section ...
98 : !> \param mm_section ...
99 : !> \param ff_type ...
100 : !> \param para_env ...
101 : !> \author teo
102 : ! **************************************************************************************************
103 36946 : SUBROUTINE read_force_field_section1(ff_section, mm_section, ff_type, para_env)
104 : TYPE(section_vals_type), POINTER :: ff_section, mm_section
105 : TYPE(force_field_type), INTENT(INOUT) :: ff_type
106 : TYPE(mp_para_env_type), POINTER :: para_env
107 :
108 : CHARACTER(LEN=default_string_length), &
109 2639 : DIMENSION(:), POINTER :: atm_names
110 : INTEGER :: nallegro, nb4, nbends, nbm, nbmhft, nbmhftd, nbonds, nchg, ndeepmd, neam, ngal, &
111 : ngal21, ngd, ngp, nimpr, nipbv, nlj, nnequip, nopbend, nquip, nshell, nsiepmann, ntab, &
112 : ntersoff, ntors, ntot, nubs, nwl
113 : LOGICAL :: explicit, unique_spline
114 : REAL(KIND=dp) :: min_eps_spline_allowed
115 : TYPE(input_info_type), POINTER :: inp_info
116 : TYPE(section_vals_type), POINTER :: tmp_section, tmp_section2
117 :
118 : INTEGER::i
119 :
120 2639 : NULLIFY (tmp_section, tmp_section2)
121 2639 : inp_info => ff_type%inp_info
122 2639 : CALL section_vals_val_get(ff_section, "PARMTYPE", i_val=ff_type%ff_type)
123 2639 : CALL section_vals_val_get(ff_section, "EI_SCALE14", r_val=ff_type%ei_scale14)
124 2639 : CALL section_vals_val_get(ff_section, "VDW_SCALE14", r_val=ff_type%vdw_scale14)
125 2639 : CALL section_vals_val_get(ff_section, "SPLINE%RCUT_NB", r_val=ff_type%rcut_nb)
126 2639 : CALL section_vals_val_get(ff_section, "SPLINE%R0_NB", r_val=ff_type%rlow_nb)
127 2639 : CALL section_vals_val_get(ff_section, "SPLINE%EPS_SPLINE", r_val=ff_type%eps_spline)
128 2639 : CALL section_vals_val_get(ff_section, "SPLINE%EMAX_SPLINE", r_val=ff_type%emax_spline)
129 2639 : CALL section_vals_val_get(ff_section, "SPLINE%EMAX_ACCURACY", r_val=ff_type%max_energy)
130 2639 : CALL section_vals_val_get(ff_section, "SPLINE%NPOINTS", i_val=ff_type%npoints)
131 2639 : CALL section_vals_val_get(ff_section, "IGNORE_MISSING_CRITICAL_PARAMS", l_val=ff_type%ignore_missing_critical)
132 2639 : CPASSERT(ff_type%max_energy <= ff_type%emax_spline)
133 : ! Read the parameter file name only if the force field type requires it..
134 3543 : SELECT CASE (ff_type%ff_type)
135 : CASE (do_ff_charmm, do_ff_amber, do_ff_g96, do_ff_g87)
136 904 : CALL section_vals_val_get(ff_section, "PARM_FILE_NAME", c_val=ff_type%ff_file_name)
137 904 : IF (TRIM(ff_type%ff_file_name) == "") &
138 0 : CPABORT("Force Field Parameter's filename is empty! Please check your input file.")
139 : CASE (do_ff_undef)
140 : ! Do Nothing
141 : CASE DEFAULT
142 2639 : CPABORT("Force field type not implemented")
143 : END SELECT
144 : ! Numerical Accuracy:
145 : ! the factors here should depend on the energy and on the shape of each potential mapped
146 : ! with splines. this would make everything un-necessarily complicated. Let's just be safe
147 : ! and assume that we cannot achieve an accuracy on the spline 2 orders of magnitude more
148 : ! than the smallest representable number (taking into account also the max_energy for the
149 : ! spline generation
150 2639 : min_eps_spline_allowed = 20.0_dp*MAX(ff_type%max_energy, 10.0_dp)*EPSILON(0.0_dp)
151 2639 : IF (ff_type%eps_spline < min_eps_spline_allowed) THEN
152 : CALL cp_warn(__LOCATION__, &
153 : "Requested spline accuracy ("//TRIM(cp_to_string(ff_type%eps_spline))//" ) "// &
154 : "is smaller than the minimum value allowed ("//TRIM(cp_to_string(min_eps_spline_allowed))// &
155 : " ) with the present machine precision ("//TRIM(cp_to_string(EPSILON(0.0_dp)))//" ). "// &
156 0 : "New EPS_SPLINE value ("//TRIM(cp_to_string(min_eps_spline_allowed))//" ). ")
157 0 : ff_type%eps_spline = min_eps_spline_allowed
158 : END IF
159 2639 : CALL section_vals_val_get(ff_section, "SHIFT_CUTOFF", l_val=ff_type%shift_cutoff)
160 2639 : CALL section_vals_val_get(ff_section, "SPLINE%UNIQUE_SPLINE", l_val=unique_spline)
161 : ! Single spline
162 2639 : potential_single_allocation = no_potential_single_allocation
163 2639 : IF (unique_spline) potential_single_allocation = do_potential_single_allocation
164 :
165 2639 : CALL section_vals_val_get(ff_section, "MULTIPLE_POTENTIAL", l_val=ff_type%multiple_potential)
166 2639 : CALL section_vals_val_get(ff_section, "DO_NONBONDED", l_val=ff_type%do_nonbonded)
167 2639 : CALL section_vals_val_get(ff_section, "DO_ELECTROSTATICS", l_val=ff_type%do_electrostatics)
168 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "NONBONDED")
169 2639 : CALL section_vals_get(tmp_section, explicit=explicit)
170 2639 : IF (explicit .AND. ff_type%do_nonbonded) THEN
171 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "LENNARD-JONES")
172 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nlj)
173 1735 : ntot = 0
174 1735 : IF (explicit) THEN
175 976 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nlj, lj_charmm=.TRUE.)
176 976 : CALL read_lj_section(inp_info%nonbonded, tmp_section2, ntot)
177 : END IF
178 :
179 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "WILLIAMS")
180 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nwl)
181 1735 : ntot = nlj
182 1735 : IF (explicit) THEN
183 373 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nwl, williams=.TRUE.)
184 373 : CALL read_wl_section(inp_info%nonbonded, tmp_section2, ntot)
185 : END IF
186 :
187 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "EAM")
188 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=neam)
189 1735 : ntot = nlj + nwl
190 1735 : IF (explicit) THEN
191 12 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + neam, eam=.TRUE.)
192 12 : CALL read_eam_section(inp_info%nonbonded, tmp_section2, ntot, para_env, mm_section)
193 : END IF
194 :
195 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "GOODWIN")
196 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ngd)
197 1735 : ntot = nlj + nwl + neam
198 1735 : IF (explicit) THEN
199 0 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ngd, goodwin=.TRUE.)
200 0 : CALL read_gd_section(inp_info%nonbonded, tmp_section2, ntot)
201 : END IF
202 :
203 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "IPBV")
204 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nipbv)
205 1735 : ntot = nlj + nwl + neam + ngd
206 1735 : IF (explicit) THEN
207 16 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nipbv, ipbv=.TRUE.)
208 16 : CALL read_ipbv_section(inp_info%nonbonded, tmp_section2, ntot)
209 : END IF
210 :
211 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "BMHFT")
212 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nbmhft)
213 1735 : ntot = nlj + nwl + neam + ngd + nipbv
214 1735 : IF (explicit) THEN
215 4 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nbmhft, bmhft=.TRUE.)
216 4 : CALL read_bmhft_section(inp_info%nonbonded, tmp_section2, ntot)
217 : END IF
218 :
219 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "BMHFTD")
220 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nbmhftd)
221 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft
222 1735 : IF (explicit) THEN
223 18 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nbmhftd, bmhftd=.TRUE.)
224 18 : CALL read_bmhftd_section(inp_info%nonbonded, tmp_section2, ntot)
225 : END IF
226 :
227 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "BUCK4RANGES")
228 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nb4)
229 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd
230 1735 : IF (explicit) THEN
231 262 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nb4, buck4r=.TRUE.)
232 262 : CALL read_b4_section(inp_info%nonbonded, tmp_section2, ntot)
233 : END IF
234 :
235 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "BUCKMORSE")
236 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nbm)
237 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4
238 1735 : IF (explicit) THEN
239 8 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nbm, buckmo=.TRUE.)
240 8 : CALL read_bm_section(inp_info%nonbonded, tmp_section2, ntot)
241 : END IF
242 :
243 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "GENPOT")
244 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ngp)
245 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm
246 1735 : IF (explicit) THEN
247 308 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ngp, gp=.TRUE.)
248 308 : CALL read_gp_section(inp_info%nonbonded, tmp_section2, ntot)
249 : END IF
250 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "TERSOFF")
251 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ntersoff)
252 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp
253 1735 : IF (explicit) THEN
254 40 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ntersoff, tersoff=.TRUE.)
255 40 : CALL read_tersoff_section(inp_info%nonbonded, tmp_section2, ntot, tmp_section2)
256 : END IF
257 :
258 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "GAL19")
259 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ngal)
260 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff
261 1735 : IF (explicit) THEN
262 1 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ngal, gal=.TRUE.)
263 1 : CALL read_gal_section(inp_info%nonbonded, tmp_section2, ntot, tmp_section2)
264 : END IF
265 :
266 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "GAL21")
267 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ngal21)
268 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal
269 1735 : IF (explicit) THEN
270 1 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ngal21, gal21=.TRUE.)
271 1 : CALL read_gal21_section(inp_info%nonbonded, tmp_section2, ntot, tmp_section2)
272 : END IF
273 :
274 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "SIEPMANN")
275 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nsiepmann)
276 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal + ngal21
277 1735 : IF (explicit) THEN
278 5 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nsiepmann, siepmann=.TRUE.)
279 5 : CALL read_siepmann_section(inp_info%nonbonded, tmp_section2, ntot, tmp_section2)
280 : END IF
281 :
282 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "quip")
283 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nquip)
284 1735 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal + ngal21 + nsiepmann
285 1735 : IF (explicit) THEN
286 2 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nquip, quip=.TRUE.)
287 2 : CALL read_quip_section(inp_info%nonbonded, tmp_section2, ntot)
288 : END IF
289 :
290 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "nequip")
291 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nnequip)
292 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal + ngal21 + nsiepmann + &
293 1735 : nquip
294 1735 : IF (explicit) THEN
295 : ! avoid repeating the nequip section for each pair
296 4 : CALL section_vals_val_get(tmp_section2, "ATOMS", c_vals=atm_names)
297 4 : nnequip = nnequip - 1 + SIZE(atm_names) + (SIZE(atm_names)*SIZE(atm_names) - SIZE(atm_names))/2
298 4 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nnequip, nequip=.TRUE.)
299 4 : CALL read_nequip_section(inp_info%nonbonded, tmp_section2, ntot)
300 : END IF
301 :
302 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "allegro")
303 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nallegro)
304 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal + ngal21 + nsiepmann + &
305 1735 : nquip + nnequip
306 1735 : IF (explicit) THEN
307 : ! avoid repeating the allegro section for each pair
308 4 : CALL section_vals_val_get(tmp_section2, "ATOMS", c_vals=atm_names)
309 4 : nallegro = nallegro - 1 + SIZE(atm_names) + (SIZE(atm_names)*SIZE(atm_names) - SIZE(atm_names))/2
310 4 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + nallegro, allegro=.TRUE.)
311 4 : CALL read_allegro_section(inp_info%nonbonded, tmp_section2, ntot)
312 : END IF
313 :
314 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "TABPOT")
315 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ntab)
316 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal + ngal21 + nsiepmann + &
317 1735 : nquip + nnequip + nallegro
318 1735 : IF (explicit) THEN
319 8 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ntab, tab=.TRUE.)
320 8 : CALL read_tabpot_section(inp_info%nonbonded, tmp_section2, ntot, para_env, mm_section)
321 : END IF
322 :
323 1735 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "DEEPMD")
324 1735 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ndeepmd)
325 : ntot = nlj + nwl + neam + ngd + nipbv + nbmhft + nbmhftd + nb4 + nbm + ngp + ntersoff + ngal + ngal21 + nsiepmann + &
326 1735 : nquip + nnequip + nallegro + ntab
327 1735 : IF (explicit) THEN
328 : ! avoid repeating the deepmd section for each pair
329 2 : CALL section_vals_val_get(tmp_section2, "ATOMS", c_vals=atm_names)
330 2 : ndeepmd = ndeepmd - 1 + SIZE(atm_names) + (SIZE(atm_names)*SIZE(atm_names) - SIZE(atm_names))/2
331 2 : CALL pair_potential_reallocate(inp_info%nonbonded, 1, ntot + ndeepmd, deepmd=.TRUE.)
332 2 : CALL read_deepmd_section(inp_info%nonbonded, tmp_section2, ntot)
333 : END IF
334 :
335 : END IF
336 :
337 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "NONBONDED14")
338 2639 : CALL section_vals_get(tmp_section, explicit=explicit)
339 2639 : IF (explicit .AND. ff_type%do_nonbonded) THEN
340 274 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "LENNARD-JONES")
341 274 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nlj)
342 274 : ntot = 0
343 274 : IF (explicit) THEN
344 12 : CALL pair_potential_reallocate(inp_info%nonbonded14, 1, ntot + nlj, lj_charmm=.TRUE.)
345 12 : CALL read_lj_section(inp_info%nonbonded14, tmp_section2, ntot)
346 : END IF
347 274 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "WILLIAMS")
348 274 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=nwl)
349 274 : ntot = nlj
350 274 : IF (explicit) THEN
351 0 : CALL pair_potential_reallocate(inp_info%nonbonded14, 1, ntot + nwl, williams=.TRUE.)
352 0 : CALL read_wl_section(inp_info%nonbonded14, tmp_section2, ntot)
353 : END IF
354 274 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "GOODWIN")
355 274 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ngd)
356 274 : ntot = nlj + nwl
357 274 : IF (explicit) THEN
358 0 : CALL pair_potential_reallocate(inp_info%nonbonded14, 1, ntot + ngd, goodwin=.TRUE.)
359 0 : CALL read_gd_section(inp_info%nonbonded14, tmp_section2, ntot)
360 : END IF
361 274 : tmp_section2 => section_vals_get_subs_vals(tmp_section, "GENPOT")
362 274 : CALL section_vals_get(tmp_section2, explicit=explicit, n_repetition=ngp)
363 274 : ntot = nlj + nwl + ngd
364 274 : IF (explicit) THEN
365 262 : CALL pair_potential_reallocate(inp_info%nonbonded14, 1, ntot + ngp, gp=.TRUE.)
366 262 : CALL read_gp_section(inp_info%nonbonded14, tmp_section2, ntot)
367 : END IF
368 : END IF
369 :
370 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "CHARGE")
371 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nchg)
372 2639 : IF (explicit) THEN
373 2079 : ntot = 0
374 2079 : CALL reallocate(inp_info%charge_atm, 1, nchg)
375 2079 : CALL reallocate(inp_info%charge, 1, nchg)
376 2079 : CALL read_chrg_section(inp_info%charge_atm, inp_info%charge, tmp_section, ntot)
377 : END IF
378 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "DIPOLE")
379 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nchg)
380 2639 : IF (explicit) THEN
381 34 : ntot = 0
382 34 : CALL reallocate(inp_info%apol_atm, 1, nchg)
383 34 : CALL reallocate(inp_info%apol, 1, nchg)
384 : CALL read_apol_section(inp_info%apol_atm, inp_info%apol, inp_info%damping_list, &
385 34 : tmp_section, ntot)
386 : END IF
387 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "QUADRUPOLE")
388 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nchg)
389 2639 : IF (explicit) THEN
390 0 : ntot = 0
391 0 : CALL reallocate(inp_info%cpol_atm, 1, nchg)
392 0 : CALL reallocate(inp_info%cpol, 1, nchg)
393 0 : CALL read_cpol_section(inp_info%cpol_atm, inp_info%cpol, tmp_section, ntot)
394 : END IF
395 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "SHELL")
396 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nshell)
397 2639 : IF (explicit) THEN
398 258 : ntot = 0
399 258 : CALL shell_p_create(inp_info%shell_list, nshell)
400 258 : CALL read_shell_section(inp_info%shell_list, tmp_section, ntot)
401 : END IF
402 :
403 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "BOND")
404 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nbonds)
405 2639 : IF (explicit) THEN
406 975 : ntot = 0
407 975 : CALL reallocate(inp_info%bond_kind, 1, nbonds)
408 975 : CALL reallocate(inp_info%bond_a, 1, nbonds)
409 975 : CALL reallocate(inp_info%bond_b, 1, nbonds)
410 975 : CALL reallocate(inp_info%bond_k, 1, 3, 1, nbonds)
411 975 : CALL reallocate(inp_info%bond_r0, 1, nbonds)
412 975 : CALL reallocate(inp_info%bond_cs, 1, nbonds)
413 : CALL read_bonds_section(inp_info%bond_kind, inp_info%bond_a, inp_info%bond_b, inp_info%bond_k, &
414 975 : inp_info%bond_r0, inp_info%bond_cs, tmp_section, ntot)
415 : END IF
416 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "BEND")
417 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nbends)
418 2639 : IF (explicit) THEN
419 939 : ntot = 0
420 939 : CALL reallocate(inp_info%bend_kind, 1, nbends)
421 939 : CALL reallocate(inp_info%bend_a, 1, nbends)
422 939 : CALL reallocate(inp_info%bend_b, 1, nbends)
423 939 : CALL reallocate(inp_info%bend_c, 1, nbends)
424 939 : CALL reallocate(inp_info%bend_k, 1, nbends)
425 939 : CALL reallocate(inp_info%bend_theta0, 1, nbends)
426 939 : CALL reallocate(inp_info%bend_cb, 1, nbends)
427 939 : CALL reallocate(inp_info%bend_r012, 1, nbends)
428 939 : CALL reallocate(inp_info%bend_r032, 1, nbends)
429 939 : CALL reallocate(inp_info%bend_kbs12, 1, nbends)
430 939 : CALL reallocate(inp_info%bend_kbs32, 1, nbends)
431 939 : CALL reallocate(inp_info%bend_kss, 1, nbends)
432 939 : IF (ASSOCIATED(inp_info%bend_legendre)) THEN
433 0 : DO i = 1, SIZE(inp_info%bend_legendre)
434 0 : IF (ASSOCIATED(inp_info%bend_legendre(i)%coeffs)) THEN
435 0 : DEALLOCATE (inp_info%bend_legendre(i)%coeffs)
436 0 : NULLIFY (inp_info%bend_legendre(i)%coeffs)
437 : END IF
438 : END DO
439 0 : DEALLOCATE (inp_info%bend_legendre)
440 : NULLIFY (inp_info%bend_legendre)
441 : END IF
442 4936 : ALLOCATE (inp_info%bend_legendre(1:nbends))
443 3058 : DO i = 1, SIZE(inp_info%bend_legendre(1:nbends))
444 2119 : NULLIFY (inp_info%bend_legendre(i)%coeffs)
445 3058 : inp_info%bend_legendre(i)%order = 0
446 : END DO
447 : CALL read_bends_section(inp_info%bend_kind, inp_info%bend_a, inp_info%bend_b, inp_info%bend_c, &
448 : inp_info%bend_k, inp_info%bend_theta0, inp_info%bend_cb, &
449 : inp_info%bend_r012, inp_info%bend_r032, inp_info%bend_kbs12, &
450 : inp_info%bend_kbs32, inp_info%bend_kss, &
451 939 : inp_info%bend_legendre, tmp_section, ntot)
452 : END IF
453 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "BEND")
454 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nubs)
455 2639 : IF (explicit) THEN
456 939 : ntot = 0
457 939 : CALL reallocate(inp_info%ub_kind, 1, nubs)
458 939 : CALL reallocate(inp_info%ub_a, 1, nubs)
459 939 : CALL reallocate(inp_info%ub_b, 1, nubs)
460 939 : CALL reallocate(inp_info%ub_c, 1, nubs)
461 939 : CALL reallocate(inp_info%ub_k, 1, 3, 1, nubs)
462 939 : CALL reallocate(inp_info%ub_r0, 1, nubs)
463 : CALL read_ubs_section(inp_info%ub_kind, inp_info%ub_a, inp_info%ub_b, inp_info%ub_c, &
464 939 : inp_info%ub_k, inp_info%ub_r0, tmp_section, ntot)
465 : END IF
466 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "TORSION")
467 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=ntors)
468 2639 : IF (explicit) THEN
469 6 : ntot = 0
470 6 : CALL reallocate(inp_info%torsion_kind, 1, ntors)
471 6 : CALL reallocate(inp_info%torsion_a, 1, ntors)
472 6 : CALL reallocate(inp_info%torsion_b, 1, ntors)
473 6 : CALL reallocate(inp_info%torsion_c, 1, ntors)
474 6 : CALL reallocate(inp_info%torsion_d, 1, ntors)
475 6 : CALL reallocate(inp_info%torsion_k, 1, ntors)
476 6 : CALL reallocate(inp_info%torsion_m, 1, ntors)
477 6 : CALL reallocate(inp_info%torsion_phi0, 1, ntors)
478 : CALL read_torsions_section(inp_info%torsion_kind, inp_info%torsion_a, inp_info%torsion_b, &
479 : inp_info%torsion_c, inp_info%torsion_d, inp_info%torsion_k, inp_info%torsion_phi0, &
480 6 : inp_info%torsion_m, tmp_section, ntot)
481 : END IF
482 :
483 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "IMPROPER")
484 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nimpr)
485 2639 : IF (explicit) THEN
486 8 : ntot = 0
487 8 : CALL reallocate(inp_info%impr_kind, 1, nimpr)
488 8 : CALL reallocate(inp_info%impr_a, 1, nimpr)
489 8 : CALL reallocate(inp_info%impr_b, 1, nimpr)
490 8 : CALL reallocate(inp_info%impr_c, 1, nimpr)
491 8 : CALL reallocate(inp_info%impr_d, 1, nimpr)
492 8 : CALL reallocate(inp_info%impr_k, 1, nimpr)
493 8 : CALL reallocate(inp_info%impr_phi0, 1, nimpr)
494 : CALL read_improper_section(inp_info%impr_kind, inp_info%impr_a, inp_info%impr_b, &
495 : inp_info%impr_c, inp_info%impr_d, inp_info%impr_k, inp_info%impr_phi0, &
496 8 : tmp_section, ntot)
497 : END IF
498 :
499 2639 : tmp_section => section_vals_get_subs_vals(ff_section, "OPBEND")
500 2639 : CALL section_vals_get(tmp_section, explicit=explicit, n_repetition=nopbend)
501 2639 : IF (explicit) THEN
502 2 : ntot = 0
503 2 : CALL reallocate(inp_info%opbend_kind, 1, nopbend)
504 2 : CALL reallocate(inp_info%opbend_a, 1, nopbend)
505 2 : CALL reallocate(inp_info%opbend_b, 1, nopbend)
506 2 : CALL reallocate(inp_info%opbend_c, 1, nopbend)
507 2 : CALL reallocate(inp_info%opbend_d, 1, nopbend)
508 2 : CALL reallocate(inp_info%opbend_k, 1, nopbend)
509 2 : CALL reallocate(inp_info%opbend_phi0, 1, nopbend)
510 : CALL read_opbend_section(inp_info%opbend_kind, inp_info%opbend_a, inp_info%opbend_b, &
511 : inp_info%opbend_c, inp_info%opbend_d, inp_info%opbend_k, inp_info%opbend_phi0, &
512 2 : tmp_section, ntot)
513 : END IF
514 :
515 2639 : END SUBROUTINE read_force_field_section1
516 :
517 : ! **************************************************************************************************
518 : !> \brief Set up of the IPBV force fields
519 : !> \param at1 ...
520 : !> \param at2 ...
521 : !> \param ipbv ...
522 : !> \author teo
523 : ! **************************************************************************************************
524 48 : SUBROUTINE set_IPBV_ff(at1, at2, ipbv)
525 : CHARACTER(LEN=*), INTENT(IN) :: at1, at2
526 : TYPE(ipbv_pot_type), POINTER :: ipbv
527 :
528 48 : IF ((at1(1:1) == 'O') .AND. (at2(1:1) == 'O')) THEN
529 16 : ipbv%rcore = 0.9_dp ! a.u.
530 16 : ipbv%m = -1.2226442563398141E+11_dp ! Kelvin/a.u.
531 16 : ipbv%b = 1.1791292385486696E+11_dp ! Hartree
532 :
533 : ! Hartree*a.u.^2
534 16 : ipbv%a(2) = 4.786380682394_dp
535 16 : ipbv%a(3) = -1543.407053545_dp
536 16 : ipbv%a(4) = 88783.31188529_dp
537 16 : ipbv%a(5) = -2361200.155376_dp
538 16 : ipbv%a(6) = 35940504.84679_dp
539 16 : ipbv%a(7) = -339762743.6358_dp
540 16 : ipbv%a(8) = 2043874926.466_dp
541 16 : ipbv%a(9) = -7654856796.383_dp
542 16 : ipbv%a(10) = 16195251405.65_dp
543 16 : ipbv%a(11) = -13140392992.18_dp
544 16 : ipbv%a(12) = -9285572894.245_dp
545 16 : ipbv%a(13) = 8756947519.029_dp
546 16 : ipbv%a(14) = 15793297761.67_dp
547 16 : ipbv%a(15) = 12917180227.21_dp
548 32 : ELSEIF (((at1(1:1) == 'O') .AND. (at2(1:1) == 'H')) .OR. &
549 : ((at1(1:1) == 'H') .AND. (at2(1:1) == 'O'))) THEN
550 16 : ipbv%rcore = 2.95_dp ! a.u.
551 :
552 16 : ipbv%m = -0.004025691139759147_dp ! Hartree/a.u.
553 16 : ipbv%b = -2.193731138097428_dp ! Hartree
554 : ! Hartree*a.u.^2
555 16 : ipbv%a(2) = -195.7716013277_dp
556 16 : ipbv%a(3) = 15343.78613395_dp
557 16 : ipbv%a(4) = -530864.4586516_dp
558 16 : ipbv%a(5) = 10707934.39058_dp
559 16 : ipbv%a(6) = -140099704.7890_dp
560 16 : ipbv%a(7) = 1250943273.785_dp
561 16 : ipbv%a(8) = -7795458330.676_dp
562 16 : ipbv%a(9) = 33955897217.31_dp
563 16 : ipbv%a(10) = -101135640744.0_dp
564 16 : ipbv%a(11) = 193107995718.7_dp
565 16 : ipbv%a(12) = -193440560940.0_dp
566 16 : ipbv%a(13) = -4224406093.918E0_dp
567 16 : ipbv%a(14) = 217192386506.5E0_dp
568 16 : ipbv%a(15) = -157581228915.5_dp
569 16 : ELSEIF ((at1(1:1) == 'H') .AND. (at2(1:1) == 'H')) THEN
570 16 : ipbv%rcore = 3.165_dp ! a.u.
571 16 : ipbv%m = 0.002639704108787555_dp ! Hartree/a.u.
572 16 : ipbv%b = -0.2735482611857583_dp ! Hartree
573 : ! Hartree*a.u.^2
574 16 : ipbv%a(2) = -26.29456010782_dp
575 16 : ipbv%a(3) = 2373.352548248_dp
576 16 : ipbv%a(4) = -93880.43551360_dp
577 16 : ipbv%a(5) = 2154624.884809_dp
578 16 : ipbv%a(6) = -31965151.34955_dp
579 16 : ipbv%a(7) = 322781785.3278_dp
580 16 : ipbv%a(8) = -2271097368.668_dp
581 16 : ipbv%a(9) = 11169163192.90_dp
582 16 : ipbv%a(10) = -37684457778.47_dp
583 16 : ipbv%a(11) = 82562104256.03_dp
584 16 : ipbv%a(12) = -100510435213.4_dp
585 16 : ipbv%a(13) = 24570342714.65E0_dp
586 16 : ipbv%a(14) = 88766181532.94E0_dp
587 16 : ipbv%a(15) = -79705131323.98_dp
588 : ELSE
589 0 : CPABORT("IPBV only for WATER")
590 : END IF
591 48 : END SUBROUTINE set_IPBV_ff
592 :
593 : ! **************************************************************************************************
594 : !> \brief Set up of the BMHFT force fields
595 : !> \param at1 ...
596 : !> \param at2 ...
597 : !> \param ft ...
598 : !> \author teo
599 : ! **************************************************************************************************
600 12 : SUBROUTINE set_BMHFT_ff(at1, at2, ft)
601 : CHARACTER(LEN=*), INTENT(IN) :: at1, at2
602 : TYPE(ft_pot_type), POINTER :: ft
603 :
604 12 : ft%b = cp_unit_to_cp2k(3.1545_dp, "angstrom^-1")
605 12 : IF ((at1(1:2) == 'NA') .AND. (at2(1:2) == 'NA')) THEN
606 4 : ft%a = cp_unit_to_cp2k(424.097_dp, "eV")
607 4 : ft%c = cp_unit_to_cp2k(1.05_dp, "eV*angstrom^6")
608 4 : ft%d = cp_unit_to_cp2k(0.499_dp, "eV*angstrom^8")
609 8 : ELSEIF (((at1(1:2) == 'NA') .AND. (at2(1:2) == 'CL')) .OR. &
610 : ((at1(1:2) == 'CL') .AND. (at2(1:2) == 'NA'))) THEN
611 4 : ft%a = cp_unit_to_cp2k(1256.31_dp, "eV")
612 4 : ft%c = cp_unit_to_cp2k(7.00_dp, "eV*angstrom^6")
613 4 : ft%d = cp_unit_to_cp2k(8.676_dp, "eV*angstrom^8")
614 4 : ELSEIF ((at1(1:2) == 'CL') .AND. (at2(1:2) == 'CL')) THEN
615 4 : ft%a = cp_unit_to_cp2k(3488.998_dp, "eV")
616 4 : ft%c = cp_unit_to_cp2k(72.50_dp, "eV*angstrom^6")
617 4 : ft%d = cp_unit_to_cp2k(145.427_dp, "eV*angstrom^8")
618 : ELSE
619 0 : CPABORT("BMHFT only for NaCl")
620 : END IF
621 :
622 12 : END SUBROUTINE set_BMHFT_ff
623 :
624 : ! **************************************************************************************************
625 : !> \brief Set up of the BMHFTD force fields
626 : !> \author Mathieu Salanne 05.2010
627 : ! **************************************************************************************************
628 0 : SUBROUTINE set_BMHFTD_ff()
629 :
630 0 : CPABORT("No default parameters present for BMHFTD")
631 :
632 0 : END SUBROUTINE set_BMHFTD_ff
633 :
634 : ! **************************************************************************************************
635 : !> \brief Reads the EAM section
636 : !> \param nonbonded ...
637 : !> \param section ...
638 : !> \param start ...
639 : !> \param para_env ...
640 : !> \param mm_section ...
641 : !> \author teo
642 : ! **************************************************************************************************
643 12 : SUBROUTINE read_eam_section(nonbonded, section, start, para_env, mm_section)
644 : TYPE(pair_potential_p_type), POINTER :: nonbonded
645 : TYPE(section_vals_type), POINTER :: section
646 : INTEGER, INTENT(IN) :: start
647 : TYPE(mp_para_env_type), POINTER :: para_env
648 : TYPE(section_vals_type), POINTER :: mm_section
649 :
650 : CHARACTER(LEN=default_string_length), &
651 12 : DIMENSION(:), POINTER :: atm_names
652 : INTEGER :: isec, n_items
653 :
654 12 : CALL section_vals_get(section, n_repetition=n_items)
655 32 : DO isec = 1, n_items
656 20 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
657 :
658 40 : nonbonded%pot(start + isec)%pot%type = ea_type
659 20 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
660 20 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
661 20 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
662 20 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
663 : CALL section_vals_val_get(section, "PARM_FILE_NAME", i_rep_section=isec, &
664 20 : c_val=nonbonded%pot(start + isec)%pot%set(1)%eam%eam_file_name)
665 20 : CALL read_eam_data(nonbonded%pot(start + isec)%pot%set(1)%eam, para_env, mm_section)
666 32 : nonbonded%pot(start + isec)%pot%rcutsq = nonbonded%pot(start + isec)%pot%set(1)%eam%acutal**2
667 : END DO
668 12 : END SUBROUTINE read_eam_section
669 :
670 : ! **************************************************************************************************
671 : !> \brief Reads the DEEPMD section
672 : !> \param nonbonded ...
673 : !> \param section ...
674 : !> \param start ...
675 : !> \author teo
676 : ! **************************************************************************************************
677 2 : SUBROUTINE read_deepmd_section(nonbonded, section, start)
678 : TYPE(pair_potential_p_type), POINTER :: nonbonded
679 : TYPE(section_vals_type), POINTER :: section
680 : INTEGER, INTENT(IN) :: start
681 :
682 : CHARACTER(LEN=default_string_length) :: deepmd_file_name
683 : CHARACTER(LEN=default_string_length), &
684 2 : DIMENSION(:), POINTER :: atm_names
685 : INTEGER :: isec, jsec, n_items
686 2 : INTEGER, DIMENSION(:), POINTER :: atm_deepmd_types
687 :
688 : n_items = 1
689 2 : isec = 1
690 2 : n_items = isec*n_items
691 2 : CALL section_vals_val_get(section, "ATOMS", c_vals=atm_names)
692 2 : CALL section_vals_val_get(section, "ATOMS_DEEPMD_TYPE", i_vals=atm_deepmd_types)
693 2 : CALL section_vals_val_get(section, "POT_FILE_NAME", c_val=deepmd_file_name)
694 :
695 4 : DO isec = 1, SIZE(atm_names)
696 6 : DO jsec = isec, SIZE(atm_names)
697 4 : nonbonded%pot(start + n_items)%pot%type = deepmd_type
698 2 : nonbonded%pot(start + n_items)%pot%at1 = atm_names(isec)
699 2 : nonbonded%pot(start + n_items)%pot%at2 = atm_names(jsec)
700 2 : CALL uppercase(nonbonded%pot(start + n_items)%pot%at1)
701 2 : CALL uppercase(nonbonded%pot(start + n_items)%pot%at2)
702 :
703 2 : nonbonded%pot(start + n_items)%pot%set(1)%deepmd%deepmd_file_name = discover_file(deepmd_file_name)
704 2 : nonbonded%pot(start + n_items)%pot%set(1)%deepmd%atom_deepmd_type = atm_deepmd_types(isec)
705 2 : nonbonded%pot(start + n_items)%pot%rcutsq = 0.0_dp
706 4 : n_items = n_items + 1
707 : END DO
708 : END DO
709 2 : END SUBROUTINE read_deepmd_section
710 :
711 : ! **************************************************************************************************
712 : !> \brief Reads the QUIP section
713 : !> \param nonbonded ...
714 : !> \param section ...
715 : !> \param start ...
716 : !> \author teo
717 : ! **************************************************************************************************
718 2 : SUBROUTINE read_quip_section(nonbonded, section, start)
719 : TYPE(pair_potential_p_type), POINTER :: nonbonded
720 : TYPE(section_vals_type), POINTER :: section
721 : INTEGER, INTENT(IN) :: start
722 :
723 : CHARACTER(LEN=default_string_length), &
724 2 : DIMENSION(:), POINTER :: args_str, atm_names
725 : INTEGER :: is, isec, n_calc_args, n_items
726 :
727 2 : CALL section_vals_get(section, n_repetition=n_items)
728 4 : DO isec = 1, n_items
729 2 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
730 :
731 4 : nonbonded%pot(start + isec)%pot%type = quip_type
732 2 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
733 2 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
734 2 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
735 2 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
736 : CALL section_vals_val_get(section, "PARM_FILE_NAME", i_rep_section=isec, &
737 2 : c_val=nonbonded%pot(start + isec)%pot%set(1)%quip%quip_file_name)
738 : CALL section_vals_val_get(section, "INIT_ARGS", i_rep_section=isec, &
739 2 : c_vals=args_str)
740 2 : nonbonded%pot(start + isec)%pot%set(1)%quip%init_args = ""
741 6 : DO is = 1, SIZE(args_str)
742 : nonbonded%pot(start + isec)%pot%set(1)%quip%init_args = &
743 : TRIM(nonbonded%pot(start + isec)%pot%set(1)%quip%init_args)// &
744 6 : " "//TRIM(args_str(is))
745 : END DO ! is
746 : CALL section_vals_val_get(section, "CALC_ARGS", i_rep_section=isec, &
747 2 : n_rep_val=n_calc_args)
748 2 : IF (n_calc_args > 0) THEN
749 : CALL section_vals_val_get(section, "CALC_ARGS", i_rep_section=isec, &
750 0 : c_vals=args_str)
751 0 : DO is = 1, SIZE(args_str)
752 : nonbonded%pot(start + isec)%pot%set(1)%quip%calc_args = &
753 : TRIM(nonbonded%pot(start + isec)%pot%set(1)%quip%calc_args)// &
754 0 : " "//TRIM(args_str(is))
755 : END DO ! is
756 : END IF
757 6 : nonbonded%pot(start + isec)%pot%rcutsq = 0.0_dp
758 : END DO
759 2 : END SUBROUTINE read_quip_section
760 :
761 : ! **************************************************************************************************
762 : !> \brief Reads the NEQUIP section
763 : !> \param nonbonded ...
764 : !> \param section ...
765 : !> \param start ...
766 : !> \author Gabriele Tocci
767 : ! **************************************************************************************************
768 4 : SUBROUTINE read_nequip_section(nonbonded, section, start)
769 : TYPE(pair_potential_p_type), POINTER :: nonbonded
770 : TYPE(section_vals_type), POINTER :: section
771 : INTEGER, INTENT(IN) :: start
772 :
773 : CHARACTER(LEN=default_string_length) :: nequip_file_name, unit_cell, &
774 : unit_coords, unit_energy, unit_forces
775 : CHARACTER(LEN=default_string_length), &
776 4 : DIMENSION(:), POINTER :: atm_names
777 : INTEGER :: isec, jsec, n_items
778 :
779 : n_items = 1
780 4 : isec = 1
781 4 : n_items = isec*n_items
782 4 : CALL section_vals_val_get(section, "ATOMS", c_vals=atm_names)
783 4 : CALL section_vals_val_get(section, "PARM_FILE_NAME", c_val=nequip_file_name)
784 4 : CALL section_vals_val_get(section, "UNIT_COORDS", c_val=unit_coords)
785 4 : CALL section_vals_val_get(section, "UNIT_ENERGY", c_val=unit_energy)
786 4 : CALL section_vals_val_get(section, "UNIT_FORCES", c_val=unit_forces)
787 4 : CALL section_vals_val_get(section, "UNIT_CELL", c_val=unit_cell)
788 :
789 12 : DO isec = 1, SIZE(atm_names)
790 24 : DO jsec = isec, SIZE(atm_names)
791 24 : nonbonded%pot(start + n_items)%pot%type = nequip_type
792 12 : nonbonded%pot(start + n_items)%pot%at1 = atm_names(isec)
793 12 : nonbonded%pot(start + n_items)%pot%at2 = atm_names(jsec)
794 12 : CALL uppercase(nonbonded%pot(start + n_items)%pot%at1)
795 12 : CALL uppercase(nonbonded%pot(start + n_items)%pot%at2)
796 12 : nonbonded%pot(start + n_items)%pot%set(1)%nequip%nequip_file_name = discover_file(nequip_file_name)
797 12 : nonbonded%pot(start + n_items)%pot%set(1)%nequip%unit_coords = unit_coords
798 12 : nonbonded%pot(start + n_items)%pot%set(1)%nequip%unit_forces = unit_forces
799 12 : nonbonded%pot(start + n_items)%pot%set(1)%nequip%unit_energy = unit_energy
800 12 : nonbonded%pot(start + n_items)%pot%set(1)%nequip%unit_cell = unit_cell
801 12 : CALL read_nequip_data(nonbonded%pot(start + n_items)%pot%set(1)%nequip)
802 12 : nonbonded%pot(start + n_items)%pot%rcutsq = nonbonded%pot(start + n_items)%pot%set(1)%nequip%rcutsq
803 20 : n_items = n_items + 1
804 : END DO
805 : END DO
806 :
807 4 : END SUBROUTINE read_nequip_section
808 :
809 : ! **************************************************************************************************
810 : !> \brief Reads the ALLEGRO section
811 : !> \param nonbonded ...
812 : !> \param section ...
813 : !> \param start ...
814 : !> \author Gabriele Tocci
815 : ! **************************************************************************************************
816 4 : SUBROUTINE read_allegro_section(nonbonded, section, start)
817 : TYPE(pair_potential_p_type), POINTER :: nonbonded
818 : TYPE(section_vals_type), POINTER :: section
819 : INTEGER, INTENT(IN) :: start
820 :
821 : CHARACTER(LEN=default_string_length) :: allegro_file_name, unit_cell, &
822 : unit_coords, unit_energy, unit_forces
823 : CHARACTER(LEN=default_string_length), &
824 4 : DIMENSION(:), POINTER :: atm_names
825 : INTEGER :: isec, jsec, n_items
826 :
827 : n_items = 1
828 4 : isec = 1
829 4 : n_items = isec*n_items
830 4 : CALL section_vals_val_get(section, "ATOMS", c_vals=atm_names)
831 4 : CALL section_vals_val_get(section, "PARM_FILE_NAME", c_val=allegro_file_name)
832 4 : CALL section_vals_val_get(section, "UNIT_COORDS", c_val=unit_coords)
833 4 : CALL section_vals_val_get(section, "UNIT_ENERGY", c_val=unit_energy)
834 4 : CALL section_vals_val_get(section, "UNIT_FORCES", c_val=unit_forces)
835 4 : CALL section_vals_val_get(section, "UNIT_CELL", c_val=unit_cell)
836 :
837 12 : DO isec = 1, SIZE(atm_names)
838 26 : DO jsec = isec, SIZE(atm_names)
839 28 : nonbonded%pot(start + n_items)%pot%type = allegro_type
840 14 : nonbonded%pot(start + n_items)%pot%at1 = atm_names(isec)
841 14 : nonbonded%pot(start + n_items)%pot%at2 = atm_names(jsec)
842 14 : CALL uppercase(nonbonded%pot(start + n_items)%pot%at1)
843 14 : CALL uppercase(nonbonded%pot(start + n_items)%pot%at2)
844 14 : nonbonded%pot(start + n_items)%pot%set(1)%allegro%allegro_file_name = discover_file(allegro_file_name)
845 14 : nonbonded%pot(start + n_items)%pot%set(1)%allegro%unit_coords = unit_coords
846 14 : nonbonded%pot(start + n_items)%pot%set(1)%allegro%unit_forces = unit_forces
847 14 : nonbonded%pot(start + n_items)%pot%set(1)%allegro%unit_energy = unit_energy
848 14 : nonbonded%pot(start + n_items)%pot%set(1)%allegro%unit_cell = unit_cell
849 14 : CALL read_allegro_data(nonbonded%pot(start + n_items)%pot%set(1)%allegro)
850 14 : nonbonded%pot(start + n_items)%pot%rcutsq = nonbonded%pot(start + n_items)%pot%set(1)%allegro%rcutsq
851 22 : n_items = n_items + 1
852 : END DO
853 : END DO
854 :
855 4 : END SUBROUTINE read_allegro_section
856 :
857 : ! **************************************************************************************************
858 : !> \brief Reads the LJ section
859 : !> \param nonbonded ...
860 : !> \param section ...
861 : !> \param start ...
862 : !> \author teo
863 : ! **************************************************************************************************
864 1006 : SUBROUTINE read_lj_section(nonbonded, section, start)
865 : TYPE(pair_potential_p_type), POINTER :: nonbonded
866 : TYPE(section_vals_type), POINTER :: section
867 : INTEGER, INTENT(IN) :: start
868 :
869 : CHARACTER(LEN=default_string_length), &
870 1006 : DIMENSION(:), POINTER :: atm_names
871 : INTEGER :: isec, n_items, n_rep
872 : REAL(KIND=dp) :: epsilon, rcut, sigma
873 :
874 1006 : CALL section_vals_get(section, n_repetition=n_items)
875 4794 : DO isec = 1, n_items
876 3788 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
877 3788 : CALL section_vals_val_get(section, "EPSILON", i_rep_section=isec, r_val=epsilon)
878 3788 : CALL section_vals_val_get(section, "SIGMA", i_rep_section=isec, r_val=sigma)
879 3788 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
880 :
881 7576 : nonbonded%pot(start + isec)%pot%type = lj_charmm_type
882 3788 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
883 3788 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
884 3788 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
885 3788 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
886 3788 : nonbonded%pot(start + isec)%pot%set(1)%lj%epsilon = epsilon
887 3788 : nonbonded%pot(start + isec)%pot%set(1)%lj%sigma6 = sigma**6
888 3788 : nonbonded%pot(start + isec)%pot%set(1)%lj%sigma12 = sigma**12
889 3788 : nonbonded%pot(start + isec)%pot%rcutsq = rcut*rcut
890 : !
891 3788 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
892 3788 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
893 2 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
894 3788 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
895 3788 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
896 12372 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
897 : END DO
898 1006 : END SUBROUTINE read_lj_section
899 :
900 : ! **************************************************************************************************
901 : !> \brief Reads the WILLIAMS section
902 : !> \param nonbonded ...
903 : !> \param section ...
904 : !> \param start ...
905 : !> \author teo
906 : ! **************************************************************************************************
907 375 : SUBROUTINE read_wl_section(nonbonded, section, start)
908 : TYPE(pair_potential_p_type), POINTER :: nonbonded
909 : TYPE(section_vals_type), POINTER :: section
910 : INTEGER, INTENT(IN) :: start
911 :
912 : CHARACTER(LEN=default_string_length), &
913 375 : DIMENSION(:), POINTER :: atm_names
914 : INTEGER :: isec, n_items, n_rep
915 : REAL(KIND=dp) :: a, b, c, rcut
916 :
917 375 : CALL section_vals_get(section, n_repetition=n_items)
918 1386 : DO isec = 1, n_items
919 1011 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
920 1011 : CALL section_vals_val_get(section, "A", i_rep_section=isec, r_val=a)
921 1011 : CALL section_vals_val_get(section, "B", i_rep_section=isec, r_val=b)
922 1011 : CALL section_vals_val_get(section, "C", i_rep_section=isec, r_val=c)
923 1011 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
924 :
925 2022 : nonbonded%pot(start + isec)%pot%type = wl_type
926 1011 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
927 1011 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
928 1011 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
929 1011 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
930 1011 : nonbonded%pot(start + isec)%pot%set(1)%willis%a = a
931 1011 : nonbonded%pot(start + isec)%pot%set(1)%willis%b = b
932 1011 : nonbonded%pot(start + isec)%pot%set(1)%willis%c = c
933 1011 : nonbonded%pot(start + isec)%pot%rcutsq = rcut*rcut
934 : !
935 1011 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
936 1011 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
937 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
938 1011 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
939 1011 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
940 3408 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
941 : END DO
942 375 : END SUBROUTINE read_wl_section
943 :
944 : ! **************************************************************************************************
945 : !> \brief Reads the GOODWIN section
946 : !> \param nonbonded ...
947 : !> \param section ...
948 : !> \param start ...
949 : !> \author teo
950 : ! **************************************************************************************************
951 0 : SUBROUTINE read_gd_section(nonbonded, section, start)
952 : TYPE(pair_potential_p_type), POINTER :: nonbonded
953 : TYPE(section_vals_type), POINTER :: section
954 : INTEGER, INTENT(IN) :: start
955 :
956 : CHARACTER(LEN=default_string_length), &
957 0 : DIMENSION(:), POINTER :: atm_names
958 : INTEGER :: isec, m, mc, n_items, n_rep
959 : REAL(KIND=dp) :: d, dc, rcut, vr0
960 :
961 0 : CALL section_vals_get(section, n_repetition=n_items)
962 0 : DO isec = 1, n_items
963 0 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
964 0 : CALL section_vals_val_get(section, "VR0", i_rep_section=isec, r_val=vr0)
965 0 : CALL section_vals_val_get(section, "D", i_rep_section=isec, r_val=d)
966 0 : CALL section_vals_val_get(section, "DC", i_rep_section=isec, r_val=dc)
967 0 : CALL section_vals_val_get(section, "M", i_rep_section=isec, i_val=m)
968 0 : CALL section_vals_val_get(section, "MC", i_rep_section=isec, i_val=mc)
969 0 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
970 :
971 0 : nonbonded%pot(start + isec)%pot%type = gw_type
972 0 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
973 0 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
974 0 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
975 0 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
976 0 : nonbonded%pot(start + isec)%pot%set(1)%goodwin%vr0 = vr0
977 0 : nonbonded%pot(start + isec)%pot%set(1)%goodwin%d = d
978 0 : nonbonded%pot(start + isec)%pot%set(1)%goodwin%dc = dc
979 0 : nonbonded%pot(start + isec)%pot%set(1)%goodwin%m = m
980 0 : nonbonded%pot(start + isec)%pot%set(1)%goodwin%mc = mc
981 0 : nonbonded%pot(start + isec)%pot%rcutsq = rcut*rcut
982 : !
983 0 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
984 0 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
985 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
986 0 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
987 0 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
988 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
989 : END DO
990 0 : END SUBROUTINE read_gd_section
991 :
992 : ! **************************************************************************************************
993 : !> \brief Reads the IPBV section
994 : !> \param nonbonded ...
995 : !> \param section ...
996 : !> \param start ...
997 : !> \author teo
998 : ! **************************************************************************************************
999 16 : SUBROUTINE read_ipbv_section(nonbonded, section, start)
1000 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1001 : TYPE(section_vals_type), POINTER :: section
1002 : INTEGER, INTENT(IN) :: start
1003 :
1004 : CHARACTER(LEN=default_string_length), &
1005 16 : DIMENSION(:), POINTER :: atm_names
1006 : INTEGER :: isec, n_items, n_rep
1007 : REAL(KIND=dp) :: rcut
1008 :
1009 16 : CALL section_vals_get(section, n_repetition=n_items)
1010 64 : DO isec = 1, n_items
1011 48 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1012 96 : nonbonded%pot(start + isec)%pot%type = ip_type
1013 48 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1014 48 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1015 48 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1016 48 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1017 : CALL set_IPBV_ff(nonbonded%pot(start + isec)%pot%at1, nonbonded%pot(start + isec)%pot%at2, &
1018 48 : nonbonded%pot(start + isec)%pot%set(1)%ipbv)
1019 48 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
1020 48 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1021 : !
1022 48 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
1023 48 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
1024 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
1025 48 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
1026 48 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
1027 112 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
1028 : END DO
1029 16 : END SUBROUTINE read_ipbv_section
1030 :
1031 : ! **************************************************************************************************
1032 : !> \brief Reads the BMHFT section
1033 : !> \param nonbonded ...
1034 : !> \param section ...
1035 : !> \param start ...
1036 : !> \author teo
1037 : ! **************************************************************************************************
1038 4 : SUBROUTINE read_bmhft_section(nonbonded, section, start)
1039 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1040 : TYPE(section_vals_type), POINTER :: section
1041 : INTEGER, INTENT(IN) :: start
1042 :
1043 : CHARACTER(LEN=default_string_length), DIMENSION(2) :: map_atoms
1044 : CHARACTER(LEN=default_string_length), &
1045 4 : DIMENSION(:), POINTER :: atm_names
1046 : INTEGER :: i, isec, n_items, n_rep
1047 : REAL(KIND=dp) :: rcut
1048 :
1049 4 : CALL section_vals_get(section, n_repetition=n_items)
1050 16 : DO isec = 1, n_items
1051 12 : CALL cite_reference(Tosi1964a)
1052 12 : CALL cite_reference(Tosi1964b)
1053 12 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1054 24 : nonbonded%pot(start + isec)%pot%type = ft_type
1055 12 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1056 12 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1057 12 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1058 12 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1059 :
1060 12 : CALL section_vals_val_get(section, "A", i_rep_section=isec, n_rep_val=i)
1061 12 : IF (i == 1) THEN
1062 : CALL section_vals_val_get(section, "A", i_rep_section=isec, &
1063 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ft%a)
1064 : CALL section_vals_val_get(section, "B", i_rep_section=isec, &
1065 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ft%b)
1066 : CALL section_vals_val_get(section, "C", i_rep_section=isec, &
1067 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ft%c)
1068 : CALL section_vals_val_get(section, "D", i_rep_section=isec, &
1069 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ft%d)
1070 : ELSE
1071 12 : CALL section_vals_val_get(section, "MAP_ATOMS", i_rep_section=isec, c_vals=atm_names)
1072 36 : map_atoms = atm_names
1073 12 : CALL uppercase(map_atoms(1))
1074 12 : CALL uppercase(map_atoms(2))
1075 12 : CALL set_BMHFT_ff(map_atoms(1), map_atoms(2), nonbonded%pot(start + isec)%pot%set(1)%ft)
1076 : END IF
1077 12 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
1078 12 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1079 : !
1080 12 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
1081 12 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
1082 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
1083 12 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
1084 12 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
1085 40 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
1086 : END DO
1087 4 : END SUBROUTINE read_bmhft_section
1088 :
1089 : ! **************************************************************************************************
1090 : !> \brief Reads the BMHFTD section
1091 : !> \param nonbonded ...
1092 : !> \param section ...
1093 : !> \param start ...
1094 : !> \author Mathieu Salanne 05.2010
1095 : ! **************************************************************************************************
1096 18 : SUBROUTINE read_bmhftd_section(nonbonded, section, start)
1097 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1098 : TYPE(section_vals_type), POINTER :: section
1099 : INTEGER, INTENT(IN) :: start
1100 :
1101 : CHARACTER(LEN=default_string_length), DIMENSION(2) :: map_atoms
1102 : CHARACTER(LEN=default_string_length), &
1103 18 : DIMENSION(:), POINTER :: atm_names
1104 : INTEGER :: i, isec, n_items, n_rep
1105 : REAL(KIND=dp) :: rcut
1106 18 : REAL(KIND=dp), DIMENSION(:), POINTER :: bd_vals
1107 :
1108 18 : NULLIFY (bd_vals)
1109 :
1110 18 : CALL section_vals_get(section, n_repetition=n_items)
1111 :
1112 84 : DO isec = 1, n_items
1113 66 : CALL cite_reference(Tosi1964a)
1114 66 : CALL cite_reference(Tosi1964b)
1115 66 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1116 132 : nonbonded%pot(start + isec)%pot%type = ftd_type
1117 66 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1118 66 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1119 66 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1120 66 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1121 :
1122 66 : CALL section_vals_val_get(section, "A", i_rep_section=isec, n_rep_val=i)
1123 66 : IF (i == 1) THEN
1124 : CALL section_vals_val_get(section, "A", i_rep_section=isec, &
1125 66 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ftd%a)
1126 : CALL section_vals_val_get(section, "B", i_rep_section=isec, &
1127 66 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ftd%b)
1128 : CALL section_vals_val_get(section, "C", i_rep_section=isec, &
1129 66 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ftd%c)
1130 : CALL section_vals_val_get(section, "D", i_rep_section=isec, &
1131 66 : r_val=nonbonded%pot(start + isec)%pot%set(1)%ftd%d)
1132 66 : CALL section_vals_val_get(section, "BD", i_rep_section=isec, r_vals=bd_vals)
1133 66 : IF (ASSOCIATED(bd_vals)) THEN
1134 66 : SELECT CASE (SIZE(bd_vals))
1135 : CASE (0)
1136 0 : CPABORT("No values specified for parameter BD in section &BMHFTD")
1137 : CASE (1)
1138 186 : nonbonded%pot(start + isec)%pot%set(1)%ftd%bd(1:2) = bd_vals(1)
1139 : CASE (2)
1140 24 : nonbonded%pot(start + isec)%pot%set(1)%ftd%bd(1:2) = bd_vals(1:2)
1141 : CASE DEFAULT
1142 66 : CPABORT("Too many values specified for parameter BD in section &BMHFTD")
1143 : END SELECT
1144 : ELSE
1145 0 : CPABORT("Parameter BD in section &BMHFTD was not specified")
1146 : END IF
1147 : ELSE
1148 0 : CALL section_vals_val_get(section, "MAP_ATOMS", i_rep_section=isec, c_vals=atm_names)
1149 0 : map_atoms = atm_names
1150 0 : CALL uppercase(map_atoms(1))
1151 0 : CALL uppercase(map_atoms(2))
1152 0 : CALL set_BMHFTD_ff()
1153 : END IF
1154 66 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
1155 66 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1156 : !
1157 66 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
1158 66 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
1159 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
1160 66 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
1161 66 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
1162 216 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
1163 : END DO
1164 18 : END SUBROUTINE read_bmhftd_section
1165 :
1166 : ! **************************************************************************************************
1167 : !> \brief Reads the Buckingham 4 Ranges potential section
1168 : !> \param nonbonded ...
1169 : !> \param section ...
1170 : !> \param start ...
1171 : !> \par History
1172 : !> MK (11.11.2010): Automatic fit of the (default) polynomial coefficients
1173 : !> \author MI,MK
1174 : ! **************************************************************************************************
1175 262 : SUBROUTINE read_b4_section(nonbonded, section, start)
1176 :
1177 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1178 : TYPE(section_vals_type), POINTER :: section
1179 : INTEGER, INTENT(IN) :: start
1180 :
1181 : CHARACTER(LEN=default_string_length), &
1182 262 : DIMENSION(:), POINTER :: atm_names
1183 : INTEGER :: i, ir, isec, n_items, n_rep, np1, np2
1184 : LOGICAL :: explicit_poly1, explicit_poly2
1185 : REAL(KIND=dp) :: a, b, c, eval_error, r1, r2, r3, rcut
1186 : REAL(KIND=dp), DIMENSION(10) :: v, x
1187 : REAL(KIND=dp), DIMENSION(10, 10) :: p, p_inv
1188 262 : REAL(KIND=dp), DIMENSION(:), POINTER :: coeff1, coeff2, list
1189 :
1190 262 : NULLIFY (coeff1)
1191 262 : NULLIFY (coeff2)
1192 :
1193 262 : CALL section_vals_get(section, n_repetition=n_items)
1194 :
1195 524 : DO isec = 1, n_items
1196 262 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1197 262 : CALL section_vals_val_get(section, "A", i_rep_section=isec, r_val=a)
1198 262 : CALL section_vals_val_get(section, "B", i_rep_section=isec, r_val=b)
1199 262 : CALL section_vals_val_get(section, "C", i_rep_section=isec, r_val=c)
1200 262 : CALL section_vals_val_get(section, "R1", i_rep_section=isec, r_val=r1)
1201 262 : CALL section_vals_val_get(section, "R2", i_rep_section=isec, r_val=r2)
1202 262 : CALL section_vals_val_get(section, "R3", i_rep_section=isec, r_val=r3)
1203 262 : CALL section_vals_val_get(section, "POLY1", explicit=explicit_poly1, n_rep_val=n_rep)
1204 : ! Check if polynomial coefficients were specified for range 2 and 3 explicitly
1205 262 : IF (explicit_poly1) THEN
1206 84 : np1 = 0
1207 168 : DO ir = 1, n_rep
1208 84 : NULLIFY (list)
1209 84 : CALL section_vals_val_get(section, "POLY1", i_rep_val=ir, r_vals=list)
1210 168 : IF (ASSOCIATED(list)) THEN
1211 84 : CALL reallocate(coeff1, 0, np1 + SIZE(list) - 1)
1212 588 : DO i = 1, SIZE(list)
1213 588 : coeff1(i + np1 - 1) = list(i)
1214 : END DO
1215 84 : np1 = np1 + SIZE(list)
1216 : END IF
1217 : END DO
1218 : END IF
1219 262 : CALL section_vals_val_get(section, "POLY2", explicit=explicit_poly2, n_rep_val=n_rep)
1220 262 : IF (explicit_poly2) THEN
1221 84 : np2 = 0
1222 168 : DO ir = 1, n_rep
1223 84 : NULLIFY (list)
1224 84 : CALL section_vals_val_get(section, "POLY2", i_rep_val=ir, r_vals=list)
1225 168 : IF (ASSOCIATED(list)) THEN
1226 84 : CALL reallocate(coeff2, 0, np2 + SIZE(list) - 1)
1227 420 : DO i = 1, SIZE(list)
1228 420 : coeff2(i + np2 - 1) = list(i)
1229 : END DO
1230 84 : np2 = np2 + SIZE(list)
1231 : END IF
1232 : END DO
1233 : END IF
1234 : ! Default is a 5th/3rd-order polynomial fit
1235 262 : IF ((.NOT. explicit_poly1) .OR. (.NOT. explicit_poly2)) THEN
1236 : ! Build matrix p and vector v to calculate the polynomial coefficients
1237 : ! in the vector x from p*x = v
1238 178 : p(:, :) = 0.0_dp
1239 : ! Row 1: Match the 5th-order polynomial and the potential at r1
1240 178 : p(1, 1) = 1.0_dp
1241 1068 : DO i = 2, 6
1242 1068 : p(1, i) = p(1, i - 1)*r1
1243 : END DO
1244 : ! Row 2: Match the first derivatives of the 5th-order polynomial and the potential at r1
1245 1068 : DO i = 2, 6
1246 1068 : p(2, i) = REAL(i - 1, KIND=dp)*p(1, i - 1)
1247 : END DO
1248 : ! Row 3: Match the second derivatives of the 5th-order polynomial and the potential at r1
1249 890 : DO i = 3, 6
1250 890 : p(3, i) = REAL(i - 1, KIND=dp)*p(2, i - 1)
1251 : END DO
1252 : ! Row 4: Match the 5th-order and the 3rd-order polynomials at r2
1253 178 : p(4, 1) = 1.0_dp
1254 1068 : DO i = 2, 6
1255 1068 : p(4, i) = p(4, i - 1)*r2
1256 : END DO
1257 178 : p(4, 7) = -1.0_dp
1258 712 : DO i = 8, 10
1259 712 : p(4, i) = p(4, i - 1)*r2
1260 : END DO
1261 : ! Row 5: Match the first derivatives of the 5th-order and the 3rd-order polynomials at r2
1262 1068 : DO i = 2, 6
1263 1068 : p(5, i) = REAL(i - 1, KIND=dp)*p(4, i - 1)
1264 : END DO
1265 712 : DO i = 8, 10
1266 712 : p(5, i) = REAL(i - 7, KIND=dp)*p(4, i - 1)
1267 : END DO
1268 : ! Row 6: Match the second derivatives of the 5th-order and the 3rd-order polynomials at r2
1269 890 : DO i = 3, 6
1270 890 : p(6, i) = REAL(i - 1, KIND=dp)*p(5, i - 1)
1271 : END DO
1272 534 : DO i = 9, 10
1273 534 : p(6, i) = REAL(i - 7, KIND=dp)*p(5, i - 1)
1274 : END DO
1275 : ! Row 7: Minimum at r2, ie. the first derivative of the 3rd-order polynomial has to be zero at r2
1276 712 : DO i = 8, 10
1277 712 : p(7, i) = -p(5, i)
1278 : END DO
1279 : ! Row 8: Match the 3rd-order polynomial and the potential at r3
1280 178 : p(8, 7) = 1.0_dp
1281 712 : DO i = 8, 10
1282 712 : p(8, i) = p(8, i - 1)*r3
1283 : END DO
1284 : ! Row 9: Match the first derivatives of the 3rd-order polynomial and the potential at r3
1285 712 : DO i = 8, 10
1286 712 : p(9, i) = REAL(i - 7, KIND=dp)*p(8, i - 1)
1287 : END DO
1288 : ! Row 10: Match the second derivatives of the 3rd-order polynomial and the potential at r3
1289 534 : DO i = 9, 10
1290 534 : p(10, i) = REAL(i - 7, KIND=dp)*p(9, i - 1)
1291 : END DO
1292 : ! Build the vector v
1293 178 : v(1) = a*EXP(-b*r1)
1294 178 : v(2) = -b*v(1)
1295 178 : v(3) = -b*v(2)
1296 890 : v(4:7) = 0.0_dp
1297 178 : v(8) = -c/p(8, 10)**2 ! = -c/r3**6
1298 178 : v(9) = -6.0_dp*v(8)/r3
1299 178 : v(10) = -7.0_dp*v(9)/r3
1300 : ! Calculate p_inv the inverse of the matrix p
1301 178 : p_inv(:, :) = 0.0_dp
1302 178 : CALL invert_matrix(p, p_inv, eval_error)
1303 178 : IF (eval_error >= 1.0E-8_dp) &
1304 : CALL cp_warn(__LOCATION__, &
1305 : "The polynomial fit for the BUCK4RANGES potential is only accurate to "// &
1306 0 : TRIM(cp_to_string(eval_error)))
1307 : ! Get the 6 coefficients of the 5th-order polynomial -> x(1:6)
1308 : ! and the 4 coefficients of the 3rd-order polynomial -> x(7:10)
1309 19758 : x(:) = MATMUL(p_inv(:, :), v(:))
1310 : ELSE
1311 84 : x(:) = 0.0_dp
1312 : END IF
1313 :
1314 262 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
1315 :
1316 524 : nonbonded%pot(start + isec)%pot%type = b4_type
1317 262 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1318 262 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1319 262 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1320 262 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1321 262 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%a = a
1322 262 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%b = b
1323 262 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%c = c
1324 262 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%r1 = r1
1325 262 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%r2 = r2
1326 262 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%r3 = r3
1327 262 : IF ((.NOT. explicit_poly1) .OR. (.NOT. explicit_poly2)) THEN
1328 178 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%npoly1 = 5
1329 1246 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%poly1(0:5) = x(1:6)
1330 178 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%npoly2 = 3
1331 890 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%poly2(0:3) = x(7:10)
1332 : ELSE
1333 84 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%npoly1 = np1 - 1
1334 84 : CPASSERT(np1 - 1 <= 10)
1335 1092 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%poly1(0:np1 - 1) = coeff1(0:np1 - 1)
1336 84 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%npoly2 = np2 - 1
1337 84 : CPASSERT(np2 - 1 <= 10)
1338 756 : nonbonded%pot(start + isec)%pot%set(1)%buck4r%poly2(0:np2 - 1) = coeff2(0:np2 - 1)
1339 : END IF
1340 262 : nonbonded%pot(start + isec)%pot%rcutsq = rcut*rcut
1341 :
1342 262 : IF (ASSOCIATED(coeff1)) THEN
1343 84 : DEALLOCATE (coeff1)
1344 : END IF
1345 262 : IF (ASSOCIATED(coeff2)) THEN
1346 84 : DEALLOCATE (coeff2)
1347 : END IF
1348 262 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
1349 262 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
1350 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
1351 262 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
1352 262 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
1353 1572 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
1354 : END DO
1355 :
1356 262 : END SUBROUTINE read_b4_section
1357 :
1358 : ! **************************************************************************************************
1359 : !> \brief Reads the GENPOT - generic potential section
1360 : !> \param nonbonded ...
1361 : !> \param section ...
1362 : !> \param start ...
1363 : !> \author Teodoro Laino - 10.2006
1364 : ! **************************************************************************************************
1365 576 : SUBROUTINE read_gp_section(nonbonded, section, start)
1366 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1367 : TYPE(section_vals_type), POINTER :: section
1368 : INTEGER, INTENT(IN) :: start
1369 :
1370 : CHARACTER(LEN=default_string_length), &
1371 576 : DIMENSION(:), POINTER :: atm_names
1372 : INTEGER :: isec, n_items, n_rep
1373 : REAL(KIND=dp) :: rcut
1374 :
1375 576 : CALL section_vals_get(section, n_repetition=n_items)
1376 3794 : DO isec = 1, n_items
1377 3218 : NULLIFY (atm_names)
1378 3218 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1379 3218 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
1380 6436 : nonbonded%pot(start + isec)%pot%type = gp_type
1381 3218 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1382 3218 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1383 3218 : nonbonded%pot(start + isec)%pot%rcutsq = rcut*rcut
1384 3218 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1385 3218 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1386 : ! Parse the genpot info
1387 : CALL get_generic_info(section, "FUNCTION", nonbonded%pot(start + isec)%pot%set(1)%gp%potential, &
1388 : nonbonded%pot(start + isec)%pot%set(1)%gp%parameters, &
1389 : nonbonded%pot(start + isec)%pot%set(1)%gp%values, &
1390 3218 : size_variables=1, i_rep_sec=isec)
1391 3218 : nonbonded%pot(start + isec)%pot%set(1)%gp%variables = nonbonded%pot(start + isec)%pot%set(1)%gp%parameters(1)
1392 : !
1393 3218 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
1394 3218 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
1395 21 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
1396 3218 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
1397 3218 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
1398 10251 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
1399 : END DO
1400 576 : END SUBROUTINE read_gp_section
1401 :
1402 : ! **************************************************************************************************
1403 : !> \brief Reads the tersoff section
1404 : !> \param nonbonded ...
1405 : !> \param section ...
1406 : !> \param start ...
1407 : !> \param tersoff_section ...
1408 : !> \author ikuo
1409 : ! **************************************************************************************************
1410 40 : SUBROUTINE read_tersoff_section(nonbonded, section, start, tersoff_section)
1411 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1412 : TYPE(section_vals_type), POINTER :: section
1413 : INTEGER, INTENT(IN) :: start
1414 : TYPE(section_vals_type), POINTER :: tersoff_section
1415 :
1416 : CHARACTER(LEN=default_string_length), &
1417 40 : DIMENSION(:), POINTER :: atm_names
1418 : INTEGER :: isec, n_items, n_rep
1419 : REAL(KIND=dp) :: rcut, rcutsq
1420 :
1421 40 : CALL section_vals_get(section, n_repetition=n_items)
1422 84 : DO isec = 1, n_items
1423 44 : CALL cite_reference(Tersoff1988)
1424 44 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1425 :
1426 88 : nonbonded%pot(start + isec)%pot%type = tersoff_type
1427 44 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1428 44 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1429 44 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1430 44 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1431 :
1432 : CALL section_vals_val_get(tersoff_section, "A", i_rep_section=isec, &
1433 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%A)
1434 : CALL section_vals_val_get(tersoff_section, "B", i_rep_section=isec, &
1435 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%B)
1436 : CALL section_vals_val_get(tersoff_section, "lambda1", i_rep_section=isec, &
1437 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%lambda1)
1438 : CALL section_vals_val_get(tersoff_section, "lambda2", i_rep_section=isec, &
1439 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%lambda2)
1440 : CALL section_vals_val_get(tersoff_section, "alpha", i_rep_section=isec, &
1441 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%alpha)
1442 : CALL section_vals_val_get(tersoff_section, "beta", i_rep_section=isec, &
1443 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%beta)
1444 : CALL section_vals_val_get(tersoff_section, "n", i_rep_section=isec, &
1445 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%n)
1446 : CALL section_vals_val_get(tersoff_section, "c", i_rep_section=isec, &
1447 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%c)
1448 : CALL section_vals_val_get(tersoff_section, "d", i_rep_section=isec, &
1449 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%d)
1450 : CALL section_vals_val_get(tersoff_section, "h", i_rep_section=isec, &
1451 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%h)
1452 : CALL section_vals_val_get(tersoff_section, "lambda3", i_rep_section=isec, &
1453 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%lambda3)
1454 : CALL section_vals_val_get(tersoff_section, "bigR", i_rep_section=isec, &
1455 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%bigR)
1456 : CALL section_vals_val_get(tersoff_section, "bigD", i_rep_section=isec, &
1457 44 : r_val=nonbonded%pot(start + isec)%pot%set(1)%tersoff%bigD)
1458 :
1459 : rcutsq = (nonbonded%pot(start + isec)%pot%set(1)%tersoff%bigR + &
1460 44 : nonbonded%pot(start + isec)%pot%set(1)%tersoff%bigD)**2
1461 44 : nonbonded%pot(start + isec)%pot%set(1)%tersoff%rcutsq = rcutsq
1462 44 : nonbonded%pot(start + isec)%pot%rcutsq = rcutsq
1463 :
1464 : ! In case it is defined override the standard specification of RCUT
1465 44 : CALL section_vals_val_get(tersoff_section, "RCUT", i_rep_section=isec, n_rep_val=n_rep)
1466 84 : IF (n_rep == 1) THEN
1467 26 : CALL section_vals_val_get(tersoff_section, "RCUT", i_rep_section=isec, r_val=rcut)
1468 26 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1469 : END IF
1470 : END DO
1471 40 : END SUBROUTINE read_tersoff_section
1472 :
1473 : ! **************************************************************************************************
1474 : !> \brief Reads the gal19 section
1475 : !> \param nonbonded ...
1476 : !> \param section ...
1477 : !> \param start ...
1478 : !> \param gal_section ...
1479 : !> \author Clabaut Paul
1480 : ! **************************************************************************************************
1481 1 : SUBROUTINE read_gal_section(nonbonded, section, start, gal_section)
1482 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1483 : TYPE(section_vals_type), POINTER :: section
1484 : INTEGER, INTENT(IN) :: start
1485 : TYPE(section_vals_type), POINTER :: gal_section
1486 :
1487 : CHARACTER(LEN=default_string_length), &
1488 1 : DIMENSION(:), POINTER :: atm_names
1489 : INTEGER :: iatom, isec, n_items, n_rep, nval
1490 : LOGICAL :: is_ok
1491 : REAL(KIND=dp) :: rcut, rval
1492 1 : REAL(KIND=dp), DIMENSION(:), POINTER :: rvalues
1493 : TYPE(cp_sll_val_type), POINTER :: list
1494 : TYPE(section_vals_type), POINTER :: subsection
1495 : TYPE(val_type), POINTER :: val
1496 :
1497 1 : CALL section_vals_get(section, n_repetition=n_items)
1498 2 : DO isec = 1, n_items
1499 1 : CALL cite_reference(Clabaut2020)
1500 1 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1501 :
1502 2 : nonbonded%pot(start + isec)%pot%type = gal_type
1503 1 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1504 1 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1505 1 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1506 1 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1507 :
1508 1 : CALL section_vals_val_get(section, "METALS", i_rep_section=isec, c_vals=atm_names)
1509 1 : nonbonded%pot(start + isec)%pot%set(1)%gal%met1 = atm_names(1)
1510 1 : nonbonded%pot(start + isec)%pot%set(1)%gal%met2 = atm_names(2)
1511 :
1512 : CALL section_vals_val_get(gal_section, "epsilon", i_rep_section=isec, &
1513 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%epsilon)
1514 : CALL section_vals_val_get(gal_section, "bxy", i_rep_section=isec, &
1515 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%bxy)
1516 : CALL section_vals_val_get(gal_section, "bz", i_rep_section=isec, &
1517 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%bz)
1518 :
1519 1 : CALL section_vals_val_get(gal_section, "r", i_rep_section=isec, r_vals=rvalues)
1520 1 : nonbonded%pot(start + isec)%pot%set(1)%gal%r1 = rvalues(1)
1521 1 : nonbonded%pot(start + isec)%pot%set(1)%gal%r2 = rvalues(2)
1522 :
1523 : CALL section_vals_val_get(gal_section, "a1", i_rep_section=isec, &
1524 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%a1)
1525 : CALL section_vals_val_get(gal_section, "a2", i_rep_section=isec, &
1526 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%a2)
1527 : CALL section_vals_val_get(gal_section, "a3", i_rep_section=isec, &
1528 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%a3)
1529 : CALL section_vals_val_get(gal_section, "a4", i_rep_section=isec, &
1530 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%a4)
1531 : CALL section_vals_val_get(gal_section, "A", i_rep_section=isec, &
1532 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%a)
1533 : CALL section_vals_val_get(gal_section, "B", i_rep_section=isec, &
1534 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%b)
1535 : CALL section_vals_val_get(gal_section, "C", i_rep_section=isec, &
1536 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal%c)
1537 1 : NULLIFY (list)
1538 1 : subsection => section_vals_get_subs_vals(section, "GCN", i_rep_section=isec)
1539 1 : CALL section_vals_val_get(subsection, "_DEFAULT_KEYWORD_", n_rep_val=nval)
1540 3 : ALLOCATE (nonbonded%pot(start + isec)%pot%set(1)%gal%gcn(nval))
1541 1 : CALL section_vals_list_get(subsection, "_DEFAULT_KEYWORD_", list=list)
1542 871 : DO iatom = 1, nval
1543 : ! we use only the first default_string_length characters of each line
1544 870 : is_ok = cp_sll_val_next(list, val)
1545 870 : CALL val_get(val, r_val=rval)
1546 : ! assign values
1547 871 : nonbonded%pot(start + isec)%pot%set(1)%gal%gcn(iatom) = rval
1548 : END DO
1549 :
1550 : CALL section_vals_val_get(gal_section, "Fit_express", i_rep_section=isec, &
1551 1 : l_val=nonbonded%pot(start + isec)%pot%set(1)%gal%express)
1552 :
1553 : ! ! In case it is defined override the standard specification of RCUT
1554 1 : CALL section_vals_val_get(gal_section, "RCUT", i_rep_section=isec, n_rep_val=n_rep)
1555 3 : IF (n_rep == 1) THEN
1556 1 : CALL section_vals_val_get(gal_section, "RCUT", i_rep_section=isec, r_val=rcut)
1557 1 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1558 1 : nonbonded%pot(start + isec)%pot%set(1)%gal%rcutsq = rcut**2
1559 : END IF
1560 : END DO
1561 1 : END SUBROUTINE read_gal_section
1562 :
1563 : ! **************************************************************************************************
1564 : !> \brief Reads the gal21 section
1565 : !> \param nonbonded ...
1566 : !> \param section ...
1567 : !> \param start ...
1568 : !> \param gal21_section ...
1569 : !> \author Clabaut Paul
1570 : ! **************************************************************************************************
1571 1 : SUBROUTINE read_gal21_section(nonbonded, section, start, gal21_section)
1572 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1573 : TYPE(section_vals_type), POINTER :: section
1574 : INTEGER, INTENT(IN) :: start
1575 : TYPE(section_vals_type), POINTER :: gal21_section
1576 :
1577 : CHARACTER(LEN=default_string_length), &
1578 1 : DIMENSION(:), POINTER :: atm_names
1579 : INTEGER :: iatom, isec, n_items, n_rep, nval
1580 : LOGICAL :: is_ok
1581 : REAL(KIND=dp) :: rcut, rval
1582 1 : REAL(KIND=dp), DIMENSION(:), POINTER :: rvalues
1583 : TYPE(cp_sll_val_type), POINTER :: list
1584 : TYPE(section_vals_type), POINTER :: subsection
1585 : TYPE(val_type), POINTER :: val
1586 :
1587 1 : CALL section_vals_get(section, n_repetition=n_items)
1588 2 : DO isec = 1, n_items
1589 1 : CALL cite_reference(Clabaut2021)
1590 1 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1591 :
1592 2 : nonbonded%pot(start + isec)%pot%type = gal21_type
1593 1 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1594 1 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1595 1 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1596 1 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1597 :
1598 1 : CALL section_vals_val_get(section, "METALS", i_rep_section=isec, c_vals=atm_names)
1599 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%met1 = atm_names(1)
1600 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%met2 = atm_names(2)
1601 :
1602 1 : CALL section_vals_val_get(gal21_section, "epsilon", i_rep_section=isec, r_vals=rvalues)
1603 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%epsilon1 = rvalues(1)
1604 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%epsilon2 = rvalues(2)
1605 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%epsilon3 = rvalues(3)
1606 :
1607 1 : CALL section_vals_val_get(gal21_section, "bxy", i_rep_section=isec, r_vals=rvalues)
1608 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%bxy1 = rvalues(1)
1609 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%bxy2 = rvalues(2)
1610 :
1611 1 : CALL section_vals_val_get(gal21_section, "bz", i_rep_section=isec, r_vals=rvalues)
1612 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%bz1 = rvalues(1)
1613 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%bz2 = rvalues(2)
1614 :
1615 1 : CALL section_vals_val_get(gal21_section, "r", i_rep_section=isec, r_vals=rvalues)
1616 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%r1 = rvalues(1)
1617 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%r2 = rvalues(2)
1618 :
1619 1 : CALL section_vals_val_get(gal21_section, "a1", i_rep_section=isec, r_vals=rvalues)
1620 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a11 = rvalues(1)
1621 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a12 = rvalues(2)
1622 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a13 = rvalues(3)
1623 :
1624 1 : CALL section_vals_val_get(gal21_section, "a2", i_rep_section=isec, r_vals=rvalues)
1625 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a21 = rvalues(1)
1626 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a22 = rvalues(2)
1627 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a23 = rvalues(3)
1628 :
1629 1 : CALL section_vals_val_get(gal21_section, "a3", i_rep_section=isec, r_vals=rvalues)
1630 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a31 = rvalues(1)
1631 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a32 = rvalues(2)
1632 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a33 = rvalues(3)
1633 :
1634 1 : CALL section_vals_val_get(gal21_section, "a4", i_rep_section=isec, r_vals=rvalues)
1635 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a41 = rvalues(1)
1636 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a42 = rvalues(2)
1637 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%a43 = rvalues(3)
1638 :
1639 1 : CALL section_vals_val_get(gal21_section, "A", i_rep_section=isec, r_vals=rvalues)
1640 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%AO1 = rvalues(1)
1641 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%AO2 = rvalues(2)
1642 :
1643 1 : CALL section_vals_val_get(gal21_section, "B", i_rep_section=isec, r_vals=rvalues)
1644 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%BO1 = rvalues(1)
1645 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%BO2 = rvalues(2)
1646 :
1647 : CALL section_vals_val_get(gal21_section, "C", i_rep_section=isec, &
1648 1 : r_val=nonbonded%pot(start + isec)%pot%set(1)%gal21%c)
1649 :
1650 1 : CALL section_vals_val_get(gal21_section, "AH", i_rep_section=isec, r_vals=rvalues)
1651 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%AH1 = rvalues(1)
1652 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%AH2 = rvalues(2)
1653 :
1654 1 : CALL section_vals_val_get(gal21_section, "BH", i_rep_section=isec, r_vals=rvalues)
1655 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%BH1 = rvalues(1)
1656 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%BH2 = rvalues(2)
1657 :
1658 1 : NULLIFY (list)
1659 1 : subsection => section_vals_get_subs_vals(section, "GCN", i_rep_section=isec)
1660 1 : CALL section_vals_val_get(subsection, "_DEFAULT_KEYWORD_", n_rep_val=nval)
1661 3 : ALLOCATE (nonbonded%pot(start + isec)%pot%set(1)%gal21%gcn(nval))
1662 1 : CALL section_vals_list_get(subsection, "_DEFAULT_KEYWORD_", list=list)
1663 871 : DO iatom = 1, nval
1664 : ! we use only the first default_string_length characters of each line
1665 870 : is_ok = cp_sll_val_next(list, val)
1666 870 : CALL val_get(val, r_val=rval)
1667 : ! assign values
1668 871 : nonbonded%pot(start + isec)%pot%set(1)%gal21%gcn(iatom) = rval
1669 : END DO
1670 :
1671 : CALL section_vals_val_get(gal21_section, "Fit_express", i_rep_section=isec, &
1672 1 : l_val=nonbonded%pot(start + isec)%pot%set(1)%gal21%express)
1673 :
1674 : ! ! In case it is defined override the standard specification of RCUT
1675 1 : CALL section_vals_val_get(gal21_section, "RCUT", i_rep_section=isec, n_rep_val=n_rep)
1676 3 : IF (n_rep == 1) THEN
1677 1 : CALL section_vals_val_get(gal21_section, "RCUT", i_rep_section=isec, r_val=rcut)
1678 1 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1679 1 : nonbonded%pot(start + isec)%pot%set(1)%gal21%rcutsq = rcut**2
1680 : END IF
1681 : END DO
1682 1 : END SUBROUTINE read_gal21_section
1683 :
1684 : ! **************************************************************************************************
1685 : !> \brief Reads the siepmann section
1686 : !> \param nonbonded ...
1687 : !> \param section ...
1688 : !> \param start ...
1689 : !> \param siepmann_section ...
1690 : !> \author Dorothea Golze
1691 : ! **************************************************************************************************
1692 5 : SUBROUTINE read_siepmann_section(nonbonded, section, start, siepmann_section)
1693 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1694 : TYPE(section_vals_type), POINTER :: section
1695 : INTEGER, INTENT(IN) :: start
1696 : TYPE(section_vals_type), POINTER :: siepmann_section
1697 :
1698 : CHARACTER(LEN=default_string_length), &
1699 5 : DIMENSION(:), POINTER :: atm_names
1700 : INTEGER :: isec, n_items, n_rep
1701 : REAL(KIND=dp) :: rcut
1702 :
1703 5 : CALL section_vals_get(section, n_repetition=n_items)
1704 10 : DO isec = 1, n_items
1705 5 : CALL cite_reference(Siepmann1995)
1706 5 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1707 :
1708 10 : nonbonded%pot(start + isec)%pot%type = siepmann_type
1709 5 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1710 5 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1711 5 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1712 5 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1713 :
1714 : CALL section_vals_val_get(siepmann_section, "B", i_rep_section=isec, &
1715 5 : r_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%B)
1716 : CALL section_vals_val_get(siepmann_section, "D", i_rep_section=isec, &
1717 5 : r_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%D)
1718 : CALL section_vals_val_get(siepmann_section, "E", i_rep_section=isec, &
1719 5 : r_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%E)
1720 : CALL section_vals_val_get(siepmann_section, "F", i_rep_section=isec, &
1721 5 : r_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%F)
1722 : CALL section_vals_val_get(siepmann_section, "beta", i_rep_section=isec, &
1723 5 : r_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%beta)
1724 : CALL section_vals_val_get(siepmann_section, "ALLOW_OH_FORMATION", i_rep_section=isec, &
1725 5 : l_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%allow_oh_formation)
1726 : CALL section_vals_val_get(siepmann_section, "ALLOW_H3O_FORMATION", i_rep_section=isec, &
1727 5 : l_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%allow_h3o_formation)
1728 : CALL section_vals_val_get(siepmann_section, "ALLOW_O_FORMATION", i_rep_section=isec, &
1729 5 : l_val=nonbonded%pot(start + isec)%pot%set(1)%siepmann%allow_o_formation)
1730 :
1731 : ! ! In case it is defined override the standard specification of RCUT
1732 5 : CALL section_vals_val_get(siepmann_section, "RCUT", i_rep_section=isec, n_rep_val=n_rep)
1733 10 : IF (n_rep == 1) THEN
1734 5 : CALL section_vals_val_get(siepmann_section, "RCUT", i_rep_section=isec, r_val=rcut)
1735 5 : nonbonded%pot(start + isec)%pot%rcutsq = rcut**2
1736 5 : nonbonded%pot(start + isec)%pot%set(1)%siepmann%rcutsq = rcut**2
1737 : END IF
1738 : END DO
1739 5 : END SUBROUTINE read_siepmann_section
1740 :
1741 : ! **************************************************************************************************
1742 : !> \brief Reads the Buckingham plus Morse potential section
1743 : !> \param nonbonded ...
1744 : !> \param section ...
1745 : !> \param start ...
1746 : !> \author MI
1747 : ! **************************************************************************************************
1748 8 : SUBROUTINE read_bm_section(nonbonded, section, start)
1749 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1750 : TYPE(section_vals_type), POINTER :: section
1751 : INTEGER, INTENT(IN) :: start
1752 :
1753 : CHARACTER(LEN=default_string_length), &
1754 8 : DIMENSION(:), POINTER :: atm_names
1755 : INTEGER :: isec, n_items, n_rep
1756 : REAL(KIND=dp) :: a1, a2, b1, b2, beta, c, d, f0, r0, rcut
1757 :
1758 8 : CALL section_vals_get(section, n_repetition=n_items)
1759 26 : DO isec = 1, n_items
1760 18 : CALL cite_reference(Yamada2000)
1761 18 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1762 18 : CALL section_vals_val_get(section, "F0", i_rep_section=isec, r_val=f0)
1763 18 : CALL section_vals_val_get(section, "A1", i_rep_section=isec, r_val=a1)
1764 18 : CALL section_vals_val_get(section, "A2", i_rep_section=isec, r_val=a2)
1765 18 : CALL section_vals_val_get(section, "B1", i_rep_section=isec, r_val=b1)
1766 18 : CALL section_vals_val_get(section, "B2", i_rep_section=isec, r_val=b2)
1767 18 : CALL section_vals_val_get(section, "C", i_rep_section=isec, r_val=c)
1768 18 : CALL section_vals_val_get(section, "D", i_rep_section=isec, r_val=d)
1769 18 : CALL section_vals_val_get(section, "R0", i_rep_section=isec, r_val=r0)
1770 18 : CALL section_vals_val_get(section, "Beta", i_rep_section=isec, r_val=beta)
1771 18 : CALL section_vals_val_get(section, "RCUT", i_rep_section=isec, r_val=rcut)
1772 :
1773 36 : nonbonded%pot(start + isec)%pot%type = bm_type
1774 18 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1775 18 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1776 18 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1777 18 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1778 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%f0 = f0
1779 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%a1 = a1
1780 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%a2 = a2
1781 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%b1 = b1
1782 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%b2 = b2
1783 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%c = c
1784 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%d = d
1785 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%r0 = r0
1786 18 : nonbonded%pot(start + isec)%pot%set(1)%buckmo%beta = beta
1787 18 : nonbonded%pot(start + isec)%pot%rcutsq = rcut*rcut
1788 : !
1789 18 : CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, n_rep_val=n_rep)
1790 18 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMIN", i_rep_section=isec, &
1791 0 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmin)
1792 18 : CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, n_rep_val=n_rep)
1793 18 : IF (n_rep == 1) CALL section_vals_val_get(section, "RMAX", i_rep_section=isec, &
1794 62 : r_val=nonbonded%pot(start + isec)%pot%set(1)%rmax)
1795 : END DO
1796 8 : END SUBROUTINE read_bm_section
1797 :
1798 : ! **************************************************************************************************
1799 : !> \brief Reads the TABPOT section
1800 : !> \param nonbonded ...
1801 : !> \param section ...
1802 : !> \param start ...
1803 : !> \param para_env ...
1804 : !> \param mm_section ...
1805 : !> \author Alex Mironenko, Da Teng
1806 : ! **************************************************************************************************
1807 8 : SUBROUTINE read_tabpot_section(nonbonded, section, start, para_env, mm_section)
1808 : TYPE(pair_potential_p_type), POINTER :: nonbonded
1809 : TYPE(section_vals_type), POINTER :: section
1810 : INTEGER, INTENT(IN) :: start
1811 : TYPE(mp_para_env_type), POINTER :: para_env
1812 : TYPE(section_vals_type), POINTER :: mm_section
1813 :
1814 : CHARACTER(LEN=default_string_length), &
1815 8 : DIMENSION(:), POINTER :: atm_names
1816 : INTEGER :: isec, n_items
1817 :
1818 8 : CALL section_vals_get(section, n_repetition=n_items)
1819 32 : DO isec = 1, n_items
1820 24 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
1821 48 : nonbonded%pot(start + isec)%pot%type = tab_type
1822 24 : nonbonded%pot(start + isec)%pot%at1 = atm_names(1)
1823 24 : nonbonded%pot(start + isec)%pot%at2 = atm_names(2)
1824 24 : CALL uppercase(nonbonded%pot(start + isec)%pot%at1)
1825 24 : CALL uppercase(nonbonded%pot(start + isec)%pot%at2)
1826 : CALL section_vals_val_get(section, "PARM_FILE_NAME", i_rep_section=isec, &
1827 24 : c_val=nonbonded%pot(start + isec)%pot%set(1)%tab%tabpot_file_name)
1828 24 : CALL read_tabpot_data(nonbonded%pot(start + isec)%pot%set(1)%tab, para_env, mm_section)
1829 32 : nonbonded%pot(start + isec)%pot%set(1)%tab%index = isec
1830 : END DO
1831 8 : END SUBROUTINE read_tabpot_section
1832 :
1833 : ! **************************************************************************************************
1834 : !> \brief Reads the CHARGE section
1835 : !> \param charge_atm ...
1836 : !> \param charge ...
1837 : !> \param section ...
1838 : !> \param start ...
1839 : !> \author teo
1840 : ! **************************************************************************************************
1841 2109 : SUBROUTINE read_chrg_section(charge_atm, charge, section, start)
1842 : CHARACTER(LEN=default_string_length), &
1843 : DIMENSION(:), POINTER :: charge_atm
1844 : REAL(KIND=dp), DIMENSION(:), POINTER :: charge
1845 : TYPE(section_vals_type), POINTER :: section
1846 : INTEGER, INTENT(IN) :: start
1847 :
1848 : CHARACTER(LEN=default_string_length) :: atm_name
1849 : INTEGER :: isec, n_items
1850 :
1851 2109 : CALL section_vals_get(section, n_repetition=n_items)
1852 7266 : DO isec = 1, n_items
1853 5157 : CALL section_vals_val_get(section, "ATOM", i_rep_section=isec, c_val=atm_name)
1854 5157 : charge_atm(start + isec) = atm_name
1855 5157 : CALL uppercase(charge_atm(start + isec))
1856 7266 : CALL section_vals_val_get(section, "CHARGE", i_rep_section=isec, r_val=charge(start + isec))
1857 : END DO
1858 2109 : END SUBROUTINE read_chrg_section
1859 :
1860 : ! **************************************************************************************************
1861 : !> \brief Reads the POLARIZABILITY section
1862 : !> \param apol_atm ...
1863 : !> \param apol ...
1864 : !> \param damping_list ...
1865 : !> \param section ...
1866 : !> \param start ...
1867 : !> \author Marcel Baer
1868 : ! **************************************************************************************************
1869 34 : SUBROUTINE read_apol_section(apol_atm, apol, damping_list, section, &
1870 : start)
1871 : CHARACTER(LEN=default_string_length), &
1872 : DIMENSION(:), POINTER :: apol_atm
1873 : REAL(KIND=dp), DIMENSION(:), POINTER :: apol
1874 : TYPE(damping_info_type), DIMENSION(:), POINTER :: damping_list
1875 : TYPE(section_vals_type), POINTER :: section
1876 : INTEGER, INTENT(IN) :: start
1877 :
1878 : CHARACTER(LEN=default_string_length) :: atm_name
1879 : INTEGER :: isec, isec_damp, n_damp, n_items, &
1880 : start_damp, tmp_damp
1881 : TYPE(section_vals_type), POINTER :: tmp_section
1882 :
1883 34 : CALL section_vals_get(section, n_repetition=n_items)
1884 34 : NULLIFY (tmp_section)
1885 34 : n_damp = 0
1886 : ! *** Counts number of DIPOLE%DAMPING sections ****
1887 102 : DO isec = 1, n_items
1888 : tmp_section => section_vals_get_subs_vals(section, "DAMPING", &
1889 68 : i_rep_section=isec)
1890 68 : CALL section_vals_get(tmp_section, n_repetition=tmp_damp)
1891 102 : n_damp = n_damp + tmp_damp
1892 :
1893 : END DO
1894 :
1895 34 : IF (n_damp > 0) THEN
1896 42 : ALLOCATE (damping_list(1:n_damp))
1897 : END IF
1898 :
1899 : ! *** Reads DIPOLE sections *****
1900 34 : start_damp = 0
1901 102 : DO isec = 1, n_items
1902 68 : CALL section_vals_val_get(section, "ATOM", i_rep_section=isec, c_val=atm_name)
1903 68 : apol_atm(start + isec) = atm_name
1904 68 : CALL uppercase(apol_atm(start + isec))
1905 68 : CALL section_vals_val_get(section, "APOL", i_rep_section=isec, r_val=apol(start + isec))
1906 :
1907 : tmp_section => section_vals_get_subs_vals(section, "DAMPING", &
1908 68 : i_rep_section=isec)
1909 68 : CALL section_vals_get(tmp_section, n_repetition=tmp_damp)
1910 80 : DO isec_damp = 1, tmp_damp
1911 12 : damping_list(start_damp + isec_damp)%atm_name1 = apol_atm(start + isec)
1912 : CALL section_vals_val_get(tmp_section, "ATOM", i_rep_section=isec_damp, &
1913 12 : c_val=atm_name)
1914 12 : damping_list(start_damp + isec_damp)%atm_name2 = atm_name
1915 12 : CALL uppercase(damping_list(start_damp + isec_damp)%atm_name2)
1916 : CALL section_vals_val_get(tmp_section, "TYPE", i_rep_section=isec_damp, &
1917 12 : c_val=atm_name)
1918 12 : damping_list(start_damp + isec_damp)%dtype = atm_name
1919 12 : CALL uppercase(damping_list(start_damp + isec_damp)%dtype)
1920 :
1921 : CALL section_vals_val_get(tmp_section, "ORDER", i_rep_section=isec_damp, &
1922 12 : i_val=damping_list(start_damp + isec_damp)%order)
1923 : CALL section_vals_val_get(tmp_section, "BIJ", i_rep_section=isec_damp, &
1924 12 : r_val=damping_list(start_damp + isec_damp)%bij)
1925 : CALL section_vals_val_get(tmp_section, "CIJ", i_rep_section=isec_damp, &
1926 80 : r_val=damping_list(start_damp + isec_damp)%cij)
1927 : END DO
1928 170 : start_damp = start_damp + tmp_damp
1929 :
1930 : END DO
1931 :
1932 34 : END SUBROUTINE read_apol_section
1933 :
1934 : ! **************************************************************************************************
1935 : !> \brief Reads the QUADRUPOLE POLARIZABILITY section
1936 : !> \param cpol_atm ...
1937 : !> \param cpol ...
1938 : !> \param section ...
1939 : !> \param start ...
1940 : !> \author Marcel Baer
1941 : ! **************************************************************************************************
1942 0 : SUBROUTINE read_cpol_section(cpol_atm, cpol, section, start)
1943 : CHARACTER(LEN=default_string_length), &
1944 : DIMENSION(:), POINTER :: cpol_atm
1945 : REAL(KIND=dp), DIMENSION(:), POINTER :: cpol
1946 : TYPE(section_vals_type), POINTER :: section
1947 : INTEGER, INTENT(IN) :: start
1948 :
1949 : CHARACTER(LEN=default_string_length) :: atm_name
1950 : INTEGER :: isec, n_items
1951 :
1952 0 : CALL section_vals_get(section, n_repetition=n_items)
1953 0 : DO isec = 1, n_items
1954 0 : CALL section_vals_val_get(section, "ATOM", i_rep_section=isec, c_val=atm_name)
1955 0 : cpol_atm(start + isec) = atm_name
1956 0 : CALL uppercase(cpol_atm(start + isec))
1957 0 : CALL section_vals_val_get(section, "CPOL", i_rep_section=isec, r_val=cpol(start + isec))
1958 : END DO
1959 0 : END SUBROUTINE read_cpol_section
1960 :
1961 : ! **************************************************************************************************
1962 : !> \brief Reads the SHELL section
1963 : !> \param shell_list ...
1964 : !> \param section ...
1965 : !> \param start ...
1966 : !> \author Marcella Iannuzzi
1967 : ! **************************************************************************************************
1968 258 : SUBROUTINE read_shell_section(shell_list, section, start)
1969 :
1970 : TYPE(shell_p_type), DIMENSION(:), POINTER :: shell_list
1971 : TYPE(section_vals_type), POINTER :: section
1972 : INTEGER, INTENT(IN) :: start
1973 :
1974 : CHARACTER(LEN=default_string_length) :: atm_name
1975 : INTEGER :: i_rep, n_rep
1976 : REAL(dp) :: ccharge, cutoff, k, maxdist, mfrac, &
1977 : scharge
1978 :
1979 258 : CALL section_vals_get(section, n_repetition=n_rep)
1980 :
1981 710 : DO i_rep = 1, n_rep
1982 : CALL section_vals_val_get(section, "_SECTION_PARAMETERS_", &
1983 452 : c_val=atm_name, i_rep_section=i_rep)
1984 452 : CALL uppercase(atm_name)
1985 452 : shell_list(start + i_rep)%atm_name = atm_name
1986 452 : CALL section_vals_val_get(section, "CORE_CHARGE", i_rep_section=i_rep, r_val=ccharge)
1987 452 : shell_list(start + i_rep)%shell%charge_core = ccharge
1988 452 : CALL section_vals_val_get(section, "SHELL_CHARGE", i_rep_section=i_rep, r_val=scharge)
1989 452 : shell_list(start + i_rep)%shell%charge_shell = scharge
1990 452 : CALL section_vals_val_get(section, "MASS_FRACTION", i_rep_section=i_rep, r_val=mfrac)
1991 452 : shell_list(start + i_rep)%shell%massfrac = mfrac
1992 452 : CALL section_vals_val_get(section, "K2_SPRING", i_rep_section=i_rep, r_val=k)
1993 452 : IF (k < 0.0_dp) THEN
1994 : CALL cp_abort(__LOCATION__, &
1995 : "An invalid value was specified for the force constant k2 of the core-shell "// &
1996 0 : "spring potential")
1997 : END IF
1998 452 : shell_list(start + i_rep)%shell%k2_spring = k
1999 452 : CALL section_vals_val_get(section, "K4_SPRING", i_rep_section=i_rep, r_val=k)
2000 452 : IF (k < 0.0_dp) THEN
2001 : CALL cp_abort(__LOCATION__, &
2002 : "An invalid value was specified for the force constant k4 of the core-shell "// &
2003 0 : "spring potential")
2004 : END IF
2005 452 : shell_list(start + i_rep)%shell%k4_spring = k
2006 452 : CALL section_vals_val_get(section, "MAX_DISTANCE", i_rep_section=i_rep, r_val=maxdist)
2007 452 : shell_list(start + i_rep)%shell%max_dist = maxdist
2008 452 : CALL section_vals_val_get(section, "SHELL_CUTOFF", i_rep_section=i_rep, r_val=cutoff)
2009 1614 : shell_list(start + i_rep)%shell%shell_cutoff = cutoff
2010 : END DO
2011 :
2012 258 : END SUBROUTINE read_shell_section
2013 :
2014 : ! **************************************************************************************************
2015 : !> \brief Reads the BONDS section
2016 : !> \param bond_kind ...
2017 : !> \param bond_a ...
2018 : !> \param bond_b ...
2019 : !> \param bond_k ...
2020 : !> \param bond_r0 ...
2021 : !> \param bond_cs ...
2022 : !> \param section ...
2023 : !> \param start ...
2024 : !> \author teo
2025 : ! **************************************************************************************************
2026 975 : SUBROUTINE read_bonds_section(bond_kind, bond_a, bond_b, bond_k, bond_r0, bond_cs, section, start)
2027 : INTEGER, DIMENSION(:), POINTER :: bond_kind
2028 : CHARACTER(LEN=default_string_length), &
2029 : DIMENSION(:), POINTER :: bond_a, bond_b
2030 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: bond_k
2031 : REAL(KIND=dp), DIMENSION(:), POINTER :: bond_r0, bond_cs
2032 : TYPE(section_vals_type), POINTER :: section
2033 : INTEGER, INTENT(IN) :: start
2034 :
2035 : CHARACTER(LEN=default_string_length), &
2036 975 : DIMENSION(:), POINTER :: atm_names
2037 : INTEGER :: isec, k, n_items
2038 975 : REAL(KIND=dp), DIMENSION(:), POINTER :: Kvals
2039 :
2040 975 : NULLIFY (Kvals, atm_names)
2041 975 : CALL section_vals_get(section, n_repetition=n_items)
2042 2826 : DO isec = 1, n_items
2043 1851 : CALL section_vals_val_get(section, "KIND", i_rep_section=isec, i_val=bond_kind(start + isec))
2044 1851 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
2045 1851 : bond_a(start + isec) = atm_names(1)
2046 1851 : bond_b(start + isec) = atm_names(2)
2047 1851 : CALL uppercase(bond_a(start + isec))
2048 1851 : CALL uppercase(bond_b(start + isec))
2049 1851 : CALL section_vals_val_get(section, "K", i_rep_section=isec, r_vals=Kvals)
2050 1851 : CPASSERT(SIZE(Kvals) <= 3)
2051 7404 : bond_k(:, start + isec) = 0.0_dp
2052 3742 : DO k = 1, SIZE(Kvals)
2053 3742 : bond_k(k, start + isec) = Kvals(k)
2054 : END DO
2055 1851 : CALL section_vals_val_get(section, "R0", i_rep_section=isec, r_val=bond_r0(start + isec))
2056 2826 : CALL section_vals_val_get(section, "CS", i_rep_section=isec, r_val=bond_cs(start + isec))
2057 : END DO
2058 975 : END SUBROUTINE read_bonds_section
2059 :
2060 : ! **************************************************************************************************
2061 : !> \brief Reads the BENDS section
2062 : !> \param bend_kind ...
2063 : !> \param bend_a ...
2064 : !> \param bend_b ...
2065 : !> \param bend_c ...
2066 : !> \param bend_k ...
2067 : !> \param bend_theta0 ...
2068 : !> \param bend_cb ...
2069 : !> \param bend_r012 ...
2070 : !> \param bend_r032 ...
2071 : !> \param bend_kbs12 ...
2072 : !> \param bend_kbs32 ...
2073 : !> \param bend_kss ...
2074 : !> \param bend_legendre ...
2075 : !> \param section ...
2076 : !> \param start ...
2077 : !> \author teo
2078 : ! **************************************************************************************************
2079 939 : SUBROUTINE read_bends_section(bend_kind, bend_a, bend_b, bend_c, bend_k, bend_theta0, bend_cb, &
2080 : bend_r012, bend_r032, bend_kbs12, bend_kbs32, bend_kss, bend_legendre, &
2081 : section, start)
2082 : INTEGER, DIMENSION(:), POINTER :: bend_kind
2083 : CHARACTER(LEN=default_string_length), &
2084 : DIMENSION(:), POINTER :: bend_a, bend_b, bend_c
2085 : REAL(KIND=dp), DIMENSION(:), POINTER :: bend_k, bend_theta0, bend_cb, bend_r012, &
2086 : bend_r032, bend_kbs12, bend_kbs32, &
2087 : bend_kss
2088 : TYPE(legendre_data_type), DIMENSION(:), POINTER :: bend_legendre
2089 : TYPE(section_vals_type), POINTER :: section
2090 : INTEGER, INTENT(IN) :: start
2091 :
2092 : CHARACTER(LEN=default_string_length), &
2093 939 : DIMENSION(:), POINTER :: atm_names
2094 : INTEGER :: isec, k, n_items, n_rep
2095 939 : REAL(KIND=dp), DIMENSION(:), POINTER :: Kvals, r_values
2096 :
2097 939 : NULLIFY (Kvals, atm_names)
2098 939 : CALL section_vals_get(section, n_repetition=n_items)
2099 3058 : bend_legendre%order = 0
2100 3058 : DO isec = 1, n_items
2101 2119 : CALL section_vals_val_get(section, "KIND", i_rep_section=isec, i_val=bend_kind(start + isec))
2102 2119 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
2103 2119 : bend_a(start + isec) = atm_names(1)
2104 2119 : bend_b(start + isec) = atm_names(2)
2105 2119 : bend_c(start + isec) = atm_names(3)
2106 2119 : CALL uppercase(bend_a(start + isec))
2107 2119 : CALL uppercase(bend_b(start + isec))
2108 2119 : CALL uppercase(bend_c(start + isec))
2109 2119 : CALL section_vals_val_get(section, "K", i_rep_section=isec, r_vals=Kvals)
2110 2119 : CPASSERT(SIZE(Kvals) == 1)
2111 2119 : bend_k(start + isec) = Kvals(1)
2112 2119 : CALL section_vals_val_get(section, "THETA0", i_rep_section=isec, r_val=bend_theta0(start + isec))
2113 2119 : CALL section_vals_val_get(section, "CB", i_rep_section=isec, r_val=bend_cb(start + isec))
2114 2119 : CALL section_vals_val_get(section, "R012", i_rep_section=isec, r_val=bend_r012(start + isec))
2115 2119 : CALL section_vals_val_get(section, "R032", i_rep_section=isec, r_val=bend_r032(start + isec))
2116 2119 : CALL section_vals_val_get(section, "KBS12", i_rep_section=isec, r_val=bend_kbs12(start + isec))
2117 2119 : CALL section_vals_val_get(section, "KBS32", i_rep_section=isec, r_val=bend_kbs32(start + isec))
2118 2119 : CALL section_vals_val_get(section, "KSS", i_rep_section=isec, r_val=bend_kss(start + isec))
2119 : ! get legendre based data
2120 2119 : CALL section_vals_val_get(section, "LEGENDRE", i_rep_section=isec, n_rep_val=n_rep)
2121 5177 : DO k = 1, n_rep
2122 2119 : CALL section_vals_val_get(section, "LEGENDRE", i_rep_val=k, r_vals=r_values, i_rep_section=isec)
2123 2119 : bend_legendre(start + isec)%order = SIZE(r_values)
2124 2119 : IF (ASSOCIATED(bend_legendre(start + isec)%coeffs)) THEN
2125 0 : DEALLOCATE (bend_legendre(start + isec)%coeffs)
2126 : END IF
2127 6357 : ALLOCATE (bend_legendre(start + isec)%coeffs(bend_legendre(start + isec)%order))
2128 10675 : bend_legendre(start + isec)%coeffs = r_values
2129 : END DO
2130 : END DO
2131 939 : END SUBROUTINE read_bends_section
2132 :
2133 : ! **************************************************************************************************
2134 : !> \brief ...
2135 : !> \param ub_kind ...
2136 : !> \param ub_a ...
2137 : !> \param ub_b ...
2138 : !> \param ub_c ...
2139 : !> \param ub_k ...
2140 : !> \param ub_r0 ...
2141 : !> \param section ...
2142 : !> \param start ...
2143 : ! **************************************************************************************************
2144 939 : SUBROUTINE read_ubs_section(ub_kind, ub_a, ub_b, ub_c, ub_k, ub_r0, section, start)
2145 : INTEGER, DIMENSION(:), POINTER :: ub_kind
2146 : CHARACTER(LEN=default_string_length), &
2147 : DIMENSION(:), POINTER :: ub_a, ub_b, ub_c
2148 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: ub_k
2149 : REAL(KIND=dp), DIMENSION(:), POINTER :: ub_r0
2150 : TYPE(section_vals_type), POINTER :: section
2151 : INTEGER, INTENT(IN) :: start
2152 :
2153 : CHARACTER(LEN=default_string_length), &
2154 939 : DIMENSION(:), POINTER :: atm_names
2155 : INTEGER :: isec, k, n_items
2156 : LOGICAL :: explicit
2157 939 : REAL(KIND=dp), DIMENSION(:), POINTER :: Kvals
2158 : TYPE(section_vals_type), POINTER :: subsection
2159 :
2160 939 : NULLIFY (atm_names)
2161 939 : CALL section_vals_get(section, n_repetition=n_items)
2162 3058 : DO isec = 1, n_items
2163 2119 : subsection => section_vals_get_subs_vals(section, "UB", i_rep_section=isec)
2164 2119 : CALL section_vals_get(subsection, explicit=explicit)
2165 3058 : IF (explicit) THEN
2166 4 : CALL section_vals_val_get(subsection, "KIND", i_val=ub_kind(start + isec))
2167 4 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
2168 4 : ub_a(start + isec) = atm_names(1)
2169 4 : ub_b(start + isec) = atm_names(2)
2170 4 : ub_c(start + isec) = atm_names(3)
2171 4 : CALL uppercase(ub_a(start + isec))
2172 4 : CALL uppercase(ub_b(start + isec))
2173 4 : CALL uppercase(ub_c(start + isec))
2174 4 : CALL section_vals_val_get(subsection, "K", r_vals=Kvals)
2175 4 : CPASSERT(SIZE(Kvals) <= 3)
2176 16 : ub_k(:, start + isec) = 0.0_dp
2177 12 : DO k = 1, SIZE(Kvals)
2178 12 : ub_k(k, start + isec) = Kvals(k)
2179 : END DO
2180 4 : CALL section_vals_val_get(subsection, "R0", r_val=ub_r0(start + isec))
2181 : END IF
2182 : END DO
2183 939 : END SUBROUTINE read_ubs_section
2184 :
2185 : ! **************************************************************************************************
2186 : !> \brief Reads the TORSIONS section
2187 : !> \param torsion_kind ...
2188 : !> \param torsion_a ...
2189 : !> \param torsion_b ...
2190 : !> \param torsion_c ...
2191 : !> \param torsion_d ...
2192 : !> \param torsion_k ...
2193 : !> \param torsion_phi0 ...
2194 : !> \param torsion_m ...
2195 : !> \param section ...
2196 : !> \param start ...
2197 : !> \author teo
2198 : ! **************************************************************************************************
2199 6 : SUBROUTINE read_torsions_section(torsion_kind, torsion_a, torsion_b, torsion_c, torsion_d, torsion_k, &
2200 : torsion_phi0, torsion_m, section, start)
2201 : INTEGER, DIMENSION(:), POINTER :: torsion_kind
2202 : CHARACTER(LEN=default_string_length), &
2203 : DIMENSION(:), POINTER :: torsion_a, torsion_b, torsion_c, &
2204 : torsion_d
2205 : REAL(KIND=dp), DIMENSION(:), POINTER :: torsion_k, torsion_phi0
2206 : INTEGER, DIMENSION(:), POINTER :: torsion_m
2207 : TYPE(section_vals_type), POINTER :: section
2208 : INTEGER, INTENT(IN) :: start
2209 :
2210 : CHARACTER(LEN=default_string_length), &
2211 6 : DIMENSION(:), POINTER :: atm_names
2212 : INTEGER :: isec, n_items
2213 :
2214 6 : NULLIFY (atm_names)
2215 6 : CALL section_vals_get(section, n_repetition=n_items)
2216 44 : DO isec = 1, n_items
2217 38 : CALL section_vals_val_get(section, "KIND", i_rep_section=isec, i_val=torsion_kind(start + isec))
2218 38 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
2219 38 : torsion_a(start + isec) = atm_names(1)
2220 38 : torsion_b(start + isec) = atm_names(2)
2221 38 : torsion_c(start + isec) = atm_names(3)
2222 38 : torsion_d(start + isec) = atm_names(4)
2223 38 : CALL uppercase(torsion_a(start + isec))
2224 38 : CALL uppercase(torsion_b(start + isec))
2225 38 : CALL uppercase(torsion_c(start + isec))
2226 38 : CALL uppercase(torsion_d(start + isec))
2227 38 : CALL section_vals_val_get(section, "K", i_rep_section=isec, r_val=torsion_k(start + isec))
2228 38 : CALL section_vals_val_get(section, "PHI0", i_rep_section=isec, r_val=torsion_phi0(start + isec))
2229 38 : CALL section_vals_val_get(section, "M", i_rep_section=isec, i_val=torsion_m(start + isec))
2230 : ! Modify parameterisation for OPLS case
2231 44 : IF (torsion_kind(start + isec) .EQ. do_ff_opls) THEN
2232 12 : IF (torsion_phi0(start + isec) .NE. 0.0_dp) THEN
2233 : CALL cp_warn(__LOCATION__, "PHI0 parameter was non-zero "// &
2234 0 : "for an OPLS-type TORSION. It will be ignored.")
2235 : END IF
2236 12 : IF (MODULO(torsion_m(start + isec), 2) .EQ. 0) THEN
2237 : ! For even M, negate the cosine using a Pi phase factor
2238 2 : torsion_phi0(start + isec) = pi
2239 : END IF
2240 : ! the K parameter appears as K/2 in the OPLS parameterisation
2241 12 : torsion_k(start + isec) = torsion_k(start + isec)*0.5_dp
2242 : END IF
2243 : END DO
2244 6 : END SUBROUTINE read_torsions_section
2245 :
2246 : ! **************************************************************************************************
2247 : !> \brief Reads the IMPROPER section
2248 : !> \param impr_kind ...
2249 : !> \param impr_a ...
2250 : !> \param impr_b ...
2251 : !> \param impr_c ...
2252 : !> \param impr_d ...
2253 : !> \param impr_k ...
2254 : !> \param impr_phi0 ...
2255 : !> \param section ...
2256 : !> \param start ...
2257 : !> \author louis vanduyfhuys
2258 : ! **************************************************************************************************
2259 8 : SUBROUTINE read_improper_section(impr_kind, impr_a, impr_b, impr_c, impr_d, impr_k, &
2260 : impr_phi0, section, start)
2261 : INTEGER, DIMENSION(:), POINTER :: impr_kind
2262 : CHARACTER(LEN=default_string_length), &
2263 : DIMENSION(:), POINTER :: impr_a, impr_b, impr_c, impr_d
2264 : REAL(KIND=dp), DIMENSION(:), POINTER :: impr_k, impr_phi0
2265 : TYPE(section_vals_type), POINTER :: section
2266 : INTEGER, INTENT(IN) :: start
2267 :
2268 : CHARACTER(LEN=default_string_length), &
2269 8 : DIMENSION(:), POINTER :: atm_names
2270 : INTEGER :: isec, n_items
2271 :
2272 8 : NULLIFY (atm_names)
2273 8 : CALL section_vals_get(section, n_repetition=n_items)
2274 16 : DO isec = 1, n_items
2275 8 : CALL section_vals_val_get(section, "KIND", i_rep_section=isec, i_val=impr_kind(start + isec))
2276 8 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
2277 8 : impr_a(start + isec) = atm_names(1)
2278 8 : impr_b(start + isec) = atm_names(2)
2279 8 : impr_c(start + isec) = atm_names(3)
2280 8 : impr_d(start + isec) = atm_names(4)
2281 8 : CALL uppercase(impr_a(start + isec))
2282 8 : CALL uppercase(impr_b(start + isec))
2283 8 : CALL uppercase(impr_c(start + isec))
2284 8 : CALL uppercase(impr_d(start + isec))
2285 8 : CALL section_vals_val_get(section, "K", i_rep_section=isec, r_val=impr_k(start + isec))
2286 16 : CALL section_vals_val_get(section, "PHI0", i_rep_section=isec, r_val=impr_phi0(start + isec))
2287 : END DO
2288 8 : END SUBROUTINE read_improper_section
2289 :
2290 : ! **************************************************************************************************
2291 : !> \brief Reads the OPBEND section
2292 : !> \param opbend_kind ...
2293 : !> \param opbend_a ...
2294 : !> \param opbend_b ...
2295 : !> \param opbend_c ...
2296 : !> \param opbend_d ...
2297 : !> \param opbend_k ...
2298 : !> \param opbend_phi0 ...
2299 : !> \param section ...
2300 : !> \param start ...
2301 : !> \author louis vanduyfhuys
2302 : ! **************************************************************************************************
2303 2 : SUBROUTINE read_opbend_section(opbend_kind, opbend_a, opbend_b, opbend_c, opbend_d, opbend_k, &
2304 : opbend_phi0, section, start)
2305 : INTEGER, DIMENSION(:), POINTER :: opbend_kind
2306 : CHARACTER(LEN=default_string_length), &
2307 : DIMENSION(:), POINTER :: opbend_a, opbend_b, opbend_c, opbend_d
2308 : REAL(KIND=dp), DIMENSION(:), POINTER :: opbend_k, opbend_phi0
2309 : TYPE(section_vals_type), POINTER :: section
2310 : INTEGER, INTENT(IN) :: start
2311 :
2312 : CHARACTER(LEN=default_string_length), &
2313 2 : DIMENSION(:), POINTER :: atm_names
2314 : INTEGER :: isec, n_items
2315 :
2316 2 : NULLIFY (atm_names)
2317 2 : CALL section_vals_get(section, n_repetition=n_items)
2318 4 : DO isec = 1, n_items
2319 2 : CALL section_vals_val_get(section, "KIND", i_rep_section=isec, i_val=opbend_kind(start + isec))
2320 2 : CALL section_vals_val_get(section, "ATOMS", i_rep_section=isec, c_vals=atm_names)
2321 2 : opbend_a(start + isec) = atm_names(1)
2322 2 : opbend_b(start + isec) = atm_names(2)
2323 2 : opbend_c(start + isec) = atm_names(3)
2324 2 : opbend_d(start + isec) = atm_names(4)
2325 2 : CALL uppercase(opbend_a(start + isec))
2326 2 : CALL uppercase(opbend_b(start + isec))
2327 2 : CALL uppercase(opbend_c(start + isec))
2328 2 : CALL uppercase(opbend_d(start + isec))
2329 2 : CALL section_vals_val_get(section, "K", i_rep_section=isec, r_val=opbend_k(start + isec))
2330 4 : CALL section_vals_val_get(section, "PHI0", i_rep_section=isec, r_val=opbend_phi0(start + isec))
2331 : END DO
2332 2 : END SUBROUTINE read_opbend_section
2333 :
2334 : ! **************************************************************************************************
2335 : !> \brief Reads the force_field input section
2336 : !> \param ff_type ...
2337 : !> \param para_env ...
2338 : !> \param mm_section ...
2339 : !> \par History
2340 : !> JGH (30.11.2001) : moved determination of setup variables to
2341 : !> molecule_input
2342 : !> \author CJM
2343 : ! **************************************************************************************************
2344 2639 : SUBROUTINE read_force_field_section(ff_type, para_env, mm_section)
2345 : TYPE(force_field_type), INTENT(INOUT) :: ff_type
2346 : TYPE(mp_para_env_type), POINTER :: para_env
2347 : TYPE(section_vals_type), POINTER :: mm_section
2348 :
2349 : TYPE(section_vals_type), POINTER :: ff_section
2350 :
2351 : NULLIFY (ff_section)
2352 2639 : ff_section => section_vals_get_subs_vals(mm_section, "FORCEFIELD")
2353 2639 : CALL read_force_field_section1(ff_section, mm_section, ff_type, para_env)
2354 2639 : END SUBROUTINE read_force_field_section
2355 :
2356 : ! **************************************************************************************************
2357 : !> \brief reads EAM potential from library
2358 : !> \param eam ...
2359 : !> \param para_env ...
2360 : !> \param mm_section ...
2361 : ! **************************************************************************************************
2362 40 : SUBROUTINE read_eam_data(eam, para_env, mm_section)
2363 : TYPE(eam_pot_type), POINTER :: eam
2364 : TYPE(mp_para_env_type), POINTER :: para_env
2365 : TYPE(section_vals_type), POINTER :: mm_section
2366 :
2367 : CHARACTER(len=*), PARAMETER :: routineN = 'read_eam_data'
2368 :
2369 : INTEGER :: handle, i, iw
2370 : TYPE(cp_logger_type), POINTER :: logger
2371 : TYPE(cp_parser_type) :: parser
2372 :
2373 20 : CALL timeset(routineN, handle)
2374 20 : NULLIFY (logger)
2375 20 : logger => cp_get_default_logger()
2376 : iw = cp_print_key_unit_nr(logger, mm_section, "PRINT%FF_INFO", &
2377 20 : extension=".mmLog")
2378 20 : IF (iw > 0) WRITE (iw, *) "Reading EAM data from: ", TRIM(eam%eam_file_name)
2379 20 : CALL parser_create(parser, TRIM(eam%eam_file_name), para_env=para_env)
2380 :
2381 20 : CALL parser_get_next_line(parser, 1)
2382 20 : IF (iw > 0) WRITE (iw, *) "Title: ", parser%input_line
2383 :
2384 20 : CALL parser_get_next_line(parser, 2)
2385 20 : READ (parser%input_line, *) eam%drar, eam%drhoar, eam%acutal, eam%npoints
2386 20 : eam%drar = cp_unit_to_cp2k(eam%drar, "angstrom")
2387 20 : eam%acutal = cp_unit_to_cp2k(eam%acutal, "angstrom")
2388 : ! Relocating arrays with the right size
2389 20 : CALL reallocate(eam%rho, 1, eam%npoints)
2390 20 : CALL reallocate(eam%rhop, 1, eam%npoints)
2391 20 : CALL reallocate(eam%rval, 1, eam%npoints)
2392 20 : CALL reallocate(eam%rhoval, 1, eam%npoints)
2393 20 : CALL reallocate(eam%phi, 1, eam%npoints)
2394 20 : CALL reallocate(eam%phip, 1, eam%npoints)
2395 20 : CALL reallocate(eam%frho, 1, eam%npoints)
2396 20 : CALL reallocate(eam%frhop, 1, eam%npoints)
2397 : ! Reading density and derivative of density (with respect to r)
2398 64020 : DO i = 1, eam%npoints
2399 64000 : CALL parser_get_next_line(parser, 1)
2400 64000 : READ (parser%input_line, *) eam%rho(i), eam%rhop(i)
2401 64000 : eam%rhop(i) = cp_unit_to_cp2k(eam%rhop(i), "angstrom^-1")
2402 64000 : eam%rval(i) = REAL(i - 1, KIND=dp)*eam%drar
2403 64020 : eam%rhoval(i) = REAL(i - 1, KIND=dp)*eam%drhoar
2404 : END DO
2405 : ! Reading pair potential PHI and its derivative (with respect to r)
2406 64020 : DO i = 1, eam%npoints
2407 64000 : CALL parser_get_next_line(parser, 1)
2408 64000 : READ (parser%input_line, *) eam%phi(i), eam%phip(i)
2409 64000 : eam%phi(i) = cp_unit_to_cp2k(eam%phi(i), "eV")
2410 64020 : eam%phip(i) = cp_unit_to_cp2k(eam%phip(i), "eV*angstrom^-1")
2411 : END DO
2412 : ! Reading embedded function and its derivative (with respect to density)
2413 64020 : DO i = 1, eam%npoints
2414 64000 : CALL parser_get_next_line(parser, 1)
2415 64000 : READ (parser%input_line, *) eam%frho(i), eam%frhop(i)
2416 64000 : eam%frho(i) = cp_unit_to_cp2k(eam%frho(i), "eV")
2417 64020 : eam%frhop(i) = cp_unit_to_cp2k(eam%frhop(i), "eV")
2418 : END DO
2419 :
2420 20 : IF (iw > 0) WRITE (iw, *) "Finished EAM data"
2421 20 : CALL parser_release(parser)
2422 20 : CALL cp_print_key_finished_output(iw, logger, mm_section, "PRINT%FF_INFO")
2423 20 : CALL timestop(handle)
2424 :
2425 60 : END SUBROUTINE read_eam_data
2426 :
2427 : ! **************************************************************************************************
2428 : !> \brief reads NequIP potential from .pth file
2429 : !> \param nequip ...
2430 : !> \author Gabriele Tocci
2431 : ! **************************************************************************************************
2432 12 : SUBROUTINE read_nequip_data(nequip)
2433 : TYPE(nequip_pot_type), POINTER :: nequip
2434 :
2435 : CHARACTER(len=*), PARAMETER :: routineN = 'read_nequip_data'
2436 :
2437 : CHARACTER(LEN=default_path_length) :: allow_tf32_str, config_str, cutoff_str
2438 : INTEGER :: handle
2439 : LOGICAL :: allow_tf32, found
2440 :
2441 12 : CALL timeset(routineN, handle)
2442 :
2443 12 : INQUIRE (FILE=nequip%nequip_file_name, EXIST=found)
2444 12 : IF (.NOT. found) THEN
2445 : CALL cp_abort(__LOCATION__, &
2446 : "Nequip model file <"//TRIM(nequip%nequip_file_name)// &
2447 0 : "> not found.")
2448 : END IF
2449 :
2450 12 : nequip%nequip_version = torch_model_read_metadata(nequip%nequip_file_name, "nequip_version")
2451 12 : cutoff_str = torch_model_read_metadata(nequip%nequip_file_name, "r_max")
2452 12 : READ (cutoff_str, *) nequip%rcutsq
2453 12 : nequip%rcutsq = cp_unit_to_cp2k(nequip%rcutsq, nequip%unit_coords)
2454 12 : nequip%rcutsq = nequip%rcutsq*nequip%rcutsq
2455 12 : nequip%unit_coords_val = cp_unit_to_cp2k(nequip%unit_coords_val, nequip%unit_coords)
2456 12 : nequip%unit_forces_val = cp_unit_to_cp2k(nequip%unit_forces_val, nequip%unit_forces)
2457 12 : nequip%unit_energy_val = cp_unit_to_cp2k(nequip%unit_energy_val, nequip%unit_energy)
2458 12 : nequip%unit_cell_val = cp_unit_to_cp2k(nequip%unit_cell_val, nequip%unit_cell)
2459 :
2460 : ! look in config which contains all the .yaml file options to see if we use float32 or float64
2461 12 : config_str = torch_model_read_metadata(nequip%nequip_file_name, "config")
2462 12 : CALL read_default_dtype(config_str, nequip%do_nequip_sp)
2463 :
2464 12 : allow_tf32_str = torch_model_read_metadata(nequip%nequip_file_name, "allow_tf32")
2465 12 : allow_tf32 = (TRIM(allow_tf32_str) == "1")
2466 12 : IF (TRIM(allow_tf32_str) /= "1" .AND. TRIM(allow_tf32_str) /= "0") THEN
2467 : CALL cp_abort(__LOCATION__, &
2468 : "The value for allow_tf32 <"//TRIM(allow_tf32_str)// &
2469 0 : "> is not supported. Check the .yaml and .pth files.")
2470 : END IF
2471 12 : CALL torch_allow_tf32(allow_tf32)
2472 :
2473 12 : CALL timestop(handle)
2474 12 : END SUBROUTINE read_nequip_data
2475 :
2476 : ! **************************************************************************************************
2477 : !> \brief reads ALLEGRO potential from .pth file
2478 : !> \param allegro ...
2479 : !> \author Gabriele Tocci
2480 : ! **************************************************************************************************
2481 14 : SUBROUTINE read_allegro_data(allegro)
2482 : TYPE(allegro_pot_type), POINTER :: allegro
2483 :
2484 : CHARACTER(len=*), PARAMETER :: routineN = 'read_allegro_data'
2485 :
2486 : CHARACTER(LEN=default_path_length) :: allow_tf32_str, config_str, cutoff_str
2487 : INTEGER :: handle
2488 : LOGICAL :: allow_tf32, found
2489 :
2490 14 : CALL timeset(routineN, handle)
2491 :
2492 14 : INQUIRE (FILE=allegro%allegro_file_name, EXIST=found)
2493 14 : IF (.NOT. found) THEN
2494 : CALL cp_abort(__LOCATION__, &
2495 : "Allegro model file <"//TRIM(allegro%allegro_file_name)// &
2496 0 : "> not found.")
2497 : END IF
2498 :
2499 14 : allegro%nequip_version = torch_model_read_metadata(allegro%allegro_file_name, "nequip_version")
2500 14 : IF (allegro%nequip_version == "") THEN
2501 : CALL cp_abort(__LOCATION__, &
2502 : "Allegro model file <"//TRIM(allegro%allegro_file_name)// &
2503 0 : "> has not been deployed; did you forget to run `nequip-deploy`?")
2504 : END IF
2505 14 : cutoff_str = torch_model_read_metadata(allegro%allegro_file_name, "r_max")
2506 14 : READ (cutoff_str, *) allegro%rcutsq
2507 14 : allegro%rcutsq = cp_unit_to_cp2k(allegro%rcutsq, allegro%unit_coords)
2508 14 : allegro%rcutsq = allegro%rcutsq*allegro%rcutsq
2509 14 : allegro%unit_coords_val = cp_unit_to_cp2k(allegro%unit_coords_val, allegro%unit_coords)
2510 14 : allegro%unit_forces_val = cp_unit_to_cp2k(allegro%unit_forces_val, allegro%unit_forces)
2511 14 : allegro%unit_energy_val = cp_unit_to_cp2k(allegro%unit_energy_val, allegro%unit_energy)
2512 14 : allegro%unit_cell_val = cp_unit_to_cp2k(allegro%unit_cell_val, allegro%unit_cell)
2513 :
2514 : ! look in config which contains all the .yaml file options to see if we use float32 or float64
2515 14 : config_str = torch_model_read_metadata(allegro%allegro_file_name, "config")
2516 14 : CALL read_default_dtype(config_str, allegro%do_allegro_sp)
2517 :
2518 14 : allow_tf32_str = torch_model_read_metadata(allegro%allegro_file_name, "allow_tf32")
2519 14 : allow_tf32 = (TRIM(allow_tf32_str) == "1")
2520 14 : IF (TRIM(allow_tf32_str) /= "1" .AND. TRIM(allow_tf32_str) /= "0") THEN
2521 : CALL cp_abort(__LOCATION__, &
2522 : "The value for allow_tf32 <"//TRIM(allow_tf32_str)// &
2523 0 : "> is not supported. Check the .yaml and .pth files.")
2524 : END IF
2525 14 : CALL torch_allow_tf32(allow_tf32)
2526 :
2527 14 : CALL timestop(handle)
2528 14 : END SUBROUTINE read_allegro_data
2529 :
2530 : ! **************************************************************************************************
2531 : !> \brief reads the default_dtype used in the Allegro/NequIP model by parsing the config file
2532 : !> \param config_str ...
2533 : !> \param do_model_sp ...
2534 : !> \author Gabriele Tocci
2535 : ! **************************************************************************************************
2536 26 : SUBROUTINE read_default_dtype(config_str, do_model_sp)
2537 :
2538 : CHARACTER(LEN=default_path_length) :: config_str
2539 : LOGICAL :: do_model_sp
2540 :
2541 : CHARACTER(len=*), PARAMETER :: routineN = 'read_default_dtype'
2542 :
2543 : INTEGER :: handle, i, idx, len_config
2544 :
2545 26 : CALL timeset(routineN, handle)
2546 :
2547 26 : len_config = LEN_TRIM(config_str)
2548 26 : idx = INDEX(config_str, "default_dtype:")
2549 26 : IF (idx /= 0) THEN
2550 26 : i = idx + 14 ! skip over "default_dtype:"
2551 52 : DO WHILE (i <= len_config .AND. config_str(i:i) == " ")
2552 52 : i = i + 1 ! skip over any whitespace
2553 : END DO
2554 :
2555 26 : IF (i > len_config) THEN
2556 : CALL cp_abort(__LOCATION__, &
2557 : "No default_dtype found, check the Nequip/Allegro .yaml or .pth files."// &
2558 0 : " Default_dtype should be either <float32> or <float64>.")
2559 26 : ELSE IF (config_str(i:i + 6) == "float32") THEN
2560 8 : do_model_sp = .TRUE.
2561 18 : ELSE IF (config_str(i:i + 6) == "float64") THEN
2562 18 : do_model_sp = .FALSE.
2563 : ELSE
2564 : CALL cp_abort(__LOCATION__, &
2565 : "The default_dtype should be either <float32> or <float64>."// &
2566 0 : " Check the NequIP/Allegro .yaml and .pth files.")
2567 : END IF
2568 : END IF
2569 :
2570 26 : CALL timestop(handle)
2571 26 : END SUBROUTINE read_default_dtype
2572 :
2573 : ! **************************************************************************************************
2574 : !> \brief reads TABPOT potential from file
2575 : !> \param tab ...
2576 : !> \param para_env ...
2577 : !> \param mm_section ...
2578 : !> \author Da Teng, Alex Mironenko
2579 : ! **************************************************************************************************
2580 48 : SUBROUTINE read_tabpot_data(tab, para_env, mm_section)
2581 : TYPE(tab_pot_type), POINTER :: tab
2582 : TYPE(mp_para_env_type), POINTER :: para_env
2583 : TYPE(section_vals_type), POINTER :: mm_section
2584 :
2585 : CHARACTER(len=*), PARAMETER :: routineN = 'read_tabpot_data'
2586 :
2587 : CHARACTER :: d1, d2
2588 : INTEGER :: d, handle, i, iw
2589 : TYPE(cp_logger_type), POINTER :: logger
2590 : TYPE(cp_parser_type) :: parser
2591 :
2592 24 : CALL timeset(routineN, handle)
2593 24 : NULLIFY (logger)
2594 24 : logger => cp_get_default_logger()
2595 : iw = cp_print_key_unit_nr(logger, mm_section, "PRINT%FF_INFO", &
2596 24 : extension=".mmLog")
2597 24 : IF (iw > 0) WRITE (iw, *) "Reading TABPOT data from: ", TRIM(tab%tabpot_file_name)
2598 24 : CALL parser_create(parser, TRIM(tab%tabpot_file_name), para_env=para_env)
2599 24 : CALL parser_get_next_line(parser, 1)
2600 24 : IF (iw > 0) WRITE (iw, *) "Title: ", parser%input_line
2601 24 : CALL parser_get_next_line(parser, 1)
2602 :
2603 : ! example format: N 1000 R 1.00 20.0
2604 : ! Assume the data is evenly spaced
2605 24 : READ (parser%input_line, *) d1, tab%npoints, d2, tab%dr, tab%rcut
2606 :
2607 : ! Relocating arrays with the right size
2608 24 : CALL reallocate(tab%r, 1, tab%npoints)
2609 24 : CALL reallocate(tab%e, 1, tab%npoints)
2610 24 : CALL reallocate(tab%f, 1, tab%npoints)
2611 :
2612 : ! Reading r, e, f
2613 21912 : DO i = 1, tab%npoints
2614 21888 : CALL parser_get_next_line(parser, 1)
2615 21888 : READ (parser%input_line, *) d, tab%r(i), tab%e(i), tab%f(i)
2616 21888 : tab%r(i) = cp_unit_to_cp2k(tab%r(i), "angstrom")
2617 21888 : tab%e(i) = cp_unit_to_cp2k(tab%e(i), "kcalmol")
2618 21912 : tab%f(i) = cp_unit_to_cp2k(tab%f(i), "kcalmol*angstrom^-1")
2619 : END DO
2620 :
2621 24 : tab%dr = tab%r(2) - tab%r(1)
2622 24 : tab%rcut = cp_unit_to_cp2k(tab%rcut, "angstrom")
2623 :
2624 24 : IF (iw > 0) WRITE (iw, *) "Finished TABPOT data"
2625 24 : CALL parser_release(parser)
2626 24 : CALL cp_print_key_finished_output(iw, logger, mm_section, "PRINT%FF_INFO")
2627 24 : CALL timestop(handle)
2628 72 : END SUBROUTINE read_tabpot_data
2629 :
2630 : END MODULE force_fields_input
|
