(git:ccc2433)
pao_linpot_full.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief Full parametrization of Fock matrix, ie. the identity parametrization.
10 !> \author Ole Schuett
11 ! **************************************************************************************************
13  USE basis_set_types, ONLY: gto_basis_set_type
14  USE kinds, ONLY: dp
15  USE qs_environment_types, ONLY: get_qs_env,&
16  qs_environment_type
17  USE qs_kind_types, ONLY: get_qs_kind,&
18  qs_kind_type
19 #include "./base/base_uses.f90"
20 
21  IMPLICIT NONE
22 
23  PRIVATE
24 
25  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pao_linpot_full'
26 
28 
29 CONTAINS
30 
31 ! **************************************************************************************************
32 !> \brief Count number of terms for given atomic kind
33 !> \param qs_env ...
34 !> \param ikind ...
35 !> \param nterms ...
36 ! **************************************************************************************************
37  SUBROUTINE linpot_full_count_terms(qs_env, ikind, nterms)
38  TYPE(qs_environment_type), POINTER :: qs_env
39  INTEGER, INTENT(IN) :: ikind
40  INTEGER, INTENT(OUT) :: nterms
41 
42  INTEGER :: n
43  TYPE(gto_basis_set_type), POINTER :: basis_set
44  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
45 
46  CALL get_qs_env(qs_env, qs_kind_set=qs_kind_set)
47  CALL get_qs_kind(qs_kind_set(ikind), basis_set=basis_set)
48  n = basis_set%nsgf
49  nterms = n + n*(n - 1)/2
50 
51  END SUBROUTINE linpot_full_count_terms
52 
53 ! **************************************************************************************************
54 !> \brief Builds potential terms
55 !> \param V_blocks ...
56 ! **************************************************************************************************
57  SUBROUTINE linpot_full_calc_terms(V_blocks)
58  REAL(dp), DIMENSION(:, :, :), INTENT(OUT) :: v_blocks
59 
60  INTEGER :: i, j, kterm, n, nterms
61 
62  n = SIZE(v_blocks, 1)
63  cpassert(SIZE(v_blocks, 2) == n)
64  nterms = SIZE(v_blocks, 3)
65 
66  v_blocks = 0.0_dp
67  kterm = 0
68  DO i = 1, n
69  DO j = i, n
70  kterm = kterm + 1
71  v_blocks(i, j, kterm) = 1.0_dp
72  v_blocks(j, i, kterm) = 1.0_dp
73  END DO
74  END DO
75 
76  cpassert(kterm == nterms)
77  END SUBROUTINE linpot_full_calc_terms
78 
79 END MODULE pao_linpot_full
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
Full parametrization of Fock matrix, ie. the identity parametrization.
subroutine, public linpot_full_calc_terms(V_blocks)
Builds potential terms.
subroutine, public linpot_full_count_terms(qs_env, ikind, nterms)
Count number of terms for given atomic kind.
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_nonbond, sab_almo, sab_kp, sab_kp_nosym, 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, ecoul_1c, rho0_s_rs, rho0_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, 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, rhs)
Get the QUICKSTEP environment.
Define the quickstep kind type and their sub types.
Definition: qs_kind_types.F:23
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.