28 USE dbcsr_api,
ONLY: dbcsr_p_type
115 #include "./base/base_uses.f90"
122 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_linres_module'
137 SUBROUTINE vcd_linres(qs_env, p_env)
138 TYPE(qs_environment_type),
POINTER :: qs_env
139 TYPE(qs_p_env_type) :: p_env
141 INTEGER :: beta, i, latom
142 LOGICAL :: mfp_is_done, mfp_repeat
143 TYPE(vcd_env_type) :: vcd_env
150 mfp_repeat = vcd_env%distributed_origin
151 mfp_is_done = .false.
153 qs_env%linres_control%linres_restart = .true.
157 DO latom = 1,
SIZE(vcd_env%dcdr_env%list_of_atoms)
158 vcd_env%dcdr_env%lambda = vcd_env%dcdr_env%list_of_atoms(latom)
165 vcd_env%dcdr_env%beta = beta
166 vcd_env%dcdr_env%deltaR(vcd_env%dcdr_env%beta, vcd_env%dcdr_env%lambda) = 1._dp
171 CALL apt_dr(qs_env, vcd_env%dcdr_env)
177 CALL apt_dv(vcd_env, qs_env)
178 CALL aat_dv(vcd_env, qs_env)
180 IF (vcd_env%do_mfp)
THEN
184 IF (.NOT. mfp_is_done .OR. mfp_repeat)
THEN
186 IF (vcd_env%origin_dependent_op_mfp)
THEN
187 cpwarn(
"Using the origin dependent MFP operator")
201 vcd_env%dcdr_env%apt_total_dcdr(:, :, vcd_env%dcdr_env%lambda) = &
202 vcd_env%dcdr_env%apt_el_dcdr(:, :, vcd_env%dcdr_env%lambda) &
203 + vcd_env%dcdr_env%apt_nuc_dcdr(:, :, vcd_env%dcdr_env%lambda)
205 vcd_env%apt_total_nvpt(:, :, vcd_env%dcdr_env%lambda) = &
206 vcd_env%apt_el_nvpt(:, :, vcd_env%dcdr_env%lambda) + vcd_env%apt_nuc_nvpt(:, :, vcd_env%dcdr_env%lambda)
208 IF (vcd_env%do_mfp) &
209 vcd_env%aat_atom_mfp(:, :, vcd_env%dcdr_env%lambda) = vcd_env%aat_atom_mfp(:, :, vcd_env%dcdr_env%lambda)*4._dp
216 END SUBROUTINE vcd_linres
228 SUBROUTINE dcdr_linres(qs_env, p_env)
229 TYPE(qs_environment_type),
POINTER :: qs_env
230 TYPE(qs_p_env_type) :: p_env
232 INTEGER :: beta, latom
233 TYPE(dcdr_env_type) :: dcdr_env
237 DO latom = 1,
SIZE(dcdr_env%list_of_atoms)
238 dcdr_env%lambda = dcdr_env%list_of_atoms(latom)
243 dcdr_env%deltaR(dcdr_env%beta, dcdr_env%lambda) = 1._dp
248 IF (.NOT. dcdr_env%localized_psi0)
THEN
249 CALL apt_dr(qs_env, dcdr_env)
250 ELSE IF (dcdr_env%localized_psi0)
THEN
256 dcdr_env%apt_total_dcdr(:, :, dcdr_env%lambda) = &
257 dcdr_env%apt_el_dcdr(:, :, dcdr_env%lambda) + dcdr_env%apt_nuc_dcdr(:, :, dcdr_env%lambda)
262 END SUBROUTINE dcdr_linres
273 TYPE(force_env_type),
POINTER :: force_env
275 CHARACTER(LEN=*),
PARAMETER :: routinen =
'linres_calculation'
278 TYPE(qs_environment_type),
POINTER :: qs_env
280 CALL timeset(routinen, handle)
284 cpassert(
ASSOCIATED(force_env))
285 cpassert(force_env%ref_count > 0)
287 SELECT CASE (force_env%in_use)
291 qs_env => force_env%qmmm_env%qs_env
293 cpabort(
"Does not recognize this force_env")
296 qs_env%linres_run = .true.
300 CALL timestop(handle)
317 TYPE(qs_environment_type),
POINTER :: qs_env
319 CHARACTER(LEN=*),
PARAMETER :: routinen =
'linres_calculation_low'
321 INTEGER :: handle, iounit
322 LOGICAL :: dcdr_present, epr_present, issc_present, &
323 lr_calculation, nmr_present, &
324 polar_present, vcd_present
325 TYPE(cp_logger_type),
POINTER :: logger
326 TYPE(dft_control_type),
POINTER :: dft_control
327 TYPE(linres_control_type),
POINTER :: linres_control
328 TYPE(qs_p_env_type) :: p_env
329 TYPE(section_vals_type),
POINTER :: lr_section, prop_section
331 CALL timeset(routinen, handle)
333 lr_calculation = .false.
334 nmr_present = .false.
335 epr_present = .false.
336 issc_present = .false.
337 polar_present = .false.
338 dcdr_present = .false.
340 NULLIFY (dft_control, linres_control, logger, prop_section, lr_section)
346 IF (lr_calculation)
THEN
347 CALL linres_init(lr_section, p_env, qs_env)
349 extension=
".linresLog")
350 CALL get_qs_env(qs_env=qs_env, dft_control=dft_control, &
351 linres_control=linres_control)
354 linres_control%property =
lr_none
362 IF (nmr_present .OR. epr_present)
THEN
363 CALL nmr_epr_linres(linres_control, qs_env, p_env, dft_control, &
364 nmr_present, epr_present, iounit)
371 IF (issc_present)
THEN
372 CALL issc_linres(linres_control, qs_env, p_env, dft_control)
378 IF (polar_present)
THEN
379 CALL polar_linres(qs_env, p_env)
386 IF (dcdr_present)
THEN
387 CALL dcdr_linres(qs_env, p_env)
394 IF (vcd_present)
THEN
395 CALL vcd_linres(qs_env, p_env)
403 WRITE (unit=iounit, fmt=
"(/,T2,A,/,T25,A,/,T2,A,/)") &
405 "ENDED LINRES CALCULATION", &
409 "PRINT%PROGRAM_RUN_INFO")
412 CALL timestop(handle)
428 SUBROUTINE linres_init(lr_section, p_env, qs_env)
430 TYPE(section_vals_type),
POINTER :: lr_section
431 TYPE(qs_p_env_type),
INTENT(OUT) :: p_env
432 TYPE(qs_environment_type),
POINTER :: qs_env
434 INTEGER :: iounit, ispin
436 TYPE(cp_logger_type),
POINTER :: logger
437 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_ks, rho_ao
438 TYPE(dft_control_type),
POINTER :: dft_control
439 TYPE(linres_control_type),
POINTER :: linres_control
440 TYPE(mo_set_type),
DIMENSION(:),
POINTER :: mos
441 TYPE(qs_rho_type),
POINTER :: rho
442 TYPE(section_vals_type),
POINTER :: loc_section
447 extension=
".linresLog")
448 NULLIFY (dft_control, linres_control, loc_section, rho, mos, matrix_ks, rho_ao)
450 ALLOCATE (linres_control)
451 CALL set_qs_env(qs_env=qs_env, linres_control=linres_control)
453 dft_control=dft_control, matrix_ks=matrix_ks, mos=mos, rho=rho)
459 l_val=linres_control%localized_psi0)
460 IF (linres_control%localized_psi0)
THEN
462 WRITE (unit=iounit, fmt=
"(/,T3,A,A)") &
463 "Localization of ground state orbitals", &
464 " before starting linear response calculation"
469 DO ispin = 1, dft_control%nspins
470 CALL calculate_density_matrix(mos(ispin), rho_ao(ispin)%matrix)
481 CALL section_vals_val_get(lr_section,
"PRECONDITIONER", i_val=linres_control%preconditioner_type)
485 WRITE (unit=iounit, fmt=
"(/,T2,A,/,T25,A,/,T2,A,/)") &
487 "START LINRES CALCULATION", &
490 WRITE (unit=iounit, fmt=
"(T2,A)") &
491 "LINRES| Properties to be calculated:"
493 IF (do_it)
WRITE (unit=iounit, fmt=
"(T62,A)")
"NMR Chemical Shift"
495 IF (do_it)
WRITE (unit=iounit, fmt=
"(T68,A)")
"EPR g Tensor"
497 IF (do_it)
WRITE (unit=iounit, fmt=
"(T43,A)")
"Indirect spin-spin coupling constants"
499 IF (do_it)
WRITE (unit=iounit, fmt=
"(T57,A)")
"Electric Polarizability"
501 IF (linres_control%localized_psi0)
WRITE (unit=iounit, fmt=
"(T2,A,T65,A)") &
502 "LINRES|",
" LOCALIZED PSI0"
504 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
505 "LINRES| Optimization algorithm",
" Conjugate Gradients"
507 SELECT CASE (linres_control%preconditioner_type)
509 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
510 "LINRES| Preconditioner",
" NONE"
512 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
513 "LINRES| Preconditioner",
" FULL_SINGLE"
515 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
516 "LINRES| Preconditioner",
" FULL_KINETIC"
518 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
519 "LINRES| Preconditioner",
" FULL_S_INVERSE"
521 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
522 "LINRES| Preconditioner",
" FULL_SINGLE_INVERSE"
524 WRITE (unit=iounit, fmt=
"(T2,A,T60,A)") &
525 "LINRES| Preconditioner",
" FULL_ALL"
527 cpabort(
"Preconditioner NYI")
530 WRITE (unit=iounit, fmt=
"(T2,A,T72,ES8.1)") &
531 "LINRES| EPS", linres_control%eps
532 WRITE (unit=iounit, fmt=
"(T2,A,T72,I8)") &
533 "LINRES| MAX_ITER", linres_control%max_iter
539 CALL p_env_create(p_env, qs_env, orthogonal_orbitals=.true., linres_control=linres_control)
544 p_env%new_preconditioner = .true.
546 "PRINT%PROGRAM_RUN_INFO")
548 END SUBROUTINE linres_init
560 SUBROUTINE nmr_epr_linres(linres_control, qs_env, p_env, dft_control, nmr_present, epr_present, iounit)
562 TYPE(linres_control_type),
POINTER :: linres_control
563 TYPE(qs_environment_type),
POINTER :: qs_env
564 TYPE(qs_p_env_type) :: p_env
565 TYPE(dft_control_type),
POINTER :: dft_control
566 LOGICAL :: nmr_present, epr_present
571 TYPE(current_env_type) :: current_env
572 TYPE(epr_env_type) :: epr_env
573 TYPE(nmr_env_type) :: nmr_env
579 IF (.NOT. linres_control%localized_psi0)
THEN
580 CALL cp_abort(__location__, &
581 "Are you sure that you want to calculate the chemical "// &
582 "shift without localized psi0?")
584 dft_control%nspins, centers_only=.true.)
586 IF (dft_control%nspins /= 2 .AND. epr_present)
THEN
587 cpabort(
"LSD is needed to perform a g tensor calculation!")
592 IF (qs_env%qmmm) do_qmmm = .true.
593 current_env%do_qmmm = do_qmmm
599 IF (current_env%all_pert_op_done)
THEN
601 IF (nmr_present)
THEN
606 IF (epr_present)
THEN
623 IF (nmr_present)
THEN
624 CALL nmr_shift(nmr_env, current_env, qs_env, ib)
628 IF (epr_present)
THEN
630 CALL epr_g_so(epr_env, current_env, qs_env, ib)
631 CALL epr_g_soo(epr_env, current_env, qs_env, ib)
636 IF (nmr_present)
THEN
642 IF (epr_present)
THEN
649 WRITE (iounit,
"(T10,A,/T20,A,/)") &
650 "CURRENT: Not all responses to perturbation operators could be calculated.", &
651 " Hence: NO nmr and NO epr possible."
657 END SUBROUTINE nmr_epr_linres
666 SUBROUTINE issc_linres(linres_control, qs_env, p_env, dft_control)
668 TYPE(linres_control_type),
POINTER :: linres_control
669 TYPE(qs_environment_type),
POINTER :: qs_env
670 TYPE(qs_p_env_type) :: p_env
671 TYPE(dft_control_type),
POINTER :: dft_control
675 TYPE(current_env_type) :: current_env
676 TYPE(issc_env_type) :: issc_env
679 IF (.NOT. linres_control%localized_psi0)
THEN
680 CALL cp_abort(__location__, &
681 "Are you sure that you want to calculate the chemical "// &
682 "shift without localized psi0?")
684 dft_control%nspins, centers_only=.true.)
689 IF (qs_env%qmmm) do_qmmm = .true.
690 current_env%do_qmmm = do_qmmm
699 DO iatom = 1, issc_env%issc_natms
709 END SUBROUTINE issc_linres
718 SUBROUTINE polar_linres(qs_env, p_env)
720 TYPE(qs_environment_type),
POINTER :: qs_env
721 TYPE(qs_p_env_type) :: p_env
729 END SUBROUTINE polar_linres
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public weber2009
integer, save, public ditler2021
integer, save, public ditler2022
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
various routines to log and control the output. The idea is that decisions about where to log should ...
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,...
Interface for the force calculations.
integer, parameter, public use_qmmm
recursive subroutine, public force_env_get(force_env, in_use, fist_env, qs_env, meta_env, fp_env, subsys, para_env, potential_energy, additional_potential, kinetic_energy, harmonic_shell, kinetic_shell, cell, sub_force_env, qmmm_env, qmmmx_env, eip_env, pwdft_env, globenv, input, force_env_section, method_name_id, root_section, mixed_env, nnp_env, embed_env)
returns various attributes about the force environment
integer, parameter, public use_qs_force
Defines the basic variable types.
integer, parameter, public dp
Calculate the derivatives of the MO coefficients wrt nuclear coordinates.
subroutine, public dcdr_env_cleanup(qs_env, dcdr_env)
Deallocate the dcdr environment.
subroutine, public dcdr_print(dcdr_env, qs_env)
Print the APT and sum rules.
subroutine, public dcdr_env_init(dcdr_env, qs_env)
Initialize the dcdr environment.
Calculate the derivatives of the MO coefficients wrt nuclear coordinates.
subroutine, public prepare_per_atom(dcdr_env, qs_env)
Prepare the environment for a choice of lambda.
subroutine, public apt_dr_localization(qs_env, dcdr_env)
Calculate atomic polar tensor using the localized dipole operator.
subroutine, public apt_dr(qs_env, dcdr_env)
Calculate atomic polar tensor.
subroutine, public dcdr_response_dr(dcdr_env, p_env, qs_env)
Get the dC/dR by solving the Sternheimer equation, using the op_dR matrix.
subroutine, public dcdr_build_op_dr(dcdr_env, qs_env)
Build the operator for the position perturbation.
collects routines that calculate density matrices
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.
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, WannierCentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Set the QUICKSTEP environment.
Chemical shift calculation by dfpt Initialization of the nmr_env, creation of the special neighbor li...
subroutine, public current_response(current_env, p_env, qs_env)
...
subroutine, public current_env_init(current_env, qs_env)
...
subroutine, public current_env_cleanup(current_env)
...
given the response wavefunctions obtained by the application of the (rxp), p, and ((dk-dl)xp) operato...
subroutine, public current_build_current(current_env, qs_env, iB)
First calculate the density matrixes, for each component of the current they are 3 because of the r d...
subroutine, public current_build_chi(current_env, qs_env, iB)
...
Calculates Nabla V_KS (local part if PSP) on the different grids.
subroutine, public epr_nablavks(epr_env, qs_env)
Evaluates Nabla V_KS on the grids.
subroutine, public epr_g_zke(epr_env, qs_env)
Calculate zke part of the g tensor.
subroutine, public epr_g_so(epr_env, current_env, qs_env, iB)
Calculates g_so.
subroutine, public epr_g_soo(epr_env, current_env, qs_env, iB)
Calculates g_soo (soft part only for now)
subroutine, public epr_ind_magnetic_field(epr_env, current_env, qs_env, iB)
...
subroutine, public epr_g_print(epr_env, qs_env)
Prints the g tensor.
g tensor calculation by dfpt Initialization of the epr_env, creation of the special neighbor lists Pe...
subroutine, public epr_env_cleanup(epr_env)
Deallocate the epr environment.
subroutine, public epr_env_init(epr_env, qs_env)
Initialize the epr environment.
Chemical shift calculation by dfpt Initialization of the issc_env, creation of the special neighbor l...
subroutine, public issc_issc(issc_env, qs_env, iatom)
...
subroutine, public issc_response(issc_env, p_env, qs_env)
Initialize the issc environment.
subroutine, public issc_env_cleanup(issc_env)
Deallocate the issc environment.
subroutine, public issc_env_init(issc_env, qs_env)
Initialize the issc environment.
subroutine, public issc_print(issc_env, qs_env)
...
localize wavefunctions linear response scf
subroutine, public linres_localize(qs_env, linres_control, nspins, centers_only)
Find the centers and spreads of the wfn, if required apply a localization algorithm.
Contains the setup for the calculation of properties by linear response by the application of second ...
subroutine, public linres_calculation(force_env)
Driver for the linear response calculatios.
subroutine, public linres_calculation_low(qs_env)
Linear response can be called as run type or as post scf calculation Initialize the perturbation envi...
from the response current density calculates the shift tensor and the susceptibility
subroutine, public nmr_shift(nmr_env, current_env, qs_env, iB)
...
subroutine, public nmr_shift_print(nmr_env, current_env, qs_env)
Shielding tensor and Chi are printed into a file if required from input It is possible to print only ...
Chemical shift calculation by dfpt Initialization of the nmr_env, creation of the special neighbor li...
subroutine, public nmr_env_cleanup(nmr_env)
Deallocate the nmr environment.
subroutine, public nmr_env_init(nmr_env, qs_env)
Initialize the nmr environment.
Calculate the operators p rxp and D needed in the optimization of the different contribution of the f...
subroutine, public current_operators(current_env, qs_env)
Calculate the first order hamiltonian applied to the ao and then apply them to the ground state orbit...
subroutine, public polar_operators(qs_env)
Calculate the dipole operator in the AO basis and its derivative wrt to MOs.
subroutine, public issc_operators(issc_env, qs_env, iatom)
...
Polarizability calculation by dfpt Initialization of the polar_env, Perturbation Hamiltonian by appli...
subroutine, public polar_polar(qs_env)
...
subroutine, public polar_print(qs_env)
Print information related to the polarisability tensor.
subroutine, public polar_env_init(qs_env)
Initialize the polar environment.
subroutine, public polar_response(p_env, qs_env)
Calculate the polarisability tensor using response theory.
Type definitiona for linear response calculations.
subroutine, public mfp_build_operator_gauge_independent(vcd_env, qs_env, alpha)
...
subroutine, public mfp_response(vcd_env, p_env, qs_env, alpha)
Get the dC/dB using the vcd_envop_dB.
subroutine, public mfp_build_operator_gauge_dependent(vcd_env, qs_env, alpha)
...
subroutine, public mfp_aat(vcd_env, qs_env)
...
Definition and initialisation of the mo data type.
Utility functions for the perturbation calculations.
subroutine, public p_env_psi0_changed(p_env, qs_env)
To be called after the value of psi0 has changed. Recalculates the quantities S_psi0 and m_epsilon.
subroutine, public p_env_create(p_env, qs_env, p1_option, p1_admm_option, orthogonal_orbitals, linres_control)
allocates and initializes the perturbation environment (no setup)
basis types for the calculation of the perturbation of density theory.
subroutine, public p_env_release(p_env)
relases the given p_env (see doc/ReferenceCounting.html)
methods of the rho structure (defined in qs_rho_types)
subroutine, public qs_rho_update_rho(rho_struct, qs_env, rho_xc_external, local_rho_set, task_list_external, task_list_external_soft, pw_env_external, para_env_external)
updates rho_r and rho_g to the rhorho_ao. if use_kinetic_energy_density also computes tau_r and tau_g...
superstucture that hold various representations of the density and keeps track of which ones are vali...
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
subroutine, public vcd_print(vcd_env, qs_env)
Print the APTs, AATs, and sum rules.
subroutine, public vcd_env_cleanup(qs_env, vcd_env)
Deallocate the vcd environment.
subroutine, public vcd_env_init(vcd_env, qs_env)
Initialize the vcd environment.
subroutine, public vcd_build_op_dv(vcd_env, qs_env)
What we are building here is the operator for the NVPT response: H0 * C1 - S0 * E0 * C1 = - op_dV lin...
subroutine, public vcd_response_dv(vcd_env, p_env, qs_env)
Get the dC/dV using the vcd_envop_dV.
subroutine, public aat_dv(vcd_env, qs_env)
Compute I_{alpha beta}^lambda = d/dV^lambda_beta <m_alpha> = d/dV^lambda_beta < r x.
subroutine, public apt_dv(vcd_env, qs_env)
Compute E_{alpha beta}^lambda = d/dV^lambda_beta <\mu_alpha> = d/dV^lambda_beta <.
subroutine, public prepare_per_atom_vcd(vcd_env, qs_env)
Initialize the matrices for the NVPT calculation.