Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : !***************************************************************************************************
9 : !> \brief Interface to the SIRIUS Library
10 : !> \par History
11 : !> 07.2018 initial create
12 : !> \author JHU
13 : !***************************************************************************************************
14 : #if defined(__SIRIUS)
15 : MODULE sirius_interface
16 : USE ISO_C_BINDING, ONLY: C_DOUBLE,&
17 : C_INT,&
18 : C_LOC
19 : USE atom_kind_orbitals, ONLY: calculate_atomic_orbitals,&
20 : gth_potential_conversion
21 : USE atom_types, ONLY: atom_gthpot_type
22 : USE atom_upf, ONLY: atom_upfpot_type
23 : USE atom_utils, ONLY: atom_local_potential
24 : USE atomic_kind_types, ONLY: atomic_kind_type,&
25 : get_atomic_kind
26 : USE cell_types, ONLY: cell_type,&
27 : real_to_scaled
28 : USE cp_log_handling, ONLY: cp_get_default_logger,&
29 : cp_logger_get_default_io_unit,&
30 : cp_logger_type
31 : USE cp_output_handling, ONLY: cp_p_file,&
32 : cp_print_key_finished_output,&
33 : cp_print_key_should_output,&
34 : cp_print_key_unit_nr
35 : USE external_potential_types, ONLY: gth_potential_type
36 : USE input_constants, ONLY: do_gapw_log
37 : USE input_cp2k_pwdft, ONLY: SIRIUS_FUNC_VDWDF,&
38 : SIRIUS_FUNC_VDWDF2,&
39 : SIRIUS_FUNC_VDWDFCX
40 : USE input_section_types, ONLY: section_vals_get,&
41 : section_vals_get_subs_vals,&
42 : section_vals_get_subs_vals2,&
43 : section_vals_type,&
44 : section_vals_val_get
45 : USE kinds, ONLY: default_string_length,&
46 : dp
47 : USE machine, ONLY: m_flush
48 : USE mathconstants, ONLY: fourpi,&
49 : gamma1
50 : USE message_passing, ONLY: mp_para_env_type
51 : USE particle_types, ONLY: particle_type
52 : USE physcon, ONLY: massunit
53 : USE pwdft_environment_types, ONLY: pwdft_energy_type,&
54 : pwdft_env_get,&
55 : pwdft_env_set,&
56 : pwdft_environment_type
57 : USE qs_grid_atom, ONLY: allocate_grid_atom,&
58 : create_grid_atom,&
59 : deallocate_grid_atom,&
60 : grid_atom_type
61 : USE qs_kind_types, ONLY: get_qs_kind,&
62 : qs_kind_type
63 : USE qs_subsys_types, ONLY: qs_subsys_get,&
64 : qs_subsys_type
65 : USE sirius, ONLY: &
66 : SIRIUS_INTEGER_ARRAY_TYPE, SIRIUS_INTEGER_TYPE, SIRIUS_LOGICAL_ARRAY_TYPE, &
67 : SIRIUS_LOGICAL_TYPE, SIRIUS_NUMBER_ARRAY_TYPE, SIRIUS_NUMBER_TYPE, &
68 : SIRIUS_STRING_ARRAY_TYPE, SIRIUS_STRING_TYPE, sirius_add_atom, sirius_add_atom_type, &
69 : sirius_add_atom_type_radial_function, sirius_add_xc_functional, sirius_context_handler, &
70 : sirius_create_context, sirius_create_ground_state, sirius_create_kset_from_grid, &
71 : sirius_finalize, sirius_find_ground_state, sirius_get_band_energies, &
72 : sirius_get_band_occupancies, sirius_get_energy, sirius_get_forces, &
73 : sirius_get_kpoint_properties, sirius_get_num_kpoints, sirius_get_parameters, &
74 : sirius_get_stress_tensor, sirius_ground_state_handler, sirius_import_parameters, &
75 : sirius_initialize, sirius_initialize_context, sirius_kpoint_set_handler, &
76 : sirius_option_get_info, sirius_option_get_section_length, sirius_option_set, &
77 : sirius_set_atom_position, sirius_set_atom_type_dion, sirius_set_atom_type_hubbard, &
78 : sirius_set_atom_type_radial_grid, sirius_set_lattice_vectors, sirius_set_mpi_grid_dims, &
79 : sirius_update_ground_state
80 : #include "./base/base_uses.f90"
81 :
82 : IMPLICIT NONE
83 :
84 : PRIVATE
85 :
86 : ! *** Global parameters ***
87 :
88 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'sirius_interface'
89 :
90 : ! *** Public subroutines ***
91 :
92 : PUBLIC :: cp_sirius_init, cp_sirius_finalize
93 : PUBLIC :: cp_sirius_create_env, cp_sirius_energy_force, cp_sirius_update_context
94 :
95 : CONTAINS
96 :
97 : !***************************************************************************************************
98 : !> \brief ...
99 : !> \param
100 : !> \par History
101 : !> 07.2018 start the Sirius library
102 : !> \author JHU
103 : ! **************************************************************************************************
104 9238 : SUBROUTINE cp_sirius_init()
105 9238 : CALL sirius_initialize(.FALSE.)
106 9238 : END SUBROUTINE cp_sirius_init
107 :
108 : !***************************************************************************************************
109 : !> \brief ...
110 : !> \param
111 : !> \par History
112 : !> 07.2018 stop the Sirius library
113 : !> \author JHU
114 : ! **************************************************************************************************
115 9238 : SUBROUTINE cp_sirius_finalize()
116 9238 : CALL sirius_finalize(.FALSE., .FALSE., .FALSE.)
117 9238 : END SUBROUTINE cp_sirius_finalize
118 :
119 : !***************************************************************************************************
120 : !> \brief ...
121 : !> \param pwdft_env ...
122 : !> \param
123 : !> \par History
124 : !> 07.2018 Create the Sirius environment
125 : !> \author JHU
126 : ! **************************************************************************************************
127 16 : SUBROUTINE cp_sirius_create_env(pwdft_env)
128 : TYPE(pwdft_environment_type), POINTER :: pwdft_env
129 : #if defined(__SIRIUS)
130 :
131 : CHARACTER(len=2) :: element_symbol
132 : CHARACTER(len=default_string_length) :: label
133 : INTEGER :: i, iatom, ibeta, ifun, ikind, iwf, j, l, &
134 : n, natom, nbeta, nkind, nmesh, &
135 : num_mag_dims, sirius_mpi_comm, vdw_func, nu, lu, output_unit
136 16 : INTEGER, DIMENSION(:), POINTER :: mpi_grid_dims
137 : INTEGER(KIND=C_INT), DIMENSION(3) :: k_grid, k_shift
138 16 : INTEGER, DIMENSION(:), POINTER :: kk
139 : LOGICAL :: up, use_ref_cell
140 : LOGICAL(4) :: use_so, use_symmetry, dft_plus_u_atom
141 16 : REAL(KIND=C_DOUBLE), ALLOCATABLE, DIMENSION(:) :: fun
142 16 : REAL(KIND=C_DOUBLE), ALLOCATABLE, DIMENSION(:, :) :: dion
143 : REAL(KIND=C_DOUBLE), DIMENSION(3) :: a1, a2, a3, v1, v2
144 : REAL(KIND=dp) :: al, angle1, angle2, cval, focc, &
145 : magnetization, mass, pf, rl, zeff, alpha_u, beta_u, &
146 : J0_u, J_u, U_u, occ_u, u_minus_J
147 16 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: beta, corden, ef, fe, locpot, rc, rp
148 : REAL(KIND=dp), DIMENSION(3) :: vr, vs, j_t
149 16 : REAL(KIND=dp), DIMENSION(:), POINTER :: density
150 16 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: wavefunction, wfninfo
151 : TYPE(atom_gthpot_type), POINTER :: gth_atompot
152 : TYPE(atom_upfpot_type), POINTER :: upf_pot
153 16 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
154 : TYPE(atomic_kind_type), POINTER :: atomic_kind
155 : TYPE(cell_type), POINTER :: my_cell
156 : TYPE(mp_para_env_type), POINTER :: para_env
157 : TYPE(grid_atom_type), POINTER :: atom_grid
158 : TYPE(gth_potential_type), POINTER :: gth_potential
159 16 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
160 16 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
161 : TYPE(qs_subsys_type), POINTER :: qs_subsys
162 : TYPE(section_vals_type), POINTER :: pwdft_section, pwdft_sub_section, &
163 : xc_fun, xc_section
164 : TYPE(sirius_context_handler) :: sctx
165 : TYPE(sirius_ground_state_handler) :: gs_handler
166 : TYPE(sirius_kpoint_set_handler) :: ks_handler
167 :
168 0 : CPASSERT(ASSOCIATED(pwdft_env))
169 :
170 16 : output_unit = cp_logger_get_default_io_unit()
171 : ! create context of simulation
172 16 : CALL pwdft_env_get(pwdft_env, para_env=para_env)
173 16 : sirius_mpi_comm = para_env%get_handle()
174 16 : CALL sirius_create_context(sirius_mpi_comm, sctx)
175 :
176 : ! the "fun" starts.
177 :
178 16 : CALL pwdft_env_get(pwdft_env=pwdft_env, pwdft_input=pwdft_section, xc_input=xc_section)
179 :
180 : CALL section_vals_val_get(pwdft_section, "ignore_convergence_failure", &
181 16 : l_val=pwdft_env%ignore_convergence_failure)
182 : ! cp2k should *have* a function that return all xc_functionals. Doing
183 : ! manually is prone to errors
184 :
185 16 : IF (ASSOCIATED(xc_section)) THEN
186 16 : ifun = 0
187 30 : DO
188 46 : ifun = ifun + 1
189 46 : xc_fun => section_vals_get_subs_vals2(xc_section, i_section=ifun)
190 46 : IF (.NOT. ASSOCIATED(xc_fun)) EXIT
191 : ! Here, we do not have to check whether the functional name starts with XC_
192 : ! because we only allow the shorter form w/o XC_
193 46 : CALL sirius_add_xc_functional(sctx, "XC_"//TRIM(xc_fun%section%name))
194 : END DO
195 : END IF
196 :
197 : ! import control section
198 16 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "control")
199 16 : IF (ASSOCIATED(pwdft_sub_section)) THEN
200 16 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "control")
201 16 : CALL section_vals_val_get(pwdft_sub_section, "mpi_grid_dims", i_vals=mpi_grid_dims)
202 : END IF
203 :
204 : ! import parameters section
205 16 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "parameters")
206 :
207 16 : IF (ASSOCIATED(pwdft_sub_section)) THEN
208 16 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "parameters")
209 16 : CALL section_vals_val_get(pwdft_sub_section, "ngridk", i_vals=kk)
210 16 : k_grid(1) = kk(1)
211 16 : k_grid(2) = kk(2)
212 16 : k_grid(3) = kk(3)
213 :
214 16 : CALL section_vals_val_get(pwdft_sub_section, "shiftk", i_vals=kk)
215 16 : k_shift(1) = kk(1)
216 16 : k_shift(2) = kk(2)
217 16 : k_shift(3) = kk(3)
218 16 : CALL section_vals_val_get(pwdft_sub_section, "num_mag_dims", i_val=num_mag_dims)
219 16 : CALL section_vals_val_get(pwdft_sub_section, "use_symmetry", l_val=use_symmetry)
220 16 : CALL section_vals_val_get(pwdft_sub_section, "so_correction", l_val=use_so)
221 :
222 : ! now check if van der walls corrections are needed
223 : vdw_func = -1
224 : #ifdef __LIBVDWXC
225 16 : CALL section_vals_val_get(pwdft_sub_section, "vdw_functional", i_val=vdw_func)
226 0 : SELECT CASE (vdw_func)
227 : CASE (SIRIUS_FUNC_VDWDF)
228 0 : CALL sirius_add_xc_functional(sctx, "XC_FUNC_VDWDF")
229 : CASE (SIRIUS_FUNC_VDWDF2)
230 0 : CALL sirius_add_xc_functional(sctx, "XC_FUNC_VDWDF2")
231 : CASE (SIRIUS_FUNC_VDWDFCX)
232 16 : CALL sirius_add_xc_functional(sctx, "XC_FUNC_VDWDF2")
233 : CASE default
234 : END SELECT
235 : #endif
236 :
237 : END IF
238 :
239 : ! import mixer section
240 16 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "mixer")
241 16 : IF (ASSOCIATED(pwdft_sub_section)) THEN
242 16 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "mixer")
243 : END IF
244 :
245 : ! import settings section
246 16 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "settings")
247 16 : IF (ASSOCIATED(pwdft_sub_section)) THEN
248 16 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "settings")
249 : END IF
250 :
251 : ! import solver section
252 16 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "iterative_solver")
253 16 : IF (ASSOCIATED(pwdft_sub_section)) THEN
254 16 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "iterative_solver")
255 : END IF
256 :
257 : #if defined(__SIRIUS_DFTD4)
258 : ! import dftd3 and dftd4 section
259 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "dftd4")
260 : IF (ASSOCIATED(pwdft_sub_section)) THEN
261 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "dftd4")
262 : END IF
263 :
264 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "dftd3")
265 : IF (ASSOCIATED(pwdft_sub_section)) THEN
266 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "dftd3")
267 : END IF
268 : #endif
269 :
270 : !
271 : ! uncomment these lines when nlcg is officially supported
272 : !
273 : #if defined(__SIRIUS_NLCG)
274 : ! import nlcg section
275 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "nlcg")
276 : IF (ASSOCIATED(pwdft_sub_section)) THEN
277 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "nlcg")
278 : END IF
279 : #endif
280 :
281 : #if defined(__SIRIUS_VCSQNM)
282 : pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "vcsqnm")
283 : IF (ASSOCIATED(pwdft_sub_section)) THEN
284 : CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "vcsqnm")
285 : END IF
286 : #endif
287 :
288 : !CALL sirius_dump_runtime_setup(sctx, "runtime.json")
289 16 : CALL sirius_import_parameters(sctx, '{}')
290 :
291 : ! lattice vectors of the unit cell should be in [a.u.] (length is in [a.u.])
292 16 : CALL pwdft_env_get(pwdft_env=pwdft_env, qs_subsys=qs_subsys)
293 16 : CALL qs_subsys_get(qs_subsys, cell=my_cell, use_ref_cell=use_ref_cell)
294 64 : a1(:) = my_cell%hmat(:, 1)
295 64 : a2(:) = my_cell%hmat(:, 2)
296 64 : a3(:) = my_cell%hmat(:, 3)
297 16 : CALL sirius_set_lattice_vectors(sctx, a1(1), a2(1), a3(1))
298 :
299 16 : IF (use_ref_cell) THEN
300 0 : CPWARN("SIRIUS| The specified CELL_REF will be ignored for PW_DFT calculations")
301 : END IF
302 :
303 : ! set up the atomic type definitions
304 : CALL qs_subsys_get(qs_subsys, &
305 : atomic_kind_set=atomic_kind_set, &
306 : qs_kind_set=qs_kind_set, &
307 16 : particle_set=particle_set)
308 16 : nkind = SIZE(atomic_kind_set)
309 38 : DO ikind = 1, nkind
310 : CALL get_atomic_kind(atomic_kind_set(ikind), &
311 22 : name=label, element_symbol=element_symbol, mass=mass)
312 22 : CALL get_qs_kind(qs_kind_set(ikind), zeff=zeff)
313 22 : NULLIFY (upf_pot, gth_potential)
314 22 : CALL get_qs_kind(qs_kind_set(ikind), upf_potential=upf_pot, gth_potential=gth_potential)
315 :
316 22 : IF (ASSOCIATED(upf_pot)) THEN
317 : CALL sirius_add_atom_type(sctx, label, fname=upf_pot%filename, &
318 : symbol=element_symbol, &
319 20 : mass=REAL(mass/massunit, KIND=C_DOUBLE))
320 :
321 2 : ELSEIF (ASSOCIATED(gth_potential)) THEN
322 : !
323 2 : NULLIFY (atom_grid)
324 2 : CALL allocate_grid_atom(atom_grid)
325 2 : nmesh = 929
326 2 : atom_grid%nr = nmesh
327 2 : CALL create_grid_atom(atom_grid, nmesh, 1, 1, 0, do_gapw_log)
328 8 : ALLOCATE (rp(nmesh), fun(nmesh))
329 2 : IF (atom_grid%rad(1) < atom_grid%rad(nmesh)) THEN
330 : up = .TRUE.
331 : ELSE
332 : up = .FALSE.
333 : END IF
334 : IF (up) THEN
335 0 : rp(1:nmesh) = atom_grid%rad(1:nmesh)
336 : ELSE
337 1860 : DO i = 1, nmesh
338 1860 : rp(i) = atom_grid%rad(nmesh - i + 1)
339 : END DO
340 : END IF
341 : ! add new atom type
342 : CALL sirius_add_atom_type(sctx, label, &
343 : zn=NINT(zeff + 0.001d0), &
344 : symbol=element_symbol, &
345 : mass=REAL(mass/massunit, KIND=C_DOUBLE), &
346 2 : spin_orbit=.FALSE.)
347 : !
348 720 : ALLOCATE (gth_atompot)
349 2 : CALL gth_potential_conversion(gth_potential, gth_atompot)
350 : ! set radial grid
351 1860 : fun(1:nmesh) = rp(1:nmesh)
352 2 : CALL sirius_set_atom_type_radial_grid(sctx, label, nmesh, fun(1))
353 : ! set beta-projectors
354 8 : ALLOCATE (ef(nmesh), beta(nmesh))
355 2 : ibeta = 0
356 10 : DO l = 0, 3
357 8 : IF (gth_atompot%nl(l) == 0) CYCLE
358 2 : rl = gth_atompot%rcnl(l)
359 : ! we need to multiply by r so that data transferred to sirius are r \beta(r) not beta(r)
360 1860 : ef(1:nmesh) = EXP(-0.5_dp*rp(1:nmesh)*rp(1:nmesh)/(rl*rl))
361 6 : DO i = 1, gth_atompot%nl(l)
362 2 : pf = rl**(l + 0.5_dp*(4._dp*i - 1._dp))
363 2 : j = l + 2*i - 1
364 2 : pf = SQRT(2._dp)/(pf*SQRT(gamma1(j)))
365 1860 : beta(:) = pf*rp**(l + 2*i - 2)*ef
366 2 : ibeta = ibeta + 1
367 1860 : fun(1:nmesh) = beta(1:nmesh)*rp(1:nmesh)
368 : CALL sirius_add_atom_type_radial_function(sctx, label, &
369 10 : "beta", fun(1), nmesh, l=l)
370 : END DO
371 : END DO
372 2 : DEALLOCATE (ef, beta)
373 2 : nbeta = ibeta
374 :
375 : ! nonlocal PP matrix elements
376 8 : ALLOCATE (dion(nbeta, nbeta))
377 6 : dion = 0.0_dp
378 10 : DO l = 0, 3
379 8 : IF (gth_atompot%nl(l) == 0) CYCLE
380 2 : ibeta = SUM(gth_atompot%nl(0:l - 1)) + 1
381 2 : i = ibeta + gth_atompot%nl(l) - 1
382 8 : dion(ibeta:i, ibeta:i) = gth_atompot%hnl(1:gth_atompot%nl(l), 1:gth_atompot%nl(l), l)
383 : END DO
384 2 : CALL sirius_set_atom_type_dion(sctx, label, nbeta, dion(1, 1))
385 2 : DEALLOCATE (dion)
386 :
387 : ! set non-linear core correction
388 2 : IF (gth_atompot%nlcc) THEN
389 0 : ALLOCATE (corden(nmesh), fe(nmesh), rc(nmesh))
390 0 : corden(:) = 0.0_dp
391 0 : n = gth_atompot%nexp_nlcc
392 0 : DO i = 1, n
393 0 : al = gth_atompot%alpha_nlcc(i)
394 0 : rc(:) = rp(:)/al
395 0 : fe(:) = EXP(-0.5_dp*rc(:)*rc(:))
396 0 : DO j = 1, gth_atompot%nct_nlcc(i)
397 0 : cval = gth_atompot%cval_nlcc(j, i)
398 0 : corden(:) = corden(:) + fe(:)*rc(:)**(2*j - 2)*cval
399 : END DO
400 : END DO
401 0 : fun(1:nmesh) = corden(1:nmesh)*rp(1:nmesh)
402 : CALL sirius_add_atom_type_radial_function(sctx, label, "ps_rho_core", &
403 0 : fun(1), nmesh)
404 0 : DEALLOCATE (corden, fe, rc)
405 : END IF
406 :
407 : ! local potential
408 6 : ALLOCATE (locpot(nmesh))
409 1860 : locpot(:) = 0.0_dp
410 2 : CALL atom_local_potential(locpot, gth_atompot, rp)
411 1860 : fun(1:nmesh) = locpot(1:nmesh)
412 : CALL sirius_add_atom_type_radial_function(sctx, label, "vloc", &
413 2 : fun(1), nmesh)
414 2 : DEALLOCATE (locpot)
415 : !
416 2 : NULLIFY (density, wavefunction, wfninfo)
417 : CALL calculate_atomic_orbitals(atomic_kind_set(ikind), qs_kind_set(ikind), &
418 : density=density, wavefunction=wavefunction, &
419 2 : wfninfo=wfninfo, agrid=atom_grid)
420 :
421 : ! set the atomic radial functions
422 6 : DO iwf = 1, SIZE(wavefunction, 2)
423 4 : focc = wfninfo(1, iwf)
424 4 : l = NINT(wfninfo(2, iwf))
425 : ! we can not easily get the principal quantum number
426 4 : nu = -1
427 4 : IF (up) THEN
428 0 : fun(1:nmesh) = wavefunction(1:nmesh, iwf)*rp(i)
429 : ELSE
430 3720 : DO i = 1, nmesh
431 3720 : fun(i) = wavefunction(nmesh - i + 1, iwf)*rp(i)
432 : END DO
433 : END IF
434 : CALL sirius_add_atom_type_radial_function(sctx, &
435 : label, "ps_atomic_wf", &
436 6 : fun(1), nmesh, l=l, occ=REAL(focc, KIND=C_DOUBLE), n=nu)
437 : END DO
438 :
439 : ! set total charge density of a free atom (to compute initial rho(r))
440 2 : IF (up) THEN
441 0 : fun(1:nmesh) = fourpi*density(1:nmesh)*atom_grid%rad(1:nmesh)**2
442 : ELSE
443 1860 : DO i = 1, nmesh
444 1860 : fun(i) = fourpi*density(nmesh - i + 1)*atom_grid%rad(nmesh - i + 1)**2
445 : END DO
446 : END IF
447 : CALL sirius_add_atom_type_radial_function(sctx, label, "ps_rho_total", &
448 2 : fun(1), nmesh)
449 :
450 2 : IF (ASSOCIATED(density)) DEALLOCATE (density)
451 2 : IF (ASSOCIATED(wavefunction)) DEALLOCATE (wavefunction)
452 2 : IF (ASSOCIATED(wfninfo)) DEALLOCATE (wfninfo)
453 :
454 2 : CALL deallocate_grid_atom(atom_grid)
455 2 : DEALLOCATE (rp, fun)
456 2 : DEALLOCATE (gth_atompot)
457 : !
458 : ELSE
459 : CALL cp_abort(__LOCATION__, &
460 0 : "CP2K/SIRIUS: atomic kind needs UPF or GTH potential definition")
461 : END IF
462 :
463 : CALL get_qs_kind(qs_kind_set(ikind), &
464 : dft_plus_u_atom=dft_plus_u_atom, &
465 : l_of_dft_plus_u=lu, &
466 : n_of_dft_plus_u=nu, &
467 : u_minus_j_target=u_minus_j, &
468 : U_of_dft_plus_u=U_u, &
469 : J_of_dft_plus_u=J_u, &
470 : alpha_of_dft_plus_u=alpha_u, &
471 : beta_of_dft_plus_u=beta_u, &
472 : J0_of_dft_plus_u=J0_u, &
473 22 : occupation_of_dft_plus_u=occ_u)
474 :
475 60 : IF (dft_plus_u_atom) THEN
476 0 : IF (nu < 1) THEN
477 0 : CPABORT("CP2K/SIRIUS (hubbard): principal quantum number not specified")
478 : END IF
479 :
480 0 : IF (lu < 0) THEN
481 0 : CPABORT("CP2K/SIRIUS (hubbard): l can not be negative.")
482 : END IF
483 :
484 0 : IF (occ_u < 0.0) THEN
485 0 : CPABORT("CP2K/SIRIUS (hubbard): the occupation number can not be negative.")
486 : END IF
487 :
488 0 : j_t(:) = 0.0
489 0 : IF (ABS(u_minus_j) < 1e-8) THEN
490 0 : j_t(1) = J_u
491 : CALL sirius_set_atom_type_hubbard(sctx, label, lu, nu, &
492 0 : occ_u, U_u, j_t, alpha_u, beta_u, J0_u)
493 : ELSE
494 : CALL sirius_set_atom_type_hubbard(sctx, label, lu, nu, &
495 0 : occ_u, u_minus_j, j_t, alpha_u, beta_u, J0_u)
496 : END IF
497 : END IF
498 :
499 : END DO
500 :
501 : ! add atoms to the unit cell
502 : ! WARNING: sirius accepts only fractional coordinates;
503 16 : natom = SIZE(particle_set)
504 50 : DO iatom = 1, natom
505 136 : vr(1:3) = particle_set(iatom)%r(1:3)
506 34 : CALL real_to_scaled(vs, vr, my_cell)
507 34 : atomic_kind => particle_set(iatom)%atomic_kind
508 34 : ikind = atomic_kind%kind_number
509 34 : CALL get_atomic_kind(atomic_kind, name=label)
510 34 : CALL get_qs_kind(qs_kind_set(ikind), zeff=zeff, magnetization=magnetization)
511 : ! angle of magnetization might come from input Atom x y z mx my mz
512 : ! or as an angle?
513 : ! Answer : SIRIUS only accept the magnetization as mx, my, mz
514 34 : IF (num_mag_dims .EQ. 3) THEN
515 2 : angle1 = 0.0_dp
516 2 : angle2 = 0.0_dp
517 2 : v1(1) = magnetization*SIN(angle1)*COS(angle2)
518 2 : v1(2) = magnetization*SIN(angle1)*SIN(angle2)
519 2 : v1(3) = magnetization*COS(angle1)
520 : ELSE
521 32 : v1 = 0._dp
522 32 : v1(3) = magnetization
523 : END IF
524 34 : v2(1:3) = vs(1:3)
525 84 : CALL sirius_add_atom(sctx, label, v2(1), v1(1))
526 : END DO
527 :
528 16 : CALL sirius_set_mpi_grid_dims(sctx, 2, mpi_grid_dims)
529 :
530 : ! initialize global variables/indices/arrays/etc. of the simulation
531 16 : CALL sirius_initialize_context(sctx)
532 :
533 : ! strictly speaking the parameter use_symmetry is initialized at the
534 : ! beginning but it does no harm to do it that way
535 16 : IF (use_symmetry) THEN
536 14 : CALL sirius_create_kset_from_grid(sctx, k_grid(1), k_shift(1), use_symmetry=.TRUE., kset_handler=ks_handler)
537 : ELSE
538 2 : CALL sirius_create_kset_from_grid(sctx, k_grid(1), k_shift(1), use_symmetry=.FALSE., kset_handler=ks_handler)
539 : END IF
540 : ! create ground-state class
541 16 : CALL sirius_create_ground_state(ks_handler, gs_handler)
542 :
543 16 : CALL pwdft_env_set(pwdft_env, sctx=sctx, gs_handler=gs_handler, ks_handler=ks_handler)
544 : #endif
545 32 : END SUBROUTINE cp_sirius_create_env
546 :
547 : !***************************************************************************************************
548 : !> \brief ...
549 : !> \param pwdft_env ...
550 : !> \param
551 : !> \par History
552 : !> 07.2018 Update the Sirius environment
553 : !> \author JHU
554 : ! **************************************************************************************************
555 16 : SUBROUTINE cp_sirius_update_context(pwdft_env)
556 : TYPE(pwdft_environment_type), POINTER :: pwdft_env
557 :
558 : INTEGER :: iatom, natom
559 : REAL(KIND=C_DOUBLE), DIMENSION(3) :: a1, a2, a3, v2
560 : REAL(KIND=dp), DIMENSION(3) :: vr, vs
561 : TYPE(cell_type), POINTER :: my_cell
562 16 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
563 : TYPE(qs_subsys_type), POINTER :: qs_subsys
564 : TYPE(sirius_context_handler) :: sctx
565 : TYPE(sirius_ground_state_handler) :: gs_handler
566 :
567 0 : CPASSERT(ASSOCIATED(pwdft_env))
568 16 : CALL pwdft_env_get(pwdft_env, sctx=sctx, gs_handler=gs_handler)
569 :
570 : ! get current positions and lattice vectors
571 16 : CALL pwdft_env_get(pwdft_env=pwdft_env, qs_subsys=qs_subsys)
572 :
573 : ! lattice vectors of the unit cell should be in [a.u.] (length is in [a.u.])
574 16 : CALL qs_subsys_get(qs_subsys, cell=my_cell)
575 64 : a1(:) = my_cell%hmat(:, 1)
576 64 : a2(:) = my_cell%hmat(:, 2)
577 64 : a3(:) = my_cell%hmat(:, 3)
578 16 : CALL sirius_set_lattice_vectors(sctx, a1(1), a2(1), a3(1))
579 :
580 : ! new atomic positions
581 16 : CALL qs_subsys_get(qs_subsys, particle_set=particle_set)
582 16 : natom = SIZE(particle_set)
583 50 : DO iatom = 1, natom
584 136 : vr(1:3) = particle_set(iatom)%r(1:3)
585 34 : CALL real_to_scaled(vs, vr, my_cell)
586 34 : v2(1:3) = vs(1:3)
587 50 : CALL sirius_set_atom_position(sctx, iatom, v2(1))
588 : END DO
589 :
590 : ! update ground-state class
591 16 : CALL sirius_update_ground_state(gs_handler)
592 :
593 16 : CALL pwdft_env_set(pwdft_env, sctx=sctx, gs_handler=gs_handler)
594 :
595 16 : END SUBROUTINE cp_sirius_update_context
596 :
597 : ! **************************************************************************************************
598 : !> \brief ...
599 : !> \param sctx ...
600 : !> \param section ...
601 : !> \param section_name ...
602 : ! **************************************************************************************************
603 80 : SUBROUTINE cp_sirius_fill_in_section(sctx, section, section_name)
604 : TYPE(sirius_context_handler), INTENT(INOUT) :: sctx
605 : TYPE(section_vals_type), POINTER :: section
606 : CHARACTER(*), INTENT(in) :: section_name
607 :
608 : CHARACTER(len=256), TARGET :: option_name
609 : CHARACTER(len=4096) :: description, usage
610 80 : CHARACTER(len=80), DIMENSION(:), POINTER :: tmp
611 : CHARACTER(len=80), TARGET :: str
612 : INTEGER :: ctype, elem, ic, j
613 80 : INTEGER, DIMENSION(:), POINTER :: ivals
614 : INTEGER, TARGET :: enum_length, ival, length, &
615 : num_possible_values, number_of_options
616 : LOGICAL :: explicit
617 80 : LOGICAL, DIMENSION(:), POINTER :: lvals
618 : LOGICAL, TARGET :: found, lval
619 80 : REAL(kind=dp), DIMENSION(:), POINTER :: rvals
620 : REAL(kind=dp), TARGET :: rval
621 :
622 80 : NULLIFY (rvals)
623 80 : NULLIFY (ivals)
624 80 : CALL sirius_option_get_section_length(section_name, number_of_options)
625 :
626 1712 : DO elem = 1, number_of_options
627 1632 : option_name = ''
628 : CALL sirius_option_get_info(section_name, &
629 : elem, &
630 : option_name, &
631 : 256, &
632 : ctype, &
633 : num_possible_values, &
634 : enum_length, &
635 : description, &
636 : 4096, &
637 : usage, &
638 1632 : 4096)
639 1712 : IF ((option_name /= 'memory_usage') .AND. (option_name /= 'xc_functionals') .AND. (option_name /= 'vk')) THEN
640 1600 : CALL section_vals_val_get(section, option_name, explicit=found)
641 1600 : IF (found) THEN
642 128 : SELECT CASE (ctype)
643 : CASE (SIRIUS_INTEGER_TYPE)
644 128 : CALL section_vals_val_get(section, option_name, i_val=ival)
645 128 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(ival))
646 : CASE (SIRIUS_NUMBER_TYPE)
647 104 : CALL section_vals_val_get(section, option_name, r_val=rval)
648 104 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(rval))
649 : CASE (SIRIUS_LOGICAL_TYPE)
650 24 : CALL section_vals_val_get(section, option_name, l_val=lval)
651 24 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(lval))
652 : CASE (SIRIUS_STRING_TYPE) ! string nightmare
653 100 : str = ''
654 100 : CALL section_vals_val_get(section, option_name, explicit=explicit, c_val=str)
655 100 : str = TRIM(ADJUSTL(str))
656 8100 : DO j = 1, LEN(str)
657 8000 : ic = ICHAR(str(j:j))
658 8100 : IF (ic >= 65 .AND. ic < 90) str(j:j) = CHAR(ic + 32)
659 : END DO
660 :
661 100 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(str), max_length=LEN_TRIM(str))
662 : CASE (SIRIUS_INTEGER_ARRAY_TYPE)
663 16 : CALL section_vals_val_get(section, option_name, i_vals=ivals)
664 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(ivals(1)), &
665 16 : max_length=num_possible_values)
666 : CASE (SIRIUS_NUMBER_ARRAY_TYPE)
667 0 : CALL section_vals_val_get(section, option_name, r_vals=rvals)
668 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(rvals(1)), &
669 0 : max_length=num_possible_values)
670 : CASE (SIRIUS_LOGICAL_ARRAY_TYPE)
671 0 : CALL section_vals_val_get(section, option_name, l_vals=lvals)
672 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(lvals(1)), &
673 0 : max_length=num_possible_values)
674 : CASE (SIRIUS_STRING_ARRAY_TYPE)
675 0 : CALL section_vals_val_get(section, option_name, explicit=explicit, n_rep_val=length)
676 372 : DO j = 1, length
677 0 : str = ''
678 0 : CALL section_vals_val_get(section, option_name, i_rep_val=j, explicit=explicit, c_vals=tmp)
679 0 : str = TRIM(ADJUSTL(tmp(j)))
680 : CALL sirius_option_set(sctx, section_name, option_name, ctype, C_LOC(str), &
681 0 : max_length=LEN_TRIM(str), append=.TRUE.)
682 : END DO
683 : CASE DEFAULT
684 : END SELECT
685 : END IF
686 : END IF
687 : END DO
688 80 : END SUBROUTINE cp_sirius_fill_in_section
689 :
690 : !***************************************************************************************************
691 : !> \brief ...
692 : !> \param pwdft_env ...
693 : !> \param calculate_forces ...
694 : !> \param calculate_stress_tensor ...
695 : !> \param
696 : !> \par History
697 : !> 07.2018 start the Sirius library
698 : !> \author JHU
699 : ! **************************************************************************************************
700 16 : SUBROUTINE cp_sirius_energy_force(pwdft_env, calculate_forces, calculate_stress_tensor)
701 : TYPE(pwdft_environment_type), INTENT(INOUT), &
702 : POINTER :: pwdft_env
703 : LOGICAL, INTENT(IN) :: calculate_forces, calculate_stress_tensor
704 :
705 : INTEGER :: iw, n1, n2
706 : LOGICAL :: do_print, gs_converged
707 : REAL(KIND=C_DOUBLE) :: etotal
708 16 : REAL(KIND=C_DOUBLE), ALLOCATABLE, DIMENSION(:, :) :: cforces
709 : REAL(KIND=C_DOUBLE), DIMENSION(3, 3) :: cstress
710 : REAL(KIND=dp), DIMENSION(3, 3) :: stress
711 16 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: forces
712 : TYPE(cp_logger_type), POINTER :: logger
713 : TYPE(pwdft_energy_type), POINTER :: energy
714 : TYPE(section_vals_type), POINTER :: print_section, pwdft_input
715 : TYPE(sirius_ground_state_handler) :: gs_handler
716 :
717 0 : CPASSERT(ASSOCIATED(pwdft_env))
718 :
719 16 : NULLIFY (logger)
720 16 : logger => cp_get_default_logger()
721 16 : iw = cp_logger_get_default_io_unit(logger)
722 :
723 16 : CALL pwdft_env_get(pwdft_env=pwdft_env, gs_handler=gs_handler)
724 16 : CALL sirius_find_ground_state(gs_handler, converged=gs_converged)
725 :
726 16 : IF (gs_converged) THEN
727 16 : IF (iw > 0) WRITE (iw, '(A)') "CP2K/SIRIUS: ground state is converged"
728 : ELSE
729 0 : IF (pwdft_env%ignore_convergence_failure) THEN
730 0 : IF (iw > 0) WRITE (iw, '(A)') "CP2K/SIRIUS Warning: ground state is not converged"
731 : ELSE
732 0 : CPABORT("CP2K/SIRIUS (ground state): SIRIUS did not converge.")
733 : END IF
734 : END IF
735 16 : IF (iw > 0) CALL m_flush(iw)
736 :
737 16 : CALL pwdft_env_get(pwdft_env=pwdft_env, energy=energy)
738 : etotal = 0.0_C_DOUBLE
739 :
740 16 : CALL sirius_get_energy(gs_handler, 'band-gap', etotal)
741 16 : energy%band_gap = etotal
742 :
743 : etotal = 0.0_C_DOUBLE
744 16 : CALL sirius_get_energy(gs_handler, 'total', etotal)
745 16 : energy%etotal = etotal
746 :
747 : ! extract entropy (TS returned by sirius is always negative, sign
748 : ! convention in QE)
749 : etotal = 0.0_C_DOUBLE
750 16 : CALL sirius_get_energy(gs_handler, 'demet', etotal)
751 16 : energy%entropy = -etotal
752 :
753 16 : IF (calculate_forces) THEN
754 12 : CALL pwdft_env_get(pwdft_env=pwdft_env, forces=forces)
755 12 : n1 = SIZE(forces, 1)
756 12 : n2 = SIZE(forces, 2)
757 :
758 48 : ALLOCATE (cforces(n2, n1))
759 132 : cforces = 0.0_C_DOUBLE
760 12 : CALL sirius_get_forces(gs_handler, 'total', cforces)
761 : ! Sirius computes the forces but cp2k use the gradient everywhere
762 : ! so a minus sign is needed.
763 : ! note also that sirius and cp2k store the forces transpose to each other
764 : ! sirius : forces(coordinates, atoms)
765 : ! cp2k : forces(atoms, coordinates)
766 138 : forces = -TRANSPOSE(cforces(:, :))
767 12 : DEALLOCATE (cforces)
768 : END IF
769 :
770 16 : IF (calculate_stress_tensor) THEN
771 0 : cstress = 0.0_C_DOUBLE
772 0 : CALL sirius_get_stress_tensor(gs_handler, 'total', cstress)
773 0 : stress(1:3, 1:3) = cstress(1:3, 1:3)
774 0 : CALL pwdft_env_set(pwdft_env=pwdft_env, stress=stress)
775 : END IF
776 :
777 16 : CALL pwdft_env_get(pwdft_env=pwdft_env, pwdft_input=pwdft_input)
778 16 : print_section => section_vals_get_subs_vals(pwdft_input, "PRINT")
779 16 : CALL section_vals_get(print_section, explicit=do_print)
780 16 : IF (do_print) THEN
781 2 : CALL cp_sirius_print_results(pwdft_env, print_section)
782 : END IF
783 32 : END SUBROUTINE cp_sirius_energy_force
784 :
785 : !***************************************************************************************************
786 : !> \brief ...
787 : !> \param pwdft_env ...
788 : !> \param print_section ...
789 : !> \param
790 : !> \par History
791 : !> 12.2019 init
792 : !> \author JHU
793 : ! **************************************************************************************************
794 2 : SUBROUTINE cp_sirius_print_results(pwdft_env, print_section)
795 : TYPE(pwdft_environment_type), INTENT(INOUT), &
796 : POINTER :: pwdft_env
797 : TYPE(section_vals_type), POINTER :: print_section
798 :
799 : CHARACTER(LEN=default_string_length) :: my_act, my_pos
800 : INTEGER :: i, ik, iounit, ispn, iterstep, iv, iw, &
801 : nbands, nhist, nkpts, nspins
802 : INTEGER(KIND=C_INT) :: cint
803 : LOGICAL :: append, dos, ionode
804 : REAL(KIND=C_DOUBLE) :: creal
805 2 : REAL(KIND=C_DOUBLE), ALLOCATABLE, DIMENSION(:) :: slist
806 : REAL(KIND=dp) :: de, e_fermi(2), emax, emin, eval
807 2 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: wkpt
808 2 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: ehist, hist, occval
809 2 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: energies, occupations
810 : TYPE(cp_logger_type), POINTER :: logger
811 : TYPE(sirius_context_handler) :: sctx
812 : TYPE(sirius_ground_state_handler) :: gs_handler
813 : TYPE(sirius_kpoint_set_handler) :: ks_handler
814 :
815 2 : NULLIFY (logger)
816 4 : logger => cp_get_default_logger()
817 : ionode = logger%para_env%is_source()
818 2 : iounit = cp_logger_get_default_io_unit(logger)
819 :
820 : ! Density of States
821 2 : dos = BTEST(cp_print_key_should_output(logger%iter_info, print_section, "DOS"), cp_p_file)
822 2 : IF (dos) THEN
823 2 : CALL pwdft_env_get(pwdft_env, ks_handler=ks_handler)
824 2 : CALL pwdft_env_get(pwdft_env, gs_handler=gs_handler)
825 2 : CALL pwdft_env_get(pwdft_env, sctx=sctx)
826 :
827 2 : CALL section_vals_val_get(print_section, "DOS%DELTA_E", r_val=de)
828 2 : CALL section_vals_val_get(print_section, "DOS%APPEND", l_val=append)
829 :
830 2 : CALL sirius_get_num_kpoints(ks_handler, cint)
831 2 : nkpts = cint
832 2 : CALL sirius_get_parameters(sctx, num_bands=cint)
833 2 : nbands = cint
834 2 : CALL sirius_get_parameters(sctx, num_spins=cint)
835 2 : nspins = cint
836 2 : e_fermi(:) = 0.0_dp
837 10 : ALLOCATE (energies(nbands, nspins, nkpts))
838 442 : energies = 0.0_dp
839 8 : ALLOCATE (occupations(nbands, nspins, nkpts))
840 442 : occupations = 0.0_dp
841 6 : ALLOCATE (wkpt(nkpts))
842 6 : ALLOCATE (slist(nbands))
843 10 : DO ik = 1, nkpts
844 8 : CALL sirius_get_kpoint_properties(ks_handler, ik, creal)
845 10 : wkpt(ik) = creal
846 : END DO
847 10 : DO ik = 1, nkpts
848 26 : DO ispn = 1, nspins
849 16 : CALL sirius_get_band_energies(ks_handler, ik, ispn, slist)
850 432 : energies(1:nbands, ispn, ik) = slist(1:nbands)
851 16 : CALL sirius_get_band_occupancies(ks_handler, ik, ispn, slist)
852 440 : occupations(1:nbands, ispn, ik) = slist(1:nbands)
853 : END DO
854 : END DO
855 442 : emin = MINVAL(energies)
856 442 : emax = MAXVAL(energies)
857 2 : nhist = NINT((emax - emin)/de) + 1
858 16 : ALLOCATE (hist(nhist, nspins), occval(nhist, nspins), ehist(nhist, nspins))
859 16110 : hist = 0.0_dp
860 16110 : occval = 0.0_dp
861 16110 : ehist = 0.0_dp
862 :
863 10 : DO ik = 1, nkpts
864 26 : DO ispn = 1, nspins
865 440 : DO i = 1, nbands
866 416 : eval = energies(i, ispn, ik) - emin
867 416 : iv = NINT(eval/de) + 1
868 416 : CPASSERT((iv > 0) .AND. (iv <= nhist))
869 416 : hist(iv, ispn) = hist(iv, ispn) + wkpt(ik)
870 432 : occval(iv, ispn) = occval(iv, ispn) + wkpt(ik)*occupations(i, ispn, ik)
871 : END DO
872 : END DO
873 : END DO
874 16110 : hist = hist/REAL(nbands, KIND=dp)
875 8054 : DO i = 1, nhist
876 24158 : ehist(i, 1:nspins) = emin + (i - 1)*de
877 : END DO
878 :
879 2 : iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
880 2 : my_act = "WRITE"
881 2 : IF (append .AND. iterstep > 1) THEN
882 0 : my_pos = "APPEND"
883 : ELSE
884 2 : my_pos = "REWIND"
885 : END IF
886 :
887 : iw = cp_print_key_unit_nr(logger, print_section, "DOS", &
888 : extension=".dos", file_position=my_pos, file_action=my_act, &
889 2 : file_form="FORMATTED")
890 2 : IF (iw > 0) THEN
891 1 : IF (nspins == 2) THEN
892 : WRITE (UNIT=iw, FMT="(T2,A,I0,A,2F12.6)") &
893 1 : "# DOS at iteration step i = ", iterstep, ", E_Fermi[a.u.] = ", e_fermi(1:2)
894 1 : WRITE (UNIT=iw, FMT="(T2,A, A)") " Energy[a.u.] Alpha_Density Occupation", &
895 2 : " Beta_Density Occupation"
896 : ELSE
897 : WRITE (UNIT=iw, FMT="(T2,A,I0,A,F12.6)") &
898 0 : "# DOS at iteration step i = ", iterstep, ", E_Fermi[a.u.] = ", e_fermi(1)
899 0 : WRITE (UNIT=iw, FMT="(T2,A)") " Energy[a.u.] Density Occupation"
900 : END IF
901 4027 : DO i = 1, nhist
902 4026 : eval = emin + (i - 1)*de
903 4027 : IF (nspins == 2) THEN
904 4026 : WRITE (UNIT=iw, FMT="(F15.8,4F15.4)") eval, hist(i, 1), occval(i, 1), &
905 8052 : hist(i, 2), occval(i, 2)
906 : ELSE
907 0 : WRITE (UNIT=iw, FMT="(F15.8,2F15.4)") eval, hist(i, 1), occval(i, 1)
908 : END IF
909 : END DO
910 : END IF
911 2 : CALL cp_print_key_finished_output(iw, logger, print_section, "DOS")
912 :
913 2 : DEALLOCATE (energies, occupations, wkpt, slist)
914 4 : DEALLOCATE (hist, occval, ehist)
915 : END IF
916 4 : END SUBROUTINE cp_sirius_print_results
917 :
918 : END MODULE sirius_interface
919 :
920 : #else
921 :
922 : !***************************************************************************************************
923 : !> \brief Empty implementation in case SIRIUS is not compiled in.
924 : !***************************************************************************************************
925 : MODULE sirius_interface
926 : USE pwdft_environment_types, ONLY: pwdft_environment_type
927 : #include "./base/base_uses.f90"
928 :
929 : IMPLICIT NONE
930 : PRIVATE
931 :
932 : PUBLIC :: cp_sirius_init, cp_sirius_finalize
933 : PUBLIC :: cp_sirius_create_env, cp_sirius_energy_force, cp_sirius_update_context
934 :
935 : CONTAINS
936 :
937 : ! **************************************************************************************************
938 : !> \brief Empty implementation in case SIRIUS is not compiled in.
939 : ! **************************************************************************************************
940 : SUBROUTINE cp_sirius_init()
941 : END SUBROUTINE cp_sirius_init
942 :
943 : ! **************************************************************************************************
944 : !> \brief Empty implementation in case SIRIUS is not compiled in.
945 : ! **************************************************************************************************
946 : SUBROUTINE cp_sirius_finalize()
947 : END SUBROUTINE cp_sirius_finalize
948 :
949 : ! **************************************************************************************************
950 : !> \brief Empty implementation in case SIRIUS is not compiled in.
951 : !> \param pwdft_env ...
952 : ! **************************************************************************************************
953 : SUBROUTINE cp_sirius_create_env(pwdft_env)
954 : TYPE(pwdft_environment_type), POINTER :: pwdft_env
955 :
956 : MARK_USED(pwdft_env)
957 : CPABORT("Sirius library is missing")
958 : END SUBROUTINE cp_sirius_create_env
959 :
960 : ! **************************************************************************************************
961 : !> \brief Empty implementation in case SIRIUS is not compiled in.
962 : !> \param pwdft_env ...
963 : !> \param calculate_forces ...
964 : !> \param calculate_stress ...
965 : ! **************************************************************************************************
966 : SUBROUTINE cp_sirius_energy_force(pwdft_env, calculate_forces, calculate_stress)
967 : TYPE(pwdft_environment_type), POINTER :: pwdft_env
968 : LOGICAL :: calculate_forces, calculate_stress
969 :
970 : MARK_USED(pwdft_env)
971 : MARK_USED(calculate_forces)
972 : MARK_USED(calculate_stress)
973 : CPABORT("Sirius library is missing")
974 : END SUBROUTINE cp_sirius_energy_force
975 :
976 : ! **************************************************************************************************
977 : !> \brief Empty implementation in case SIRIUS is not compiled in.
978 : !> \param pwdft_env ...
979 : ! **************************************************************************************************
980 : SUBROUTINE cp_sirius_update_context(pwdft_env)
981 : TYPE(pwdft_environment_type), POINTER :: pwdft_env
982 :
983 : MARK_USED(pwdft_env)
984 : CPABORT("Sirius library is missing")
985 : END SUBROUTINE cp_sirius_update_context
986 :
987 : END MODULE sirius_interface
988 :
989 : #endif
|