#include <init.hpp>
Inheritance diagram for q3::UltrasoftPlanewaveEffectivePotential:Public Member Functions | |
| UltrasoftPlanewaveEffectivePotential ()=default | |
| void | add_per_atom_d_terms (const std::vector< std::vector< Real > > &terms_by_atom) |
| void | add_per_atom_d_terms_for_spin (int spin, const std::vector< std::vector< Real > > &terms) |
| const std::vector< std::vector< Complex > > & | cached_structure_factors () const |
| void | compute (const Crystal &crystal_structure, const PlanewaveBasis &planewave_basis, const UltrasoftPlanewaveBasis &ultrasoft_pseudopotential_basis, const std::vector< Complex > &total_potential_g) |
| void | compute_by_spin (const Crystal &crystal_structure, const PlanewaveBasis &planewave_basis, const UltrasoftPlanewaveBasis &ultrasoft_pseudopotential_basis, const std::vector< std::vector< Complex > > &total_potential_g_by_spin) |
| Real | element (int atom_index, int projector_index_i, int projector_index_j, int spin=0) const |
| bool | is_initialized () const |
| void | precompute_structure_factors (const Crystal &crystal_structure, const PlanewaveBasis &planewave_basis) |
| Public Member Functions inherited from q3::NormConservingPlanewaveEffectivePotential | |
| virtual | ~NormConservingPlanewaveEffectivePotential () override=default |
| PlanewaveEffectivePotential ()=default | |
| PlanewaveEffectivePotential (const PlanewaveBasis &basis, int number_of_spin_components) | |
| Public Member Functions inherited from q3::PlanewaveEffectivePotential | |
| PlanewaveEffectivePotential ()=default | |
| PlanewaveEffectivePotential (const PlanewaveBasis &basis, int number_of_spin_components) | |
| const PlanewaveBasis & | basis () const |
| std::vector< Complex > & | charge_density_reciprocal () |
| const std::vector< Complex > & | charge_density_reciprocal () const |
| std::vector< Complex > & | exchange_correlation_reciprocal_space_data () |
| const std::vector< Complex > & | exchange_correlation_reciprocal_space_data () const |
| int | grid_size () const |
| void | initialize_kinetic_potential () |
| std::vector< std::vector< Real > > & | kinetic_potential_data () |
| const std::vector< std::vector< Real > > & | kinetic_potential_data () const |
| void | log_summary () const |
| std::vector< Real > & | real_space_by_spin_index (int spin_index) |
| const std::vector< Real > & | real_space_by_spin_index (int spin_index) const |
| void | to_real_space () override |
| void | to_reciprocal_space () override |
| Public Member Functions inherited from q3::EffectivePotential | |
| virtual | ~EffectivePotential ()=default |
| Public Member Functions inherited from q3::SpinResolvedEntity | |
| virtual | ~SpinResolvedEntity ()=default |
| bool | is_spin_polarized () const |
| int | number_of_spin_components () const |
| Public Member Functions inherited from q3::FFTransformable< SpinResolvedPotentialGrid, PotentialReciprocalCoefficients > | |
| virtual | ~FFTransformable ()=default |
| SpinResolvedPotentialGrid & | real_space_representation () |
| const SpinResolvedPotentialGrid & | real_space_representation () const |
| PotentialReciprocalCoefficients & | reciprocal_space_representation () |
| const PotentialReciprocalCoefficients & | reciprocal_space_representation () const |
Protected Attributes | |
| std::vector< std::vector< Real > > | effective_d_matrix_ |
| std::vector< std::vector< std::vector< Real > > > | effective_d_matrix_by_spin_ |
One D per spin channel; empty when the spin-free compute was used. | |
| Protected Attributes inherited from q3::SpinResolvedEntity | |
| int | number_of_spin_components_ = 1 |
| Protected Attributes inherited from q3::FFTransformable< SpinResolvedPotentialGrid, PotentialReciprocalCoefficients > | |
| SpinResolvedPotentialGrid | real_space_data_ |
| PotentialReciprocalCoefficients | reciprocal_space_data_ |
Additional Inherited Members | |
| Public Attributes inherited from q3::PlanewaveEffectivePotential | |
| Real | exchange_correlation_energy = 0.0 |
| Real | exchange_correlation_potential_energy = 0.0 |
| Real | hartree_energy = 0.0 |
|
default |
|
inline |
Add per-atom terms that a layered flavor contributes to D.
PAW's D is this D plus a one-center term:
D_ij = D⁰_ij + ∫ v_eff Q_ij + (D¹ − D̃¹)_Hxc
└─ dataset ─┘ └─ here ─┘ └─ passed in ─┘
The first two are what compute() built and are never recomputed for PAW — the kinetic and core-reference pieces baked into PP_DIJ/DION at generation time would be double-counted if they were. So the one-center term is a decoration of this result, and this is the seam it arrives through: the correction depends on the on-site occupancies, which only the SCF loop knows, so it is passed in rather than computed here.
Expressed as a method rather than a subclass override deliberately. The iterators hold this type by value, so a PawPlanewaveEffectivePotential could not actually be substituted for it without pointer-ifying every holder — a subclass that cannot be substituted is worse than a named operation on the object that owns D.
Call after compute(), once per SCF step. An empty entry means "no terms for this atom", so a mixed PAW/ultrasoft cell needs no special case.
A size mismatch means the caller's channel expansion disagrees with the projector mapping table's, which would otherwise write plausible numbers into the wrong elements. There is no partial-overlap reading of these matrices that is meaningful, so refusing is the only safe response.
Definition at line 347 of file init.hpp.
References effective_d_matrix_.
Referenced by add_per_atom_d_terms_for_spin().
|
inline |
Add per-atom terms to one spin channel's D (PAW's one-center term).
Refuses a channel that was never computed rather than returning quietly. A silent no-op would leave that channel's Hamiltonian carrying the other channel's one-center D — via element()'s fallback — while the double-counting term still subtracted its own, so the SCF would converge to a stationary point of something that is not the energy it prints.
Definition at line 88 of file init.hpp.
References add_per_atom_d_terms(), effective_d_matrix_, and effective_d_matrix_by_spin_.
|
inline |
Return the cached structure factors (computed on first call to compute()). Useful for sharing with augmentation charge code.
Definition at line 431 of file init.hpp.
Referenced by q3::UltrasoftPlanewavePseudopotentialKohnShamSCFIterator::add_augmentation_to_charge_density().
|
inline |
Phase 1: For each US species, precompute aux_ij(G) = conj(V(G)) * Q_ij(G) This is done once per species, not once per atom.
Correct formula: D_ij = ∫ Q(r−R)·V(r) dr = Σ_G Q(G)·conj(V(G))·e^{−iGR}·Ω [For real V: V(−G)=V*(G), so the integral picks up V*(G).] Taking conj(Q) instead of conj(V) gives the same result only when Im(V(G)) = 0 — true for symmetric structures but NOT for displaced atoms.
A species contributes ∫v_eff·Q_ij here exactly when it has augmentation charges — which PAW does. Skipping on != ultrasoft would leave a PAW atom's D at D_ion alone, a large and entirely silent error.
D_ij = ∫ Q(r−R)·V(r) dr = Σ_G Q(G)·V*(G)·e^{−iGR}·Ω Must conjugate V, NOT Q (Q may also be complex due to (−i)^l factor).
Phase 2: Compute D_eff for all atoms using batched BLAS zgemm per species. Instead of 54 individual zgemv calls, we group atoms by species and compute all integrals in a single zgemm: result = aux^T × sf_matrix. aux is (n_g × n_pairs), sf_matrix is (n_g × n_atoms_of_species), result is (n_pairs × n_atoms_of_species).
zgemm: result = aux^T × sf_matrix aux is (n_g × n_pairs) column-major, sf_matrix is (n_g × n_atoms) column-major result = aux^T × sf = (n_pairs × n_g) × (n_g × n_atoms) = (n_pairs × n_atoms)
Unpack results into per-atom D-matrices For gamma-only: multiply by 2 (half-sphere → full sphere), then subtract the G=0 contribution once (it was double-counted).
If gamma-only, precompute the G=0 correction: G=0 term = Re(conj(V(0))*Q_ij(0)*sf(0)) = Re(conj(V(0))*Q_ij(0)*1) sf(G=0) = exp(i·0·τ) = 1 for all atoms. At G=0, V(0) is real and Q(0) is real (l=0 only), so conj(V(0))=V(0). This correction is the same for all atoms of this species.
Definition at line 100 of file init.hpp.
References q3::Crystal::cell_volume(), effective_d_matrix_, q3::PlanewaveBasis::first_nonzero_g_index(), q3::PlanewaveBasis::is_gamma_only(), q3::Crystal::number_of_atoms(), q3::PlanewaveBasis::number_of_g_vectors(), precompute_structure_factors(), q3::NormConservingPlanewaveBasis::projector_mapping_table(), q3::NormConservingPlanewaveBasis::pseudopotential_set(), q3::UltrasoftPlanewaveBasis::q_function_g_data(), q3::Crystal::species_type_index(), and q3::blas::zgemm_().
Referenced by compute_by_spin().
|
inline |
Compute the effective D-matrix for all atoms from the current V_eff(G).
Uses the pre-computed Q_ij(G) table from the ultrasoft basis to avoid per-G interpolation. The computation reduces to: D_eff(i,j) = D_ion(i,j) + Re(Σ_G conj(V(G)) * Q_ij(G) * e^{-iG·τ}) * Ω
| crystal_structure | Crystal geometry |
| planewave_basis | Plane-wave basis (G-vectors) |
| ultrasoft_pseudopotential_basis | Ultrasoft basis (pre-computed Q(G), projector mappings) |
| total_potential_g | Total effective potential in G-space Compute D for every spin channel. |
D_ij = D_ion + ∫v_eff·Q_ij, and v_eff is spin-dependent, so D is too. Calling this with one potential reproduces the unpolarized result exactly; calling it with two gives each channel its own D, which is what a magnetic system needs and what a shared, spin-averaged D silently denies it.
Leave effective_d_matrix_ holding channel 0 so anything still reading the spin-free accessor sees a real D rather than whatever the last loop iteration left behind.
Definition at line 57 of file init.hpp.
References compute(), effective_d_matrix_, and effective_d_matrix_by_spin_.
Referenced by q3::UltrasoftPlanewavePseudopotentialKohnShamSCFIterator::prepare_hamiltonian().
|
inline |
Read an element of the effective D-matrix.
Falls back to the spin-free matrix when no per-spin one was computed, which is every unpolarized run and every path that still calls the single-potential compute.
Definition at line 377 of file init.hpp.
References effective_d_matrix_, and effective_d_matrix_by_spin_.
Referenced by q3::UltrasoftPlanewaveKohnShamQuantumSystem::d_matrix_element().
|
inline |
Definition at line 402 of file init.hpp.
References effective_d_matrix_.
Referenced by q3::UltrasoftPlanewaveKohnShamQuantumSystem::d_matrix_element().
|
inline |
Precompute and cache structure factors for all atoms. These are constant during SCF since atomic positions don't change.
Definition at line 410 of file init.hpp.
References q3::compute_phase_factor(), q3::Crystal::convert_from_crystal_to_cartesian(), q3::PlanewaveBasis::g_vectors(), q3::Crystal::number_of_atoms(), q3::PlanewaveBasis::number_of_g_vectors(), and q3::constants::two_pi.
Referenced by compute().
|
protected |
Effective D-matrix: per-atom × (projector_channels²).
Protected rather than private because PAW's D is this D plus a one-center term — the seam named in the A3 addendum. A subclass adds to it in place; nothing else may.
Definition at line 443 of file init.hpp.
Referenced by add_per_atom_d_terms(), add_per_atom_d_terms_for_spin(), compute(), compute_by_spin(), element(), and is_initialized().
|
protected |
One D per spin channel; empty when the spin-free compute was used.
Definition at line 446 of file init.hpp.
Referenced by add_per_atom_d_terms_for_spin(), compute_by_spin(), and element().
1.9.8