16 USE iso_c_binding,
ONLY: c_double,&
55 pwdft_environment_type
65 sirius_integer_array_type, sirius_integer_type, sirius_logical_array_type, &
66 sirius_logical_type, sirius_number_array_type, sirius_number_type, &
67 sirius_string_array_type, sirius_string_type, sirius_add_atom, sirius_add_atom_type, &
68 sirius_add_atom_type_radial_function, sirius_add_xc_functional, sirius_context_handler, &
69 sirius_create_context, sirius_create_ground_state, sirius_create_kset_from_grid, &
70 sirius_finalize, sirius_find_ground_state, sirius_get_band_energies, &
71 sirius_get_band_occupancies, sirius_get_energy, sirius_get_forces, &
72 sirius_get_kpoint_properties, sirius_get_num_kpoints, sirius_get_parameters, &
73 sirius_get_stress_tensor, sirius_ground_state_handler, sirius_import_parameters, &
74 sirius_initialize, sirius_initialize_context, sirius_kpoint_set_handler, &
75 sirius_option_get_info, sirius_option_get_section_length, sirius_option_set, &
76 sirius_set_atom_position, sirius_set_atom_type_dion, sirius_set_atom_type_hubbard, &
77 sirius_set_atom_type_radial_grid, sirius_set_lattice_vectors, sirius_set_mpi_grid_dims, &
78 sirius_update_ground_state
79 #include "./base/base_uses.f90"
87 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'sirius_interface'
104 CALL sirius_initialize(.false.)
115 CALL sirius_finalize(.false., .false., .false.)
127 TYPE(pwdft_environment_type),
POINTER :: pwdft_env
128 #if defined(__SIRIUS)
130 CHARACTER(len=2) :: element_symbol
131 CHARACTER(len=default_string_length) :: label
132 INTEGER :: i, iatom, ibeta, ifun, ikind, iwf, j, l, &
133 n, natom, nbeta, nkind, nmesh, &
134 num_mag_dims, sirius_mpi_comm, vdw_func, nu, lu
135 INTEGER,
DIMENSION(:),
POINTER :: mpi_grid_dims
136 INTEGER(KIND=C_INT),
DIMENSION(3) :: k_grid, k_shift
137 INTEGER,
DIMENSION(:),
POINTER :: kk
138 LOGICAL :: up, use_ref_cell
139 LOGICAL(4) :: use_so, use_symmetry, dft_plus_u_atom
140 REAL(KIND=c_double),
ALLOCATABLE,
DIMENSION(:) :: fun
141 REAL(KIND=c_double),
ALLOCATABLE,
DIMENSION(:, :) :: dion
142 REAL(KIND=c_double),
DIMENSION(3) :: a1, a2, a3, v1, v2
143 REAL(KIND=
dp) :: al, angle1, angle2, cval, focc, &
144 magnetization, mass, pf, rl, zeff, alpha_u, beta_u, &
145 j0_u, j_u, u_u, occ_u, u_minus_j
146 REAL(KIND=
dp),
ALLOCATABLE,
DIMENSION(:) :: beta, corden, ef, fe, locpot, rc, rp
147 REAL(KIND=
dp),
DIMENSION(3) :: vr, vs
148 REAL(KIND=
dp),
DIMENSION(:),
POINTER :: density
149 REAL(KIND=
dp),
DIMENSION(:, :),
POINTER :: wavefunction, wfninfo
150 TYPE(atom_gthpot_type),
POINTER :: gth_atompot
151 TYPE(atom_upfpot_type),
POINTER :: upf_pot
152 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
153 TYPE(atomic_kind_type),
POINTER :: atomic_kind
154 TYPE(cell_type),
POINTER :: my_cell
155 TYPE(mp_para_env_type),
POINTER :: para_env
156 TYPE(grid_atom_type),
POINTER :: atom_grid
157 TYPE(gth_potential_type),
POINTER :: gth_potential
158 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
159 TYPE(qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
160 TYPE(qs_subsys_type),
POINTER :: qs_subsys
161 TYPE(section_vals_type),
POINTER :: pwdft_section, pwdft_sub_section, &
163 TYPE(sirius_context_handler) :: sctx
164 TYPE(sirius_ground_state_handler) :: gs_handler
165 TYPE(sirius_kpoint_set_handler) :: ks_handler
167 cpassert(
ASSOCIATED(pwdft_env))
171 sirius_mpi_comm = para_env%get_handle()
172 CALL sirius_create_context(sirius_mpi_comm, sctx)
176 CALL pwdft_env_get(pwdft_env=pwdft_env, pwdft_input=pwdft_section, xc_input=xc_section)
181 IF (
ASSOCIATED(xc_section))
THEN
186 IF (.NOT.
ASSOCIATED(xc_fun))
EXIT
189 CALL sirius_add_xc_functional(sctx,
"XC_"//trim(xc_fun%section%name))
195 IF (
ASSOCIATED(pwdft_sub_section))
THEN
196 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section,
"control")
203 IF (
ASSOCIATED(pwdft_sub_section))
THEN
204 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section,
"parameters")
222 SELECT CASE (vdw_func)
224 CALL sirius_add_xc_functional(sctx,
"XC_FUNC_VDWDF")
226 CALL sirius_add_xc_functional(sctx,
"XC_FUNC_VDWDF2")
228 CALL sirius_add_xc_functional(sctx,
"XC_FUNC_VDWDF2")
237 IF (
ASSOCIATED(pwdft_sub_section))
THEN
238 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section,
"mixer")
243 IF (
ASSOCIATED(pwdft_sub_section))
THEN
244 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section,
"settings")
249 IF (
ASSOCIATED(pwdft_sub_section))
THEN
250 CALL cp_sirius_fill_in_section(sctx, pwdft_sub_section,
"iterative_solver")
264 CALL sirius_import_parameters(sctx,
'{}')
268 CALL qs_subsys_get(qs_subsys, cell=my_cell, use_ref_cell=use_ref_cell)
269 a1(:) = my_cell%hmat(:, 1)
270 a2(:) = my_cell%hmat(:, 2)
271 a3(:) = my_cell%hmat(:, 3)
272 CALL sirius_set_lattice_vectors(sctx, a1(1), a2(1), a3(1))
274 IF (use_ref_cell)
THEN
275 cpwarn(
"SIRIUS| The specified CELL_REF will be ignored for PW_DFT calculations")
280 atomic_kind_set=atomic_kind_set, &
281 qs_kind_set=qs_kind_set, &
282 particle_set=particle_set)
283 nkind =
SIZE(atomic_kind_set)
286 name=label, element_symbol=element_symbol, mass=mass)
288 NULLIFY (upf_pot, gth_potential)
289 CALL get_qs_kind(qs_kind_set(ikind), upf_potential=upf_pot, gth_potential=gth_potential)
291 IF (
ASSOCIATED(upf_pot))
THEN
292 CALL sirius_add_atom_type(sctx, label, fname=upf_pot%filename, &
293 symbol=element_symbol, &
294 mass=real(mass/
massunit, kind=c_double))
296 ELSEIF (
ASSOCIATED(gth_potential))
THEN
303 ALLOCATE (rp(nmesh), fun(nmesh))
304 IF (atom_grid%rad(1) < atom_grid%rad(nmesh))
THEN
310 rp(1:nmesh) = atom_grid%rad(1:nmesh)
313 rp(i) = atom_grid%rad(nmesh - i + 1)
317 CALL sirius_add_atom_type(sctx, label, &
318 zn=nint(zeff + 0.001d0), &
319 symbol=element_symbol, &
320 mass=real(mass/
massunit, kind=c_double), &
323 ALLOCATE (gth_atompot)
326 fun(1:nmesh) = rp(1:nmesh)
327 CALL sirius_set_atom_type_radial_grid(sctx, label, nmesh, fun(1))
329 ALLOCATE (ef(nmesh), beta(nmesh))
332 IF (gth_atompot%nl(l) == 0) cycle
333 rl = gth_atompot%rcnl(l)
335 ef(1:nmesh) = exp(-0.5_dp*rp(1:nmesh)*rp(1:nmesh)/(rl*rl))
336 DO i = 1, gth_atompot%nl(l)
337 pf = rl**(l + 0.5_dp*(4._dp*i - 1._dp))
339 pf = sqrt(2._dp)/(pf*sqrt(
gamma1(j)))
340 beta(:) = pf*rp**(l + 2*i - 2)*ef
342 fun(1:nmesh) = beta(1:nmesh)*rp(1:nmesh)
343 CALL sirius_add_atom_type_radial_function(sctx, label, &
344 "beta", fun(1), nmesh, l=l)
347 DEALLOCATE (ef, beta)
351 ALLOCATE (dion(nbeta, nbeta))
354 IF (gth_atompot%nl(l) == 0) cycle
355 ibeta = sum(gth_atompot%nl(0:l - 1)) + 1
356 i = ibeta + gth_atompot%nl(l) - 1
357 dion(ibeta:i, ibeta:i) = gth_atompot%hnl(1:gth_atompot%nl(l), 1:gth_atompot%nl(l), l)
359 CALL sirius_set_atom_type_dion(sctx, label, nbeta, dion(1, 1))
363 IF (gth_atompot%nlcc)
THEN
364 ALLOCATE (corden(nmesh), fe(nmesh), rc(nmesh))
366 n = gth_atompot%nexp_nlcc
368 al = gth_atompot%alpha_nlcc(i)
370 fe(:) = exp(-0.5_dp*rc(:)*rc(:))
371 DO j = 1, gth_atompot%nct_nlcc(i)
372 cval = gth_atompot%cval_nlcc(j, i)
373 corden(:) = corden(:) + fe(:)*rc(:)**(2*j - 2)*cval
376 fun(1:nmesh) = corden(1:nmesh)*rp(1:nmesh)
377 CALL sirius_add_atom_type_radial_function(sctx, label,
"ps_rho_core", &
379 DEALLOCATE (corden, fe, rc)
383 ALLOCATE (locpot(nmesh))
386 fun(1:nmesh) = locpot(1:nmesh)
387 CALL sirius_add_atom_type_radial_function(sctx, label,
"vloc", &
391 NULLIFY (density, wavefunction, wfninfo)
393 density=density, wavefunction=wavefunction, &
394 wfninfo=wfninfo, agrid=atom_grid)
397 DO iwf = 1,
SIZE(wavefunction, 2)
398 focc = wfninfo(1, iwf)
399 l = nint(wfninfo(2, iwf))
403 fun(1:nmesh) = wavefunction(1:nmesh, iwf)*rp(i)
406 fun(i) = wavefunction(nmesh - i + 1, iwf)*rp(i)
409 CALL sirius_add_atom_type_radial_function(sctx, &
410 label,
"ps_atomic_wf", &
411 fun(1), nmesh, l=l, occ=real(focc, kind=c_double), n=nu)
416 fun(1:nmesh) =
fourpi*density(1:nmesh)*atom_grid%rad(1:nmesh)**2
419 fun(i) =
fourpi*density(nmesh - i + 1)*atom_grid%rad(nmesh - i + 1)**2
422 CALL sirius_add_atom_type_radial_function(sctx, label,
"ps_rho_total", &
425 IF (
ASSOCIATED(density))
DEALLOCATE (density)
426 IF (
ASSOCIATED(wavefunction))
DEALLOCATE (wavefunction)
427 IF (
ASSOCIATED(wfninfo))
DEALLOCATE (wfninfo)
431 DEALLOCATE (gth_atompot)
434 CALL cp_abort(__location__, &
435 "CP2K/SIRIUS: atomic kind needs UPF or GTH potential definition")
439 dft_plus_u_atom=dft_plus_u_atom, &
440 l_of_dft_plus_u=lu, &
441 n_of_dft_plus_u=nu, &
442 u_minus_j_target=u_minus_j, &
443 u_of_dft_plus_u=u_u, &
444 j_of_dft_plus_u=j_u, &
445 alpha_of_dft_plus_u=alpha_u, &
446 beta_of_dft_plus_u=beta_u, &
447 j0_of_dft_plus_u=j0_u, &
448 occupation_of_dft_plus_u=occ_u)
450 IF (dft_plus_u_atom)
THEN
452 cpabort(
"CP2K/SIRIUS (hubbard): principal quantum number not specified")
456 cpabort(
"CP2K/SIRIUS (hubbard): l can not be negative.")
459 IF (occ_u < 0.0)
THEN
460 cpabort(
"CP2K/SIRIUS (hubbard): the occupation number can not be negative.")
463 IF (abs(u_minus_j) < 1e-8)
THEN
464 CALL sirius_set_atom_type_hubbard(sctx, label, lu, nu, &
465 occ_u, u_u, j_u, alpha_u, beta_u, j0_u)
467 CALL sirius_set_atom_type_hubbard(sctx, label, lu, nu, &
468 occ_u, u_minus_j, 0.0_dp, alpha_u, beta_u, j0_u)
476 natom =
SIZE(particle_set)
478 vr(1:3) = particle_set(iatom)%r(1:3)
480 atomic_kind => particle_set(iatom)%atomic_kind
481 ikind = atomic_kind%kind_number
483 CALL get_qs_kind(qs_kind_set(ikind), zeff=zeff, magnetization=magnetization)
487 IF (num_mag_dims .EQ. 3)
THEN
490 v1(1) = magnetization*sin(angle1)*cos(angle2)
491 v1(2) = magnetization*sin(angle1)*sin(angle2)
492 v1(3) = magnetization*cos(angle1)
495 v1(3) = magnetization
498 CALL sirius_add_atom(sctx, label, v2(1), v1(1))
501 CALL sirius_set_mpi_grid_dims(sctx, 2, mpi_grid_dims)
504 CALL sirius_initialize_context(sctx)
508 IF (use_symmetry)
THEN
509 CALL sirius_create_kset_from_grid(sctx, k_grid(1), k_shift(1), use_symmetry=.true., kset_handler=ks_handler)
511 CALL sirius_create_kset_from_grid(sctx, k_grid(1), k_shift(1), use_symmetry=.false., kset_handler=ks_handler)
514 CALL sirius_create_ground_state(ks_handler, gs_handler)
516 CALL pwdft_env_set(pwdft_env, sctx=sctx, gs_handler=gs_handler, ks_handler=ks_handler)
529 TYPE(pwdft_environment_type),
POINTER :: pwdft_env
531 INTEGER :: iatom, natom
532 REAL(KIND=c_double),
DIMENSION(3) :: a1, a2, a3, v2
533 REAL(KIND=
dp),
DIMENSION(3) :: vr, vs
534 TYPE(cell_type),
POINTER :: my_cell
535 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
536 TYPE(qs_subsys_type),
POINTER :: qs_subsys
537 TYPE(sirius_context_handler) :: sctx
538 TYPE(sirius_ground_state_handler) :: gs_handler
540 cpassert(
ASSOCIATED(pwdft_env))
541 CALL pwdft_env_get(pwdft_env, sctx=sctx, gs_handler=gs_handler)
548 a1(:) = my_cell%hmat(:, 1)
549 a2(:) = my_cell%hmat(:, 2)
550 a3(:) = my_cell%hmat(:, 3)
551 CALL sirius_set_lattice_vectors(sctx, a1(1), a2(1), a3(1))
555 natom =
SIZE(particle_set)
557 vr(1:3) = particle_set(iatom)%r(1:3)
560 CALL sirius_set_atom_position(sctx, iatom, v2(1))
564 CALL sirius_update_ground_state(gs_handler)
566 CALL pwdft_env_set(pwdft_env, sctx=sctx, gs_handler=gs_handler)
576 SUBROUTINE cp_sirius_fill_in_section(sctx, section, section_name)
577 TYPE(sirius_context_handler),
INTENT(INOUT) :: sctx
578 TYPE(section_vals_type),
POINTER :: section
579 CHARACTER(*),
INTENT(in) :: section_name
581 CHARACTER(len=256),
TARGET :: option_name
582 CHARACTER(len=4096) :: description, usage
583 CHARACTER(len=80),
DIMENSION(:),
POINTER :: tmp
584 CHARACTER(len=80),
TARGET :: str
585 INTEGER :: ctype, elem, ic, j
586 INTEGER,
DIMENSION(:),
POINTER :: ivals
587 INTEGER,
TARGET :: enum_length, ival, length, &
588 num_possible_values, number_of_options
590 LOGICAL,
DIMENSION(:),
POINTER :: lvals
591 LOGICAL,
TARGET :: found, lval
592 REAL(kind=
dp),
DIMENSION(:),
POINTER :: rvals
593 REAL(kind=
dp),
TARGET :: rval
597 CALL sirius_option_get_section_length(section_name, number_of_options)
599 DO elem = 1, number_of_options
601 CALL sirius_option_get_info(section_name, &
606 num_possible_values, &
612 IF ((option_name /=
'memory_usage') .AND. (option_name /=
'xc_functionals') .AND. (option_name /=
'vk'))
THEN
616 CASE (sirius_integer_type)
618 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(ival))
619 CASE (sirius_number_type)
621 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(rval))
622 CASE (sirius_logical_type)
624 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(lval))
625 CASE (sirius_string_type)
628 str = trim(adjustl(str))
631 IF (ic >= 65 .AND. ic < 90) str(j:j) = char(ic + 32)
634 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(str), max_length=len_trim(str))
635 CASE (sirius_integer_array_type)
637 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(ivals(1)), &
638 max_length=num_possible_values)
639 CASE (sirius_number_array_type)
641 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(rvals(1)), &
642 max_length=num_possible_values)
643 CASE (sirius_logical_array_type)
645 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(lvals(1)), &
646 max_length=num_possible_values)
647 CASE (sirius_string_array_type)
652 str = trim(adjustl(tmp(j)))
653 CALL sirius_option_set(sctx, section_name, option_name, ctype, c_loc(str), &
654 max_length=len_trim(str), append=.true.)
661 END SUBROUTINE cp_sirius_fill_in_section
674 TYPE(pwdft_environment_type),
INTENT(INOUT), &
676 LOGICAL,
INTENT(IN) :: calculate_forces, calculate_stress_tensor
680 REAL(KIND=c_double) :: etotal
681 REAL(KIND=c_double),
ALLOCATABLE,
DIMENSION(:, :) :: cforces
682 REAL(KIND=c_double),
DIMENSION(3, 3) :: cstress
683 REAL(KIND=
dp),
DIMENSION(3, 3) :: stress
684 REAL(KIND=
dp),
DIMENSION(:, :),
POINTER :: forces
685 TYPE(pwdft_energy_type),
POINTER :: energy
686 TYPE(section_vals_type),
POINTER :: print_section, pwdft_input
687 TYPE(sirius_ground_state_handler) :: gs_handler
689 cpassert(
ASSOCIATED(pwdft_env))
691 CALL pwdft_env_get(pwdft_env=pwdft_env, gs_handler=gs_handler)
692 CALL sirius_find_ground_state(gs_handler)
694 etotal = 0.0_c_double
696 CALL sirius_get_energy(gs_handler,
'band-gap', etotal)
697 energy%band_gap = etotal
699 etotal = 0.0_c_double
700 CALL sirius_get_energy(gs_handler,
'total', etotal)
701 energy%etotal = etotal
704 etotal = 0.0_c_double
705 CALL sirius_get_energy(gs_handler,
'demet', etotal)
706 energy%entropy = etotal
708 IF (calculate_forces)
THEN
713 ALLOCATE (cforces(n2, n1))
714 cforces = 0.0_c_double
715 CALL sirius_get_forces(gs_handler,
'total', cforces)
721 forces = -transpose(cforces(:, :))
725 IF (calculate_stress_tensor)
THEN
726 cstress = 0.0_c_double
727 CALL sirius_get_stress_tensor(gs_handler,
'total', cstress)
728 stress(1:3, 1:3) = cstress(1:3, 1:3)
732 CALL pwdft_env_get(pwdft_env=pwdft_env, pwdft_input=pwdft_input)
736 CALL cp_sirius_print_results(pwdft_env, print_section)
750 SUBROUTINE cp_sirius_print_results(pwdft_env, print_section)
751 TYPE(pwdft_environment_type),
INTENT(INOUT), &
753 TYPE(section_vals_type),
POINTER :: print_section
755 CHARACTER(LEN=default_string_length) :: my_act, my_pos
756 INTEGER :: i, ik, iounit, ispn, iterstep, iv, iw, &
757 nbands, nhist, nkpts, nspins
758 INTEGER(KIND=C_INT) :: cint
759 LOGICAL :: append, dos, ionode
760 REAL(KIND=c_double) :: creal
761 REAL(KIND=c_double),
ALLOCATABLE,
DIMENSION(:) :: slist
762 REAL(KIND=
dp) :: de, e_fermi(2), emax, emin, eval
763 REAL(KIND=
dp),
ALLOCATABLE,
DIMENSION(:) :: wkpt
764 REAL(KIND=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: ehist, hist, occval
765 REAL(KIND=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: energies, occupations
766 TYPE(cp_logger_type),
POINTER :: logger
767 TYPE(sirius_context_handler) :: sctx
768 TYPE(sirius_ground_state_handler) :: gs_handler
769 TYPE(sirius_kpoint_set_handler) :: ks_handler
773 ionode = logger%para_env%is_source()
786 CALL sirius_get_num_kpoints(ks_handler, cint)
788 CALL sirius_get_parameters(sctx, num_bands=cint)
790 CALL sirius_get_parameters(sctx, num_spins=cint)
793 ALLOCATE (energies(nbands, nspins, nkpts))
795 ALLOCATE (occupations(nbands, nspins, nkpts))
797 ALLOCATE (wkpt(nkpts))
798 ALLOCATE (slist(nbands))
800 CALL sirius_get_kpoint_properties(ks_handler, ik, creal)
805 CALL sirius_get_band_energies(ks_handler, ik, ispn, slist)
806 energies(1:nbands, ispn, ik) = slist(1:nbands)
807 CALL sirius_get_band_occupancies(ks_handler, ik, ispn, slist)
808 occupations(1:nbands, ispn, ik) = slist(1:nbands)
811 emin = minval(energies)
812 emax = maxval(energies)
813 nhist = nint((emax - emin)/de) + 1
814 ALLOCATE (hist(nhist, nspins), occval(nhist, nspins), ehist(nhist, nspins))
822 eval = energies(i, ispn, ik) - emin
823 iv = nint(eval/de) + 1
824 cpassert((iv > 0) .AND. (iv <= nhist))
825 hist(iv, ispn) = hist(iv, ispn) + wkpt(ik)
826 occval(iv, ispn) = occval(iv, ispn) + wkpt(ik)*occupations(i, ispn, ik)
830 hist = hist/real(nbands, kind=
dp)
832 ehist(i, 1:nspins) = emin + (i - 1)*de
835 iterstep = logger%iter_info%iteration(logger%iter_info%n_rlevel)
837 IF (append .AND. iterstep > 1)
THEN
844 extension=
".dos", file_position=my_pos, file_action=my_act, &
845 file_form=
"FORMATTED")
847 IF (nspins == 2)
THEN
848 WRITE (unit=iw, fmt=
"(T2,A,I0,A,2F12.6)") &
849 "# DOS at iteration step i = ", iterstep,
", E_Fermi[a.u.] = ", e_fermi(1:2)
850 WRITE (unit=iw, fmt=
"(T2,A, A)")
" Energy[a.u.] Alpha_Density Occupation", &
851 " Beta_Density Occupation"
853 WRITE (unit=iw, fmt=
"(T2,A,I0,A,F12.6)") &
854 "# DOS at iteration step i = ", iterstep,
", E_Fermi[a.u.] = ", e_fermi(1)
855 WRITE (unit=iw, fmt=
"(T2,A)")
" Energy[a.u.] Density Occupation"
858 eval = emin + (i - 1)*de
859 IF (nspins == 2)
THEN
860 WRITE (unit=iw, fmt=
"(F15.8,4F15.4)") eval, hist(i, 1), occval(i, 1), &
861 hist(i, 2), occval(i, 2)
863 WRITE (unit=iw, fmt=
"(F15.8,2F15.4)") eval, hist(i, 1), occval(i, 1)
869 DEALLOCATE (energies, occupations, wkpt, slist)
870 DEALLOCATE (hist, occval, ehist)
872 END SUBROUTINE cp_sirius_print_results
883 #include "./base/base_uses.f90"
910 TYPE(pwdft_environment_type),
POINTER :: pwdft_env
913 cpabort(
"Sirius library is missing")
923 TYPE(pwdft_environment_type),
POINTER :: pwdft_env
924 LOGICAL :: calculate_forces, calculate_stress
927 mark_used(calculate_forces)
928 mark_used(calculate_stress)
929 cpabort(
"Sirius library is missing")
937 TYPE(pwdft_environment_type),
POINTER :: pwdft_env
940 cpabort(
"Sirius library is missing")
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.
Routines that process Quantum Espresso UPF files.
Some basic routines for atomic calculations.
pure subroutine, public atom_local_potential(locpot, gthpot, rr)
...
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.
subroutine, public real_to_scaled(s, r, cell)
Transform real to scaled cell coordinates. s=h_inv*r.
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.
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
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:
real(kind=dp), parameter, public massunit
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, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, 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_r3d_rs_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_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.