95 SUBROUTINE put_rho0_on_grid(qs_env, rho0, tot_rs_int, my_pools, my_rs_grids, my_rs_descs)
99 REAL(kind=
dp),
INTENT(OUT) :: tot_rs_int
103 OPTIONAL,
POINTER :: my_rs_grids
105 OPTIONAL,
POINTER :: my_rs_descs
107 CHARACTER(LEN=*),
PARAMETER :: routinen =
'put_rho0_on_grid'
109 INTEGER :: auxbas_grid, handle, iat, iatom, igrid, &
110 ikind, ithread, j, l0_ikind, lmax0, &
111 nat, nch_ik, nch_max, npme
112 INTEGER,
DIMENSION(:),
POINTER :: atom_list, cores
114 REAL(kind=
dp) :: eps_rho_rspace, rpgf0, zet0
115 REAL(kind=
dp),
DIMENSION(3) :: ra
116 REAL(kind=
dp),
DIMENSION(:),
POINTER :: qlm_c
117 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: pab
130 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
137 CALL timeset(routinen, handle)
139 NULLIFY (atomic_kind_set, qs_kind_set, cores, pab, qlm_c)
141 NULLIFY (dft_control, pw_env, particle_set, para_env, cell, rho0_s_gs, rho0_s_rs)
142 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control, &
143 particle_set=particle_set, &
144 atomic_kind_set=atomic_kind_set, &
145 qs_kind_set=qs_kind_set, &
147 pw_env=pw_env, cell=cell)
148 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
150 NULLIFY (descs, pw_pools)
151 CALL pw_env_get(pw_env=pw_env, rs_descs=descs, rs_grids=grids, pw_pools=pw_pools)
152 auxbas_grid = pw_env%auxbas_grid
154 NULLIFY (rho0_s_gs, rho0_s_rs)
156 zet0_h=zet0, igrid_zet0_s=igrid, &
157 rho0_s_gs=rho0_s_gs, &
161 NULLIFY (rs_grid, desc, pw_pool)
163 IF (
PRESENT(my_pools))
THEN
164 desc => my_rs_descs(igrid)%rs_desc
165 rs_grid => my_rs_grids(igrid)
166 pw_pool => my_pools(igrid)%pool
168 desc => descs(igrid)%rs_desc
169 rs_grid => grids(igrid)
170 pw_pool => pw_pools(igrid)%pool
173 cpassert(
ASSOCIATED(desc))
174 cpassert(
ASSOCIATED(pw_pool))
176 IF (igrid /= auxbas_grid)
THEN
177 CALL pw_pool%create_pw(coeff_rspace)
178 CALL pw_pool%create_pw(coeff_gspace)
184 ALLOCATE (pab(nch_max, 1))
186 DO ikind = 1,
SIZE(atomic_kind_set)
187 CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom_list, natom=nat)
188 CALL get_qs_kind(qs_kind_set(ikind), paw_atom=paw_atom)
190 IF (.NOT. paw_atom .AND. dft_control%qs_control%gapw_control%nopaw_as_gpw) cycle
192 CALL get_rho0_mpole(rho0_mpole=rho0, ikind=ikind, l0_ikind=l0_ikind, &
203 iatom = atom_list(iat)
204 ra(:) =
pbc(particle_set(iatom)%r, cell)
205 IF (rs_grid%desc%parallel .AND. .NOT. rs_grid%desc%distributed)
THEN
207 IF (
modulo(nat, rs_grid%desc%group_size) == rs_grid%desc%my_pos)
THEN
222 iatom = atom_list(iat)
226 pab(1:nch_ik, 1) = qlm_c(1:nch_ik)
228 ra(:) =
pbc(particle_set(iatom)%r, cell)
231 l0_ikind, zet0, 0, 0, 0.0_dp, 0, &
232 ra, (/0.0_dp, 0.0_dp, 0.0_dp/), 1.0_dp, pab, 0, 0, &
234 use_subpatch=.true., subpatch_pattern=0)
240 IF (
ASSOCIATED(cores))
THEN
246 IF (igrid /= auxbas_grid)
THEN
250 CALL pw_axpy(coeff_gspace, rho0_s_gs)
254 CALL pw_pool%give_back_pw(coeff_rspace)
255 CALL pw_pool%give_back_pw(coeff_gspace)
258 CALL pw_pool%create_pw(rho0_r_tmp)
265 CALL pw_pool%give_back_pw(rho0_r_tmp)
270 CALL timestop(handle)
341 local_rho_set_2nd, atener, kforce, my_pools, my_rs_descs)
346 LOGICAL,
INTENT(IN) :: calculate_forces
347 TYPE(
local_rho_type),
OPTIONAL,
POINTER :: local_rho_set, local_rho_set_2nd
348 REAL(kind=
dp),
DIMENSION(:),
OPTIONAL :: atener
349 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: kforce
353 OPTIONAL,
POINTER :: my_rs_descs
355 CHARACTER(LEN=*),
PARAMETER :: routinen =
'integrate_vhg0_rspace'
357 INTEGER :: auxbas_grid, bo(2), handle, i, iat, iatom, ic, icg, ico, ig1, ig2, igrid, ii, &
358 ikind, ipgf1, ipgf2, is, iset1, iset2, iso, iso1, iso2, ispin, j, l0_ikind, llmax, &
359 llmax_nuc, lmax0, lshell, lx, ly, lz, m1, m2, max_iso_not0_local, max_s_harm, &
360 max_s_harm_nuc, maxl, maxl_nuc, maxso, maxso_nuc, mepos, n1, n2, nat, nch_ik, nch_max, &
361 ncurr, nset, nset_nuc, nsotot, nsotot_nuc, nspins, num_pe
362 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: cg_n_list
363 INTEGER,
ALLOCATABLE,
DIMENSION(:, :, :) :: cg_list
364 INTEGER,
DIMENSION(:),
POINTER :: atom_list, lmax, lmax_nuc, lmin, &
365 lmin_nuc, npgf, npgf_nuc
366 LOGICAL :: grid_distributed, paw_atom, use_virial
367 REAL(kind=
dp) :: eps_rho_rspace, force_tmp(3), fscale, &
369 REAL(kind=
dp),
DIMENSION(3, 3) :: my_virial_a, my_virial_b
370 REAL(kind=
dp),
DIMENSION(:),
POINTER :: hab_sph, norm_l, qlm
371 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: hab, hdab_sph, intloc, intloc_nuc, pab
372 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: a_hdab_sph, hdab, qlm_gg
373 REAL(kind=
dp),
DIMENSION(:, :, :, :),
POINTER :: a_hdab
387 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
393 TYPE(
rho_atom_coeff),
DIMENSION(:),
POINTER :: int_local_h, int_local_s
400 CALL timeset(routinen, handle)
411 NULLIFY (atomic_kind_set, qs_kind_set, dft_control, particle_set)
412 NULLIFY (cell, force, pw_env, rho0_mpole, rho_atom_set, rhoz_cneo_set)
415 atomic_kind_set=atomic_kind_set, &
416 qs_kind_set=qs_kind_set, &
418 dft_control=dft_control, &
419 force=force, pw_env=pw_env, &
420 rho0_mpole=rho0_mpole, &
421 rho_atom_set=rho_atom_set, &
422 rhoz_cneo_set=rhoz_cneo_set, &
423 particle_set=particle_set, &
426 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
428 nspins = dft_control%nspins
443 IF (
PRESENT(local_rho_set)) &
444 CALL get_local_rho(local_rho_set, rho0_mpole=rho0_mpole, rho_atom_set=rho_atom_set, &
445 rhoz_cneo_set=rhoz_cneo_set)
449 IF (
PRESENT(local_rho_set_2nd))
THEN
450 CALL get_local_rho(local_rho_set_2nd, rho_atom_set=rho_atom_set)
453 zet0_h=zet0, igrid_zet0_s=igrid, &
457 NULLIFY (rs_descs, pw_pools)
458 cpassert(
ASSOCIATED(pw_env))
459 CALL pw_env_get(pw_env, rs_descs=rs_descs, pw_pools=pw_pools)
462 auxbas_grid = pw_env%auxbas_grid
466 IF (
PRESENT(my_pools))
THEN
467 rs_desc => my_rs_descs(igrid)%rs_desc
468 pw_pool => my_pools(igrid)%pool
470 rs_desc => rs_descs(igrid)%rs_desc
471 pw_pool => pw_pools(igrid)%pool
474 CALL pw_pool%create_pw(coeff_gspace)
475 CALL pw_pool%create_pw(coeff_rspace)
477 IF (igrid /= auxbas_grid)
THEN
478 pw_aux => pw_pools(auxbas_grid)%pool
479 CALL pw_aux%create_pw(coeff_gaux)
481 CALL pw_copy(coeff_gaux, coeff_gspace)
483 CALL pw_aux%give_back_pw(coeff_gaux)
484 CALL pw_aux%create_pw(coeff_raux)
485 fscale = coeff_rspace%pw_grid%dvol/coeff_raux%pw_grid%dvol
487 CALL pw_aux%give_back_pw(coeff_raux)
490 IF (coeff_gspace%pw_grid%spherical)
THEN
494 CALL pw_copy(v_rspace, coeff_rspace)
497 CALL pw_pool%give_back_pw(coeff_gspace)
504 CALL pw_pool%give_back_pw(coeff_rspace)
508 eps_rho_rspace = dft_control%qs_control%eps_rho_rspace
512 NULLIFY (hab, hab_sph, hdab, hdab_sph, pab, a_hdab, a_hdab_sph)
518 CALL reallocate(a_hdab, 1, 3, 1, 3, 1, nch_max, 1, 1)
519 CALL reallocate(a_hdab_sph, 1, 3, 1, 3, 1, nch_max)
524 grid_distributed = rs_v%desc%distributed
527 IF (
PRESENT(kforce))
THEN
531 DO ikind = 1,
SIZE(atomic_kind_set, 1)
532 NULLIFY (basis_1c_set, atom_list, harmonics, cneo_potential)
533 CALL get_atomic_kind(atomic_kind_set(ikind), atom_list=atom_list, natom=nat)
535 basis_set=basis_1c_set, basis_type=
"GAPW_1C", &
537 harmonics=harmonics, &
538 cneo_potential=cneo_potential)
540 IF (.NOT. paw_atom) cycle
542 NULLIFY (qlm_gg, lmax, npgf)
544 l0_ikind=l0_ikind, qlm_gg=qlm_gg, &
548 lmax=lmax, lmin=lmin, &
549 maxso=maxso, maxl=maxl, &
550 nset=nset, npgf=npgf)
553 ALLOCATE (intloc(nsotot, nsotot))
559 max_s_harm = harmonics%max_s_harm
560 llmax = harmonics%llmax
566 IF (
ASSOCIATED(cneo_potential))
THEN
567 NULLIFY (nuc_basis_set)
569 basis_set=nuc_basis_set, &
573 lmax=lmax_nuc, lmin=lmin_nuc, &
574 maxso=maxso_nuc, maxl=maxl_nuc, &
575 nset=nset_nuc, npgf=npgf_nuc)
576 nsotot_nuc = maxso_nuc*nset_nuc
577 ALLOCATE (intloc_nuc(nsotot_nuc, nsotot_nuc))
578 max_s_harm_nuc = cneo_potential%harmonics%max_s_harm
579 llmax_nuc = cneo_potential%harmonics%llmax
582 ALLOCATE (cg_list(2,
nsoset(max(maxl, maxl_nuc))**2, &
583 max(max_s_harm, max_s_harm_nuc)), &
584 cg_n_list(max(max_s_harm, max_s_harm_nuc)))
586 num_pe = para_env%num_pe
587 mepos = para_env%mepos
590 IF (.NOT. grid_distributed .AND. j /= mepos) cycle
592 DO iat = bo(1), bo(2)
593 iatom = atom_list(iat)
594 ra(:) =
pbc(particle_set(iatom)%r, cell)
597 CALL get_rho0_mpole(rho0_mpole=rho0_mpole, iat=iatom, qlm_tot=qlm)
609 l0_ikind, zet0, 0, 0, 0.0_dp, 0, &
610 ra, (/0.0_dp, 0.0_dp, 0.0_dp/), rs_v, &
611 hab, pab, o1=0, o2=0, &
613 calculate_forces=calculate_forces, &
614 use_virial=use_virial, my_virial_a=my_virial_a, my_virial_b=my_virial_b, &
615 hdab=hdab, a_hdab=a_hdab, use_subpatch=.true., subpatch_pattern=0)
618 DO lshell = 0, l0_ikind
619 DO is = 1,
nso(lshell)
620 iso = is +
nsoset(lshell - 1)
621 hab_sph(iso) = 0.0_dp
622 hdab_sph(1:3, iso) = 0.0_dp
623 a_hdab_sph(1:3, 1:3, iso) = 0.0_dp
624 DO ic = 1,
nco(lshell)
625 ico = ic +
ncoset(lshell - 1)
629 hab_sph(iso) = hab_sph(iso) + &
633 IF (calculate_forces)
THEN
634 hdab_sph(1:3, iso) = hdab_sph(1:3, iso) + &
642 a_hdab_sph(i, ii, iso) = a_hdab_sph(i, ii, iso) + &
645 a_hdab(i, ii, ico, 1)
659 CALL get_none0_cg_list(harmonics%my_CG, lmin(iset1), lmax(iset1), lmin(iset2), lmax(iset2), &
660 max_s_harm, llmax, cg_list, cg_n_list, max_iso_not0_local)
662 DO ipgf1 = 1, npgf(iset1)
664 DO ipgf2 = 1, npgf(iset2)
666 DO iso = 1, min(
nsoset(l0_ikind), max_iso_not0_local)
667 DO icg = 1, cg_n_list(iso)
668 iso1 = cg_list(1, icg, iso)
669 iso2 = cg_list(2, icg, iso)
671 ig1 = iso1 + n1*(ipgf1 - 1) + m1
672 ig2 = iso2 + n2*(ipgf2 - 1) + m2
674 intloc(ig1, ig2) = intloc(ig1, ig2) + qlm_gg(ig1, ig2, iso)*hab_sph(iso)
686 IF (
ASSOCIATED(cneo_potential))
THEN
689 DO iset1 = 1, nset_nuc
690 n1 =
nsoset(lmax_nuc(iset1))
692 DO iset2 = 1, nset_nuc
693 n2 =
nsoset(lmax_nuc(iset2))
695 lmax_nuc(iset1), lmin_nuc(iset2), lmax_nuc(iset2), &
696 max_s_harm_nuc, llmax_nuc, cg_list, cg_n_list, &
698 DO ipgf1 = 1, npgf_nuc(iset1)
699 DO ipgf2 = 1, npgf_nuc(iset2)
701 DO iso = 1, min(
nsoset(l0_ikind), max_iso_not0_local)
702 DO icg = 1, cg_n_list(iso)
703 iso1 = cg_list(1, icg, iso)
704 iso2 = cg_list(2, icg, iso)
706 ig1 = iso1 + n1*(ipgf1 - 1) + m1
707 ig2 = iso2 + n2*(ipgf2 - 1) + m2
709 intloc_nuc(ig1, ig2) = intloc_nuc(ig1, ig2) - cneo_potential%zeff* &
710 cneo_potential%Qlm_gg(ig1, ig2, iso)*hab_sph(iso)
723 IF (grid_distributed)
THEN
725 CALL para_env%sum(intloc)
726 IF (
ASSOCIATED(cneo_potential))
THEN
727 CALL para_env%sum(intloc_nuc)
732 rho_atom => rho_atom_set(iatom)
733 CALL get_rho_atom(rho_atom=rho_atom, ga_vlocal_gb_h=int_local_h, ga_vlocal_gb_s=int_local_s)
735 int_local_h(ispin)%r_coef = int_local_h(ispin)%r_coef + intloc
736 int_local_s(ispin)%r_coef = int_local_s(ispin)%r_coef + intloc
738 IF (
ASSOCIATED(cneo_potential))
THEN
739 rhoz_cneo => rhoz_cneo_set(iatom)
740 rhoz_cneo%ga_Vlocal_gb_h = rhoz_cneo%ga_Vlocal_gb_h + intloc_nuc
741 rhoz_cneo%ga_Vlocal_gb_s = rhoz_cneo%ga_Vlocal_gb_s + intloc_nuc
745 IF (
PRESENT(atener))
THEN
746 DO iso = 1,
nsoset(l0_ikind)
747 atener(iatom) = atener(iatom) + 0.5_dp*qlm(iso)*hab_sph(iso)
751 IF (calculate_forces)
THEN
752 force_tmp(1:3) = 0.0_dp
753 DO iso = 1,
nsoset(l0_ikind)
754 force_tmp(1) = force_tmp(1) + qlm(iso)*hdab_sph(1, iso)
755 force_tmp(2) = force_tmp(2) + qlm(iso)*hdab_sph(2, iso)
756 force_tmp(3) = force_tmp(3) + qlm(iso)*hdab_sph(3, iso)
758 force(ikind)%g0s_Vh_elec(1:3, iat) = force(ikind)%g0s_Vh_elec(1:3, iat) + fscale*force_tmp(1:3)
762 DO iso = 1,
nsoset(l0_ikind)
766 virial%pv_gapw(i, ii) = virial%pv_gapw(i, ii) + fscale*qlm(iso)*a_hdab_sph(i, ii, iso)
767 virial%pv_virial(i, ii) = virial%pv_virial(i, ii) + fscale*qlm(iso)*a_hdab_sph(i, ii, iso)
777 IF (
ASSOCIATED(intloc_nuc))
DEALLOCATE (intloc_nuc)
778 DEALLOCATE (cg_list, cg_n_list)
784 DEALLOCATE (hab, hdab, hab_sph, hdab_sph, pab, a_hdab, a_hdab_sph)
786 CALL timestop(handle)
788 IF (rho0_mpole%do_cneo)
THEN
789 IF (
PRESENT(kforce))
THEN
791 rhoz_cneo_set, kforce)
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
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.