(git:15c1bfc)
Loading...
Searching...
No Matches
sirius_interface.F
Go to the documentation of this file.
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)
16 USE iso_c_binding, ONLY: c_double,&
17 c_int,&
18 c_loc
22 USE atom_upf, ONLY: atom_upfpot_type
26 USE cell_types, ONLY: cell_type,&
31 USE cp_output_handling, ONLY: cp_p_file,&
45 USE kinds, ONLY: default_string_length,&
46 dp
47 USE machine, ONLY: m_flush
48 USE mathconstants, ONLY: fourpi,&
49 gamma1
52 USE physcon, ONLY: massunit
61 USE qs_kind_types, ONLY: get_qs_kind,&
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
94
95CONTAINS
96
97!***************************************************************************************************
98!> \brief ...
99!> \param
100!> \par History
101!> 07.2018 start the Sirius library
102!> \author JHU
103! **************************************************************************************************
104 SUBROUTINE cp_sirius_init()
105 CALL sirius_initialize(.false.)
106 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 SUBROUTINE cp_sirius_finalize()
116 CALL sirius_finalize(.false., .false., .false.)
117 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 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, ii, jj, iatom, ibeta, ifun, ikind, iwf, j, l, &
134 n, ns, natom, nbeta, nbs, nkind, nmesh, &
135 num_mag_dims, sirius_mpi_comm, vdw_func, nu, lu, output_unit
136 INTEGER, DIMENSION(:), POINTER :: mpi_grid_dims
137 INTEGER(KIND=C_INT), DIMENSION(3) :: k_grid, k_shift
138 INTEGER, DIMENSION(:), POINTER :: kk
139 LOGICAL :: up, use_ref_cell
140 LOGICAL(4) :: use_so, use_symmetry, dft_plus_u_atom
141 REAL(KIND=c_double), ALLOCATABLE, DIMENSION(:) :: fun
142 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, vnlp, vnlm
147 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: beta, corden, ef, fe, locpot, rc, rp
148 REAL(KIND=dp), DIMENSION(3) :: vr, vs, j_t
149 REAL(KIND=dp), DIMENSION(:), POINTER :: density
150 REAL(KIND=dp), DIMENSION(:, :), POINTER :: wavefunction, wfninfo
151 TYPE(atom_gthpot_type), POINTER :: gth_atompot
152 TYPE(atom_upfpot_type), POINTER :: upf_pot
153 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 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
160 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 cpassert(ASSOCIATED(pwdft_env))
169
170 output_unit = cp_logger_get_default_io_unit()
171 ! create context of simulation
172 CALL pwdft_env_get(pwdft_env, para_env=para_env)
173 sirius_mpi_comm = para_env%get_handle()
174 CALL sirius_create_context(sirius_mpi_comm, sctx)
175
176! the "fun" starts.
177
178 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 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 IF (ASSOCIATED(xc_section)) THEN
186 ifun = 0
187 DO
188 ifun = ifun + 1
189 xc_fun => section_vals_get_subs_vals2(xc_section, i_section=ifun)
190 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 CALL sirius_add_xc_functional(sctx, "XC_"//trim(xc_fun%section%name))
194 END DO
195 END IF
196
197 ! import control section
198 pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "control")
199 IF (ASSOCIATED(pwdft_sub_section)) THEN
200 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "control")
201 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 pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "parameters")
206
207 IF (ASSOCIATED(pwdft_sub_section)) THEN
208 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "parameters")
209 CALL section_vals_val_get(pwdft_sub_section, "ngridk", i_vals=kk)
210 k_grid(1) = kk(1)
211 k_grid(2) = kk(2)
212 k_grid(3) = kk(3)
213
214 CALL section_vals_val_get(pwdft_sub_section, "shiftk", i_vals=kk)
215 k_shift(1) = kk(1)
216 k_shift(2) = kk(2)
217 k_shift(3) = kk(3)
218 CALL section_vals_val_get(pwdft_sub_section, "num_mag_dims", i_val=num_mag_dims)
219 CALL section_vals_val_get(pwdft_sub_section, "use_symmetry", l_val=use_symmetry)
220 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 CALL section_vals_val_get(pwdft_sub_section, "vdw_functional", i_val=vdw_func)
226 SELECT CASE (vdw_func)
227 CASE (sirius_func_vdwdf)
228 CALL sirius_add_xc_functional(sctx, "XC_FUNC_VDWDF")
229 CASE (sirius_func_vdwdf2)
230 CALL sirius_add_xc_functional(sctx, "XC_FUNC_VDWDF2")
232 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 pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "mixer")
241 IF (ASSOCIATED(pwdft_sub_section)) THEN
242 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "mixer")
243 END IF
244
245! import settings section
246 pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "settings")
247 IF (ASSOCIATED(pwdft_sub_section)) THEN
248 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section, "settings")
249 END IF
250
251 ! import solver section
252 pwdft_sub_section => section_vals_get_subs_vals(pwdft_section, "iterative_solver")
253 IF (ASSOCIATED(pwdft_sub_section)) THEN
254 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 CALL sirius_import_parameters(sctx, '{}')
290
291! lattice vectors of the unit cell should be in [a.u.] (length is in [a.u.])
292 CALL pwdft_env_get(pwdft_env=pwdft_env, qs_subsys=qs_subsys)
293 CALL qs_subsys_get(qs_subsys, cell=my_cell, use_ref_cell=use_ref_cell)
294 a1(:) = my_cell%hmat(:, 1)
295 a2(:) = my_cell%hmat(:, 2)
296 a3(:) = my_cell%hmat(:, 3)
297 CALL sirius_set_lattice_vectors(sctx, a1(1), a2(1), a3(1))
298
299 IF (use_ref_cell) THEN
300 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 particle_set=particle_set)
308 nkind = SIZE(atomic_kind_set)
309 DO ikind = 1, nkind
310 CALL get_atomic_kind(atomic_kind_set(ikind), &
311 name=label, element_symbol=element_symbol, mass=mass)
312 CALL get_qs_kind(qs_kind_set(ikind), zeff=zeff)
313 NULLIFY (upf_pot, gth_potential)
314 CALL get_qs_kind(qs_kind_set(ikind), upf_potential=upf_pot, gth_potential=gth_potential)
315
316 IF (ASSOCIATED(upf_pot)) THEN
317 CALL sirius_add_atom_type(sctx, label, fname=upf_pot%filename, &
318 symbol=element_symbol, &
319 mass=real(mass/massunit, kind=c_double))
320
321 ELSEIF (ASSOCIATED(gth_potential)) THEN
322!
323 NULLIFY (atom_grid)
324 CALL allocate_grid_atom(atom_grid)
325 nmesh = 929
326 atom_grid%nr = nmesh
327 CALL create_grid_atom(atom_grid, nmesh, 1, 1, 0, do_gapw_log)
328 ALLOCATE (rp(nmesh), fun(nmesh))
329 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 rp(1:nmesh) = atom_grid%rad(1:nmesh)
336 ELSE
337 DO i = 1, nmesh
338 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 spin_orbit=gth_potential%soc)
347!
348 ALLOCATE (gth_atompot)
349 CALL gth_potential_conversion(gth_potential, gth_atompot)
350! set radial grid
351 fun(1:nmesh) = rp(1:nmesh)
352 CALL sirius_set_atom_type_radial_grid(sctx, label, nmesh, fun(1))
353! set beta-projectors
354! GTH SOC uses the same projectors, SIRIUS can use the same or different projectors for l+1/2, l-1/2 (l > 0 l+1/2 l < 0 l-/2 )
355 ALLOCATE (ef(nmesh), beta(nmesh))
356 ibeta = 0
357 DO l = 0, 3
358 IF (gth_atompot%nl(l) == 0) cycle
359 rl = gth_atompot%rcnl(l)
360! we need to multiply by r so that data transferred to sirius are r \beta(r) not beta(r)
361 ef(1:nmesh) = exp(-0.5_dp*rp(1:nmesh)*rp(1:nmesh)/(rl*rl))
362 DO i = 1, gth_atompot%nl(l)
363 pf = rl**(l + 0.5_dp*(4._dp*i - 1._dp))
364 j = l + 2*i - 1
365 pf = sqrt(2._dp)/(pf*sqrt(gamma1(j)))
366 beta(:) = pf*rp**(l + 2*i - 2)*ef
367 ibeta = ibeta + 1
368 fun(1:nmesh) = beta(1:nmesh)*rp(1:nmesh)
369 CALL sirius_add_atom_type_radial_function(sctx, label, &
370 "beta", fun(1), nmesh, l=l)
371 ! we double the number of beta projectors for SO and l>0
372 IF (gth_atompot%soc .AND. l /= 0) THEN
373 CALL sirius_add_atom_type_radial_function(sctx, label, &
374 "beta", fun(1), nmesh, l=-l)
375 END IF
376 END DO
377 END DO
378 DEALLOCATE (ef, beta)
379 nbeta = ibeta
380
381! nonlocal PP matrix elements
382 IF (gth_atompot%soc) THEN
383 nbs = 2*nbeta - gth_atompot%nl(0)
384 ALLOCATE (dion(nbs, nbs))
385 ELSE
386 ALLOCATE (dion(nbeta, nbeta))
387 END IF
388 dion = 0.0_dp
389 IF (gth_atompot%soc) THEN
390 ns = gth_atompot%nl(0)
391 IF (ns /= 0) THEN
392 dion(1:ns, 1:ns) = gth_atompot%hnl(1:ns, 1:ns, 0)
393 END IF
394 DO l = 1, 3
395 IF (gth_atompot%nl(l) == 0) cycle
396 DO i = 1, gth_atompot%nl(l)
397 ii = ns + 2*sum(gth_atompot%nl(1:l - 1))
398 ii = ii + 2*(i - 1) + 1
399 DO j = 1, gth_atompot%nl(l)
400 jj = ns + 2*sum(gth_atompot%nl(1:l - 1))
401 jj = jj + 2*(j - 1) + 1
402 vnlp = gth_atompot%hnl(i, j, l) + 0.5_dp*l*gth_atompot%knl(i, j, l)
403 vnlm = gth_atompot%hnl(i, j, l) - 0.5_dp*(l + 1)*gth_atompot%knl(i, j, l)
404 dion(ii, jj) = vnlp
405 dion(ii + 1, jj + 1) = vnlm
406 END DO
407 END DO
408 END DO
409 CALL sirius_set_atom_type_dion(sctx, label, nbs, dion(1, 1))
410 ELSE
411 DO l = 0, 3
412 IF (gth_atompot%nl(l) == 0) cycle
413 ibeta = sum(gth_atompot%nl(0:l - 1)) + 1
414 i = ibeta + gth_atompot%nl(l) - 1
415 dion(ibeta:i, ibeta:i) = gth_atompot%hnl(1:gth_atompot%nl(l), 1:gth_atompot%nl(l), l)
416 END DO
417 CALL sirius_set_atom_type_dion(sctx, label, nbeta, dion(1, 1))
418 END IF
419
420 DEALLOCATE (dion)
421
422! set non-linear core correction
423 IF (gth_atompot%nlcc) THEN
424 ALLOCATE (corden(nmesh), fe(nmesh), rc(nmesh))
425 corden(:) = 0.0_dp
426 n = gth_atompot%nexp_nlcc
427 DO i = 1, n
428 al = gth_atompot%alpha_nlcc(i)
429 rc(:) = rp(:)/al
430 fe(:) = exp(-0.5_dp*rc(:)*rc(:))
431 DO j = 1, gth_atompot%nct_nlcc(i)
432 cval = gth_atompot%cval_nlcc(j, i)
433 corden(:) = corden(:) + fe(:)*rc(:)**(2*j - 2)*cval
434 END DO
435 END DO
436 fun(1:nmesh) = corden(1:nmesh)*rp(1:nmesh)
437 CALL sirius_add_atom_type_radial_function(sctx, label, "ps_rho_core", &
438 fun(1), nmesh)
439 DEALLOCATE (corden, fe, rc)
440 END IF
441
442! local potential
443 ALLOCATE (locpot(nmesh))
444 locpot(:) = 0.0_dp
445 CALL atom_local_potential(locpot, gth_atompot, rp)
446 fun(1:nmesh) = locpot(1:nmesh)
447 CALL sirius_add_atom_type_radial_function(sctx, label, "vloc", &
448 fun(1), nmesh)
449 DEALLOCATE (locpot)
450!
451 NULLIFY (density, wavefunction, wfninfo)
452 CALL calculate_atomic_orbitals(atomic_kind_set(ikind), qs_kind_set(ikind), &
453 density=density, wavefunction=wavefunction, &
454 wfninfo=wfninfo, agrid=atom_grid)
455
456! set the atomic radial functions
457 DO iwf = 1, SIZE(wavefunction, 2)
458 focc = wfninfo(1, iwf)
459 l = nint(wfninfo(2, iwf))
460! we can not easily get the principal quantum number
461 nu = -1
462 IF (up) THEN
463 fun(1:nmesh) = wavefunction(1:nmesh, iwf)*rp(i)
464 ELSE
465 DO i = 1, nmesh
466 fun(i) = wavefunction(nmesh - i + 1, iwf)*rp(i)
467 END DO
468 END IF
469 CALL sirius_add_atom_type_radial_function(sctx, &
470 label, "ps_atomic_wf", &
471 fun(1), nmesh, l=l, occ=real(focc, kind=c_double), n=nu)
472 END DO
473
474! set total charge density of a free atom (to compute initial rho(r))
475 IF (up) THEN
476 fun(1:nmesh) = fourpi*density(1:nmesh)*atom_grid%rad(1:nmesh)**2
477 ELSE
478 DO i = 1, nmesh
479 fun(i) = fourpi*density(nmesh - i + 1)*atom_grid%rad(nmesh - i + 1)**2
480 END DO
481 END IF
482 CALL sirius_add_atom_type_radial_function(sctx, label, "ps_rho_total", &
483 fun(1), nmesh)
484
485 IF (ASSOCIATED(density)) DEALLOCATE (density)
486 IF (ASSOCIATED(wavefunction)) DEALLOCATE (wavefunction)
487 IF (ASSOCIATED(wfninfo)) DEALLOCATE (wfninfo)
488
489 CALL deallocate_grid_atom(atom_grid)
490 DEALLOCATE (rp, fun)
491 DEALLOCATE (gth_atompot)
492!
493 ELSE
494 CALL cp_abort(__location__, &
495 "CP2K/SIRIUS: atomic kind needs UPF or GTH potential definition")
496 END IF
497
498 CALL get_qs_kind(qs_kind_set(ikind), &
499 dft_plus_u_atom=dft_plus_u_atom, &
500 l_of_dft_plus_u=lu, &
501 n_of_dft_plus_u=nu, &
502 u_minus_j_target=u_minus_j, &
503 u_of_dft_plus_u=u_u, &
504 j_of_dft_plus_u=j_u, &
505 alpha_of_dft_plus_u=alpha_u, &
506 beta_of_dft_plus_u=beta_u, &
507 j0_of_dft_plus_u=j0_u, &
508 occupation_of_dft_plus_u=occ_u)
509
510 IF (dft_plus_u_atom) THEN
511 IF (nu < 1) THEN
512 cpabort("CP2K/SIRIUS (hubbard): principal quantum number not specified")
513 END IF
514
515 IF (lu < 0) THEN
516 cpabort("CP2K/SIRIUS (hubbard): l can not be negative.")
517 END IF
518
519 IF (occ_u < 0.0) THEN
520 cpabort("CP2K/SIRIUS (hubbard): the occupation number can not be negative.")
521 END IF
522
523 j_t(:) = 0.0
524 IF (abs(u_minus_j) < 1e-8) THEN
525 j_t(1) = j_u
526 CALL sirius_set_atom_type_hubbard(sctx, label, lu, nu, &
527 occ_u, u_u, j_t, alpha_u, beta_u, j0_u)
528 ELSE
529 CALL sirius_set_atom_type_hubbard(sctx, label, lu, nu, &
530 occ_u, u_minus_j, j_t, alpha_u, beta_u, j0_u)
531 END IF
532 END IF
533
534 END DO
535
536! add atoms to the unit cell
537! WARNING: sirius accepts only fractional coordinates;
538 natom = SIZE(particle_set)
539 DO iatom = 1, natom
540 vr(1:3) = particle_set(iatom)%r(1:3)
541 CALL real_to_scaled(vs, vr, my_cell)
542 atomic_kind => particle_set(iatom)%atomic_kind
543 ikind = atomic_kind%kind_number
544 CALL get_atomic_kind(atomic_kind, name=label)
545 CALL get_qs_kind(qs_kind_set(ikind), zeff=zeff, magnetization=magnetization)
546! angle of magnetization might come from input Atom x y z mx my mz
547! or as an angle?
548! Answer : SIRIUS only accept the magnetization as mx, my, mz
549 IF (num_mag_dims .EQ. 3) THEN
550 angle1 = 0.0_dp
551 angle2 = 0.0_dp
552 v1(1) = magnetization*sin(angle1)*cos(angle2)
553 v1(2) = magnetization*sin(angle1)*sin(angle2)
554 v1(3) = magnetization*cos(angle1)
555 ELSE
556 v1 = 0._dp
557 v1(3) = magnetization
558 END IF
559 v2(1:3) = vs(1:3)
560 CALL sirius_add_atom(sctx, label, v2(1), v1(1))
561 END DO
562
563 CALL sirius_set_mpi_grid_dims(sctx, 2, mpi_grid_dims)
564
565! initialize global variables/indices/arrays/etc. of the simulation
566 CALL sirius_initialize_context(sctx)
567
568 ! strictly speaking the parameter use_symmetry is initialized at the
569 ! beginning but it does no harm to do it that way
570 IF (use_symmetry) THEN
571 CALL sirius_create_kset_from_grid(sctx, k_grid(1), k_shift(1), use_symmetry=.true., kset_handler=ks_handler)
572 ELSE
573 CALL sirius_create_kset_from_grid(sctx, k_grid(1), k_shift(1), use_symmetry=.false., kset_handler=ks_handler)
574 END IF
575! create ground-state class
576 CALL sirius_create_ground_state(ks_handler, gs_handler)
577
578 CALL pwdft_env_set(pwdft_env, sctx=sctx, gs_handler=gs_handler, ks_handler=ks_handler)
579#endif
580 END SUBROUTINE cp_sirius_create_env
581
582!***************************************************************************************************
583!> \brief ...
584!> \param pwdft_env ...
585!> \param
586!> \par History
587!> 07.2018 Update the Sirius environment
588!> \author JHU
589! **************************************************************************************************
590 SUBROUTINE cp_sirius_update_context(pwdft_env)
591 TYPE(pwdft_environment_type), POINTER :: pwdft_env
592
593 INTEGER :: iatom, natom
594 REAL(KIND=c_double), DIMENSION(3) :: a1, a2, a3, v2
595 REAL(KIND=dp), DIMENSION(3) :: vr, vs
596 TYPE(cell_type), POINTER :: my_cell
597 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
598 TYPE(qs_subsys_type), POINTER :: qs_subsys
599 TYPE(sirius_context_handler) :: sctx
600 TYPE(sirius_ground_state_handler) :: gs_handler
601
602 cpassert(ASSOCIATED(pwdft_env))
603 CALL pwdft_env_get(pwdft_env, sctx=sctx, gs_handler=gs_handler)
604
605! get current positions and lattice vectors
606 CALL pwdft_env_get(pwdft_env=pwdft_env, qs_subsys=qs_subsys)
607
608! lattice vectors of the unit cell should be in [a.u.] (length is in [a.u.])
609 CALL qs_subsys_get(qs_subsys, cell=my_cell)
610 a1(:) = my_cell%hmat(:, 1)
611 a2(:) = my_cell%hmat(:, 2)
612 a3(:) = my_cell%hmat(:, 3)
613 CALL sirius_set_lattice_vectors(sctx, a1(1), a2(1), a3(1))
614
615! new atomic positions
616 CALL qs_subsys_get(qs_subsys, particle_set=particle_set)
617 natom = SIZE(particle_set)
618 DO iatom = 1, natom
619 vr(1:3) = particle_set(iatom)%r(1:3)
620 CALL real_to_scaled(vs, vr, my_cell)
621 v2(1:3) = vs(1:3)
622 CALL sirius_set_atom_position(sctx, iatom, v2(1))
623 END DO
624
625! update ground-state class
626 CALL sirius_update_ground_state(gs_handler)
627
628 CALL pwdft_env_set(pwdft_env, sctx=sctx, gs_handler=gs_handler)
629
630 END SUBROUTINE cp_sirius_update_context
631
632! **************************************************************************************************
633!> \brief ...
634!> \param sctx ...
635!> \param section ...
636!> \param section_name ...
637! **************************************************************************************************
638 SUBROUTINE cp_sirius_fill_in_section(sctx, section, section_name)
639 TYPE(sirius_context_handler), INTENT(INOUT) :: sctx
640 TYPE(section_vals_type), POINTER :: section
641 CHARACTER(*), INTENT(in) :: section_name
642
643 CHARACTER(len=256), TARGET :: option_name
644 CHARACTER(len=4096) :: description, usage
645 CHARACTER(len=80), DIMENSION(:), POINTER :: tmp
646 CHARACTER(len=80), TARGET :: str
647 INTEGER :: ctype, elem, ic, j
648 INTEGER, DIMENSION(:), POINTER :: ivals
649 INTEGER, TARGET :: enum_length, ival, length, &
650 num_possible_values, number_of_options
651 LOGICAL :: explicit
652 LOGICAL, DIMENSION(:), POINTER :: lvals
653 LOGICAL, TARGET :: found, lval
654 REAL(kind=dp), DIMENSION(:), POINTER :: rvals
655 REAL(kind=dp), TARGET :: rval
656
657 NULLIFY (rvals)
658 NULLIFY (ivals)
659 CALL sirius_option_get_section_length(section_name, number_of_options)
660
661 DO elem = 1, number_of_options
662 option_name = ''
663 CALL sirius_option_get_info(section_name, &
664 elem, &
665 option_name, &
666 256, &
667 ctype, &
668 num_possible_values, &
669 enum_length, &
670 description, &
671 4096, &
672 usage, &
673 4096)
674 IF ((option_name /= 'memory_usage') .AND. (option_name /= 'xc_functionals') .AND. (option_name /= 'vk')) THEN
675 CALL section_vals_val_get(section, option_name, explicit=found)
676 IF (found) THEN
677 SELECT CASE (ctype)
678 CASE (sirius_integer_type)
679 CALL section_vals_val_get(section, option_name, i_val=ival)
680 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(ival))
681 CASE (sirius_number_type)
682 CALL section_vals_val_get(section, option_name, r_val=rval)
683 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(rval))
684 CASE (sirius_logical_type)
685 CALL section_vals_val_get(section, option_name, l_val=lval)
686 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(lval))
687 CASE (sirius_string_type) ! string nightmare
688 str = ''
689 CALL section_vals_val_get(section, option_name, explicit=explicit, c_val=str)
690 str = trim(adjustl(str))
691 DO j = 1, len(str)
692 ic = ichar(str(j:j))
693 IF (ic >= 65 .AND. ic < 90) str(j:j) = char(ic + 32)
694 END DO
695
696 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(str), max_length=len_trim(str))
697 CASE (sirius_integer_array_type)
698 CALL section_vals_val_get(section, option_name, i_vals=ivals)
699 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(ivals(1)), &
700 max_length=num_possible_values)
701 CASE (sirius_number_array_type)
702 CALL section_vals_val_get(section, option_name, r_vals=rvals)
703 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(rvals(1)), &
704 max_length=num_possible_values)
705 CASE (sirius_logical_array_type)
706 CALL section_vals_val_get(section, option_name, l_vals=lvals)
707 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(lvals(1)), &
708 max_length=num_possible_values)
709 CASE (sirius_string_array_type)
710 CALL section_vals_val_get(section, option_name, explicit=explicit, n_rep_val=length)
711 DO j = 1, length
712 str = ''
713 CALL section_vals_val_get(section, option_name, i_rep_val=j, explicit=explicit, c_vals=tmp)
714 str = trim(adjustl(tmp(j)))
715 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(str), &
716 max_length=len_trim(str), append=.true.)
717 END DO
718 CASE DEFAULT
719 END SELECT
720 END IF
721 END IF
722 END DO
723 END SUBROUTINE cp_sirius_fill_in_section
724
725!***************************************************************************************************
726!> \brief ...
727!> \param pwdft_env ...
728!> \param calculate_forces ...
729!> \param calculate_stress_tensor ...
730!> \param
731!> \par History
732!> 07.2018 start the Sirius library
733!> \author JHU
734! **************************************************************************************************
735 SUBROUTINE cp_sirius_energy_force(pwdft_env, calculate_forces, calculate_stress_tensor)
736 TYPE(pwdft_environment_type), INTENT(INOUT), &
737 POINTER :: pwdft_env
738 LOGICAL, INTENT(IN) :: calculate_forces, calculate_stress_tensor
739
740 INTEGER :: iw, n1, n2
741 LOGICAL :: do_print, gs_converged
742 REAL(KIND=c_double) :: etotal
743 REAL(KIND=c_double), ALLOCATABLE, DIMENSION(:, :) :: cforces
744 REAL(KIND=c_double), DIMENSION(3, 3) :: cstress
745 REAL(KIND=dp), DIMENSION(3, 3) :: stress
746 REAL(KIND=dp), DIMENSION(:, :), POINTER :: forces
747 TYPE(cp_logger_type), POINTER :: logger
748 TYPE(pwdft_energy_type), POINTER :: energy
749 TYPE(section_vals_type), POINTER :: print_section, pwdft_input
750 TYPE(sirius_ground_state_handler) :: gs_handler
751
752 cpassert(ASSOCIATED(pwdft_env))
753
754 NULLIFY (logger)
755 logger => cp_get_default_logger()
757
758 CALL pwdft_env_get(pwdft_env=pwdft_env, gs_handler=gs_handler)
759 CALL sirius_find_ground_state(gs_handler, converged=gs_converged)
760
761 IF (gs_converged) THEN
762 IF (iw > 0) WRITE (iw, '(A)') "CP2K/SIRIUS: ground state is converged"
763 ELSE
764 IF (pwdft_env%ignore_convergence_failure) THEN
765 IF (iw > 0) WRITE (iw, '(A)') "CP2K/SIRIUS Warning: ground state is not converged"
766 ELSE
767 cpabort("CP2K/SIRIUS (ground state): SIRIUS did not converge.")
768 END IF
769 END IF
770 IF (iw > 0) CALL m_flush(iw)
771
772 CALL pwdft_env_get(pwdft_env=pwdft_env, energy=energy)
773 etotal = 0.0_c_double
774
775 CALL sirius_get_energy(gs_handler, 'band-gap', etotal)
776 energy%band_gap = etotal
777
778 etotal = 0.0_c_double
779 CALL sirius_get_energy(gs_handler, 'total', etotal)
780 energy%etotal = etotal
781
782 ! extract entropy (TS returned by sirius is always negative, sign
783 ! convention in QE)
784 etotal = 0.0_c_double
785 CALL sirius_get_energy(gs_handler, 'demet', etotal)
786 energy%entropy = -etotal
787
788 IF (calculate_forces) THEN
789 CALL pwdft_env_get(pwdft_env=pwdft_env, forces=forces)
790 n1 = SIZE(forces, 1)
791 n2 = SIZE(forces, 2)
792
793 ALLOCATE (cforces(n2, n1))
794 cforces = 0.0_c_double
795 CALL sirius_get_forces(gs_handler, 'total', cforces)
796 ! Sirius computes the forces but cp2k use the gradient everywhere
797 ! so a minus sign is needed.
798 ! note also that sirius and cp2k store the forces transpose to each other
799 ! sirius : forces(coordinates, atoms)
800 ! cp2k : forces(atoms, coordinates)
801 forces = -transpose(cforces(:, :))
802 DEALLOCATE (cforces)
803 END IF
804
805 IF (calculate_stress_tensor) THEN
806 cstress = 0.0_c_double
807 CALL sirius_get_stress_tensor(gs_handler, 'total', cstress)
808 stress(1:3, 1:3) = cstress(1:3, 1:3)
809 CALL pwdft_env_set(pwdft_env=pwdft_env, stress=stress)
810 END IF
811
812 CALL pwdft_env_get(pwdft_env=pwdft_env, pwdft_input=pwdft_input)
813 print_section => section_vals_get_subs_vals(pwdft_input, "PRINT")
814 CALL section_vals_get(print_section, explicit=do_print)
815 IF (do_print) THEN
816 CALL cp_sirius_print_results(pwdft_env, print_section)
817 END IF
818 END SUBROUTINE cp_sirius_energy_force
819
820!***************************************************************************************************
821!> \brief ...
822!> \param pwdft_env ...
823!> \param print_section ...
824!> \param
825!> \par History
826!> 12.2019 init
827!> \author JHU
828! **************************************************************************************************
829 SUBROUTINE cp_sirius_print_results(pwdft_env, print_section)
830 TYPE(pwdft_environment_type), INTENT(INOUT), &
831 POINTER :: pwdft_env
832 TYPE(section_vals_type), POINTER :: print_section
833
834 CHARACTER(LEN=default_string_length) :: my_act, my_pos
835 INTEGER :: i, ik, iounit, ispn, iterstep, iv, iw, &
836 nbands, nhist, nkpts, nspins
837 INTEGER(KIND=C_INT) :: cint
838 LOGICAL :: append, dos, ionode
839 REAL(KIND=c_double) :: creal
840 REAL(KIND=c_double), ALLOCATABLE, DIMENSION(:) :: slist
841 REAL(KIND=dp) :: de, e_fermi(2), emax, emin, eval
842 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: wkpt
843 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: ehist, hist, occval
844 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: energies, occupations
845 TYPE(cp_logger_type), POINTER :: logger
846 TYPE(sirius_context_handler) :: sctx
847 TYPE(sirius_ground_state_handler) :: gs_handler
848 TYPE(sirius_kpoint_set_handler) :: ks_handler
849
850 NULLIFY (logger)
851 logger => cp_get_default_logger()
852 ionode = logger%para_env%is_source()
853 iounit = cp_logger_get_default_io_unit(logger)
854
855 ! Density of States
856 dos = btest(cp_print_key_should_output(logger%iter_info, print_section, "DOS"), cp_p_file)
857 IF (dos) THEN
858 CALL pwdft_env_get(pwdft_env, ks_handler=ks_handler)
859 CALL pwdft_env_get(pwdft_env, gs_handler=gs_handler)
860 CALL pwdft_env_get(pwdft_env, sctx=sctx)
861
862 CALL section_vals_val_get(print_section, "DOS%DELTA_E", r_val=de)
863 CALL section_vals_val_get(print_section, "DOS%APPEND", l_val=append)
864
865 CALL sirius_get_num_kpoints(ks_handler, cint)
866 nkpts = cint
867 CALL sirius_get_parameters(sctx, num_bands=cint)
868 nbands = cint
869 CALL sirius_get_parameters(sctx, num_spins=cint)
870 nspins = cint
871 e_fermi(:) = 0.0_dp
872 ALLOCATE (energies(nbands, nspins, nkpts))
873 energies = 0.0_dp
874 ALLOCATE (occupations(nbands, nspins, nkpts))
875 occupations = 0.0_dp
876 ALLOCATE (wkpt(nkpts))
877 ALLOCATE (slist(nbands))
878 DO ik = 1, nkpts
879 CALL sirius_get_kpoint_properties(ks_handler, ik, creal)
880 wkpt(ik) = creal
881 END DO
882 DO ik = 1, nkpts
883 DO ispn = 1, nspins
884 CALL sirius_get_band_energies(ks_handler, ik, ispn, slist)
885 energies(1:nbands, ispn, ik) = slist(1:nbands)
886 CALL sirius_get_band_occupancies(ks_handler, ik, ispn, slist)
887 occupations(1:nbands, ispn, ik) = slist(1:nbands)
888 END DO
889 END DO
890 emin = minval(energies)
891 emax = maxval(energies)
892 nhist = nint((emax - emin)/de) + 1
893 ALLOCATE (hist(nhist, nspins), occval(nhist, nspins), ehist(nhist, nspins))
894 hist = 0.0_dp
895 occval = 0.0_dp
896 ehist = 0.0_dp
897
898 DO ik = 1, nkpts
899 DO ispn = 1, nspins
900 DO i = 1, nbands
901 eval = energies(i, ispn, ik) - emin
902 iv = nint(eval/de) + 1
903 cpassert((iv > 0) .AND. (iv <= nhist))
904 hist(iv, ispn) = hist(iv, ispn) + wkpt(ik)
905 occval(iv, ispn) = occval(iv, ispn) + wkpt(ik)*occupations(i, ispn, ik)
906 END DO
907 END DO
908 END DO
909 hist = hist/real(nbands, kind=dp)
910 DO i = 1, nhist
911 ehist(i, 1:nspins) = emin + (i - 1)*de
912 END DO
913
914 iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
915 my_act = "WRITE"
916 IF (append .AND. iterstep > 1) THEN
917 my_pos = "APPEND"
918 ELSE
919 my_pos = "REWIND"
920 END IF
921
922 iw = cp_print_key_unit_nr(logger, print_section, "DOS", &
923 extension=".dos", file_position=my_pos, file_action=my_act, &
924 file_form="FORMATTED")
925 IF (iw > 0) THEN
926 IF (nspins == 2) THEN
927 WRITE (unit=iw, fmt="(T2,A,I0,A,2F12.6)") &
928 "# DOS at iteration step i = ", iterstep, ", E_Fermi[a.u.] = ", e_fermi(1:2)
929 WRITE (unit=iw, fmt="(T2,A, A)") " Energy[a.u.] Alpha_Density Occupation", &
930 " Beta_Density Occupation"
931 ELSE
932 WRITE (unit=iw, fmt="(T2,A,I0,A,F12.6)") &
933 "# DOS at iteration step i = ", iterstep, ", E_Fermi[a.u.] = ", e_fermi(1)
934 WRITE (unit=iw, fmt="(T2,A)") " Energy[a.u.] Density Occupation"
935 END IF
936 DO i = 1, nhist
937 eval = emin + (i - 1)*de
938 IF (nspins == 2) THEN
939 WRITE (unit=iw, fmt="(F15.8,4F15.4)") eval, hist(i, 1), occval(i, 1), &
940 hist(i, 2), occval(i, 2)
941 ELSE
942 WRITE (unit=iw, fmt="(F15.8,2F15.4)") eval, hist(i, 1), occval(i, 1)
943 END IF
944 END DO
945 END IF
946 CALL cp_print_key_finished_output(iw, logger, print_section, "DOS")
947
948 DEALLOCATE (energies, occupations, wkpt, slist)
949 DEALLOCATE (hist, occval, ehist)
950 END IF
951 END SUBROUTINE cp_sirius_print_results
952
953END MODULE sirius_interface
954
955#else
956
957!***************************************************************************************************
958!> \brief Empty implementation in case SIRIUS is not compiled in.
959!***************************************************************************************************
962#include "./base/base_uses.f90"
963
964 IMPLICIT NONE
965 PRIVATE
966
969
970CONTAINS
971
972! **************************************************************************************************
973!> \brief Empty implementation in case SIRIUS is not compiled in.
974! **************************************************************************************************
975 SUBROUTINE cp_sirius_init()
976 END SUBROUTINE cp_sirius_init
977
978! **************************************************************************************************
979!> \brief Empty implementation in case SIRIUS is not compiled in.
980! **************************************************************************************************
982 END SUBROUTINE cp_sirius_finalize
983
984! **************************************************************************************************
985!> \brief Empty implementation in case SIRIUS is not compiled in.
986!> \param pwdft_env ...
987! **************************************************************************************************
988 SUBROUTINE cp_sirius_create_env(pwdft_env)
989 TYPE(pwdft_environment_type), POINTER :: pwdft_env
990
991 mark_used(pwdft_env)
992 cpabort("Sirius library is missing")
993 END SUBROUTINE cp_sirius_create_env
994
995! **************************************************************************************************
996!> \brief Empty implementation in case SIRIUS is not compiled in.
997!> \param pwdft_env ...
998!> \param calculate_forces ...
999!> \param calculate_stress ...
1000! **************************************************************************************************
1001 SUBROUTINE cp_sirius_energy_force(pwdft_env, calculate_forces, calculate_stress)
1002 TYPE(pwdft_environment_type), POINTER :: pwdft_env
1003 LOGICAL :: calculate_forces, calculate_stress
1004
1005 mark_used(pwdft_env)
1006 mark_used(calculate_forces)
1007 mark_used(calculate_stress)
1008 cpabort("Sirius library is missing")
1009 END SUBROUTINE cp_sirius_energy_force
1010
1011! **************************************************************************************************
1012!> \brief Empty implementation in case SIRIUS is not compiled in.
1013!> \param pwdft_env ...
1014! **************************************************************************************************
1015 SUBROUTINE cp_sirius_update_context(pwdft_env)
1016 TYPE(pwdft_environment_type), POINTER :: pwdft_env
1017
1018 mark_used(pwdft_env)
1019 cpabort("Sirius library is missing")
1020 END SUBROUTINE cp_sirius_update_context
1021
1022END MODULE sirius_interface
1023
1024#endif
calculate the orbitals for a given atomic kind type
subroutine, public calculate_atomic_orbitals(atomic_kind, qs_kind, agrid, iunit, pmat, fmat, density, wavefunction, wfninfo, confine, xc_section, nocc)
...
subroutine, public gth_potential_conversion(gth_potential, gth_atompot)
...
Define the atom type and its sub types.
Definition atom_types.F:15
Routines that process Quantum Espresso UPF files.
Definition atom_upf.F:14
Some basic routines for atomic calculations.
Definition atom_utils.F:15
pure subroutine, public atom_local_potential(locpot, gthpot, rr)
...
Definition atom_utils.F:799
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Handles all functions related to the CELL.
Definition cell_types.F:15
subroutine, public real_to_scaled(s, r, cell)
Transform real to scaled cell coordinates. s=h_inv*r.
Definition cell_types.F:486
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
integer, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
Definition of the atomic potential types.
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_gapw_log
integer, parameter, public sirius_func_vdwdfcx
integer, parameter, public sirius_func_vdwdf2
integer, parameter, public sirius_func_vdwdf
objects that represent the structure of input sections and the data contained in an input section
type(section_vals_type) function, pointer, public section_vals_get_subs_vals2(section_vals, i_section, i_rep_section)
returns the values of the n-th non default subsection (null if no such section exists (not so many no...
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition machine.F:131
Definition of mathematical constants and functions.
real(kind=dp), dimension(0:maxfac), parameter, public gamma1
real(kind=dp), parameter, public fourpi
Interface to the message passing library MPI.
Define the data structure for the particle information.
Definition of physical constants:
Definition physcon.F:68
real(kind=dp), parameter, public massunit
Definition physcon.F:141
The type definitions for the PWDFT environment.
subroutine, public pwdft_env_get(pwdft_env, pwdft_input, force_env_input, xc_input, cp_subsys, qs_subsys, para_env, energy, forces, stress, sctx, gs_handler, ks_handler)
Returns various attributes of the pwdft environment.
subroutine, public pwdft_env_set(pwdft_env, pwdft_input, force_env_input, xc_input, qs_subsys, cp_subsys, para_env, energy, forces, stress, sctx, gs_handler, ks_handler)
Sets various attributes of the pwdft environment.
subroutine, public deallocate_grid_atom(grid_atom)
Deallocate a Gaussian-type orbital (GTO) basis set data set.
subroutine, public allocate_grid_atom(grid_atom)
Initialize components of the grid_atom_type structure.
subroutine, public create_grid_atom(grid_atom, nr, na, llmax, ll, quadrature)
...
Define the quickstep kind type and their sub types.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
types that represent a quickstep subsys
subroutine, public qs_subsys_get(subsys, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, cp_subsys, nelectron_total, nelectron_spin)
...
Interface to the SIRIUS Library.
subroutine, public cp_sirius_update_context(pwdft_env)
Empty implementation in case SIRIUS is not compiled in.
subroutine, public cp_sirius_init()
Empty implementation in case SIRIUS is not compiled in.
subroutine, public cp_sirius_energy_force(pwdft_env, calculate_forces, calculate_stress)
Empty implementation in case SIRIUS is not compiled in.
subroutine, public cp_sirius_finalize()
Empty implementation in case SIRIUS is not compiled in.
subroutine, public cp_sirius_create_env(pwdft_env)
Empty implementation in case SIRIUS is not compiled in.
Provides all information about a pseudopotential.
Definition atom_types.F:98
Provides all information about an atomic kind.
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment
Provides all information about a quickstep kind.