#include <init.hpp>

Inheritance diagram for q3::UltrasoftPlanewaveBasis:

Public Member Functions
	UltrasoftPlanewaveBasis ()=default

	~UltrasoftPlanewaveBasis () override=default

void	add_augmentation_charges (const Crystal &crystal_structure, const Matrix< Real > &projector_overlap_coefficients, std::vector< Complex > &charge_density_g, const std::vector< std::vector< Complex > > *cached_structure_factors=nullptr) const

Real	augmentation_overlap_integral (int species_type, int projector_a, int projector_b) const

const AugmentationQTable &	augmentation_q_table () const

Complex	compute_q_g_single (int species_type_index, int projector_index_i, int projector_index_j, const Vector3D &g_cartesian, Real two_pi_over_lattice_parameter) const

const GauntCoefficientTable &	gaunt_table () const

bool	has_q_functions_g (int species_type_index) const

void	init (const Crystal &crystal_structure, const PlanewaveBasis &planewave_basis, const PseudopotentialSet &pseudopotential_set) override

Real	interpolate_qfunc (int species_type_index, int projector_pair_index, Real q_value) const

bool	is_ultrasoft () const override

const Complex *	q_function_g_data (int species_type_index, int channel_pair_index, int number_of_g_vectors) const

Complex	q_function_g_value (int species_type_index, int channel_pair_index, int g_index, int number_of_g_vectors) const

const SphericalHarmonicsTable &	spherical_harmonics_table () const

Public Member Functions inherited from q3::NormConservingPlanewaveBasis
	NormConservingPlanewaveBasis ()=default

virtual	~NormConservingPlanewaveBasis ()=default

const std::vector< Real > &	atomic_charge_density_g (int species_type) const

Real	bare_d_matrix_element (int species_type, int projector_a, int projector_b) const

int	beta_angular_momentum (int species_type, int projector_index) const

const BetaFunctionTable &	beta_function_table () const

const std::vector< Complex > &	core_charge_density_reciprocal () const

const CoreChargeDensityTable &	core_charge_density_table () const

bool	has_nonlinear_core_correction () const

Real	interpolate_beta (int species_type_index, int beta_index, Real q_value) const

Real	interpolate_core_charge_density (int element_type, Real q_value) const

Real	interpolate_core_charge_density_derivative (int element_type, Real q_value) const

Real	interpolate_dbeta (int species_type_index, int beta_index, Real q_value) const

Real	interpolate_local_potential (int element_type, Real q_value) const

Real	interpolate_local_potential_derivative (int element_type, Real q_value) const

Real	interpolate_rho_core (int element_type, Real q_value) const

const std::vector< Real > &	local_potential_g (int species_type_index) const

const LocalPotentialTable &	local_potential_table () const

int	number_of_projector_channels (int species_type) const

int	number_of_projectors (int species_type) const

const PlanewaveBasis &	planewave_basis () const

int	projector_angular_momentum (int species_type, int projector_index) const

int	projector_combined_lm_index (int species_type, int projector_index) const

int	projector_magnetic_quantum_number (int species_type, int projector_index) const

const ProjectorMappingTable &	projector_mapping_table () const

int	projector_to_beta_index (int species_type, int projector_index) const

const PseudopotentialSet &	pseudopotential_set () const

const std::vector< Complex > &	total_local_potential_g () const

int	total_number_of_projectors () const

Protected Member Functions
void	init_ultrasoft (const Crystal &crystal_structure, const PlanewaveBasis &planewave_basis)

void	precompute_q_functions_g (const PlanewaveBasis &planewave_basis)

Protected Member Functions inherited from q3::NormConservingPlanewaveBasis
void	compute_mesh_truncation_sizes (int number_of_species_types)

Static Protected Member Functions
static int	projector_degree_for_gaunt_table (const PseudopotentialSet &pseudopotential_set, int maximum_angular_momentum_q)

Additional Inherited Members
Protected Attributes inherited from q3::NormConservingPlanewaveBasis
std::vector< std::vector< Real > >	atomic_charge_density_g_

std::unique_ptr< BetaFunctionTable >	beta_function_table_

Real	cell_volume_ {1.0}

CoreChargeData	core_charge_data_

CoreChargeDensityTable	core_charge_density_table_

std::vector< std::vector< Real > >	local_potential_g_per_species_

LocalPotentialTable	local_potential_table_

std::vector< int >	mesh_sizes_

const PlanewaveBasis *	planewave_basis_ptr_ = nullptr

ProjectorMappingTable	projector_mapping_table_

PseudopotentialSet	pseudopotential_set_

std::vector< Complex >	total_local_potential_g_

Detailed Description

Definition at line 41 of file init.hpp.

Constructor & Destructor Documentation

◆ UltrasoftPlanewaveBasis()

q3::UltrasoftPlanewaveBasis::UltrasoftPlanewaveBasis ( )

default

◆ ~UltrasoftPlanewaveBasis()

q3::UltrasoftPlanewaveBasis::~UltrasoftPlanewaveBasis ( )

overridedefault

Member Function Documentation

◆ add_augmentation_charges()

void q3::UltrasoftPlanewaveBasis::add_augmentation_charges	(	const Crystal &	crystal_structure,
		const Matrix< Real > &	projector_overlap_coefficients,
		std::vector< Complex > &	charge_density_g,
		const std::vector< std::vector< Complex > > *	cached_structure_factors = `nullptr`
	)		const

inline

Add ultrasoft augmentation charges to the charge density in G-space.

Uses the pre-computed Q_ij(G) table for O(atoms × pairs × G) cost, eliminating the per-G interpolation that dominated the previous implementation.

Parameters

crystal_structure	Crystal geometry
projector_overlap_coefficients	⟨ψ\|β⟩⟨β\|ψ⟩ summed over bands
charge_density_g	Charge density in G-space (modified in place)

Skip species that have no Q(G) table, exactly as the augmentation force kernel in forces/.../ultrasoft/internal/augmentation.hpp already does — the two assemble the same A_I(G) and must agree on which atoms contribute to it.

In a mixed cell precompute_q_functions_g leaves the entry for a norm-conserving species empty, but ρ_ij is accumulated over every atom, so its ⟨ψ|β_i⟩⟨β_j|ψ⟩ are nonzero and this loop reached q_function_g_data for it anyway. That takes the address of an element of an empty vector and then reads number_of_g_vectors complex values through it — out of bounds on every SCF iteration, silently corrupting the density rather than failing outright. A norm-conserving species has no augmentation charge to add, so contributing nothing is also the right answer.

Accumulate ρ_ij-weighted Q_ij(G) over all channel pairs. Exploit Q_ij(G) = Q_ji(G) symmetry: loop only over upper triangle. For i < j: use (ρ_ij + ρ_ji) * Q_ij(G) in a single pass. For i == j: use ρ_ii * Q_ii(G) directly.

Definition at line 178 of file init.hpp.

References q3::compute_phase_factor(), q3::Crystal::convert_from_crystal_to_cartesian(), q3::PlanewaveBasis::g_vectors(), has_q_functions_g(), q3::Crystal::number_of_atoms(), q3::ProjectorMappingTable::number_of_channels(), q3::PlanewaveBasis::number_of_g_vectors(), q3::NormConservingPlanewaveBasis::number_of_projector_channels(), q3::NormConservingPlanewaveBasis::planewave_basis_ptr_, q3::NormConservingPlanewaveBasis::projector_mapping_table_, q_function_g_data(), q3::Crystal::species_type_index(), and q3::constants::two_pi.

◆ augmentation_overlap_integral()

Real q3::UltrasoftPlanewaveBasis::augmentation_overlap_integral	(	int	species_type,
		int	projector_a,
		int	projector_b
	)		const

inline

Augmentation overlap integral q_{ij} = ∫ Q_ij(r) dr, or zero when the species carries no augmentation charges.

The gate belongs here, not in each caller. This class is instantiated whenever any species needs augmentation, so a mixed cell — one ultrasoft or PAW species alongside one norm-conserving species — brings its norm-conserving species through every method on it. as_ultrasoft is a dynamic_cast whose null check is an assert, so asking a norm-conserving species for its q_ij aborts a Debug build and dereferences null in the Release builds users actually run.

q_ij = 0 is the physically correct answer, not merely a safe one: a norm-conserving dataset has no augmentation charges, so it contributes nothing to S = 1 + Σ q_ij |β_i⟩⟨β_j|. Its projectors still participate in the nonlocal potential; only the augmentation term vanishes.

Callers used to carry this test themselves and mostly did not — the whole-system has_ultrasoft() flag reads as if it settled the question per species, and it does not. hubbard_augmentation_overlap and the force kernel in forces/.../nonlocal.hpp grew their own copies for this exact reason; a mixed US+NC cell still aborted in the wavefunction initializer's overlap, which had none.

Parameters

species_type	Species index
projector_a	Beta projector index i (radial)
projector_b	Beta projector index j (radial)

Definition at line 123 of file init.hpp.

References q3::PseudopotentialSet::as_ultrasoft(), q3::UltrasoftPseudopotential::augmentation_overlap_integrals, and q3::NormConservingPlanewaveBasis::pseudopotential_set_.

Referenced by q3::hubbard_augmentation_overlap().

◆ augmentation_q_table()

const AugmentationQTable & q3::UltrasoftPlanewaveBasis::augmentation_q_table ( ) const

inline

Access the augmentation Q-function table.

Definition at line 62 of file init.hpp.

◆ compute_q_g_single()

Complex q3::UltrasoftPlanewaveBasis::compute_q_g_single	(	int	species_type_index,
		int	projector_index_i,
		int	projector_index_j,
		const Vector3D &	g_cartesian,
		Real	two_pi_over_lattice_parameter
	)		const

inline

Compute the augmentation Q-function Q_ij(G) for a single G-vector.

Parameters

species_type_index	Species index
projector_index_i	Projector channel index i
projector_index_j	Projector channel index j
g_cartesian	G-vector in Cartesian reciprocal units (2π/a)
two_pi_over_lattice_parameter	Conversion factor 2π/alat

Returns: Complex Q_ij(G) value

Same per-species gate as augmentation_overlap_integral: Q_ij(G) is identically zero for a species with no augmentation charges, and reaching for its (nonexistent) table is the Release null dereference.

Definition at line 288 of file init.hpp.

References q3::ProjectorMappingTable::angular_momentum(), q3::PseudopotentialSet::as_ultrasoft(), q3::UltrasoftPseudopotential::augmentation_overlap_integrals, q3::GauntCoefficientTable::coefficients, q3::ProjectorMappingTable::combined_lm_index(), q3::compute_real_spherical_harmonic(), q3::AugmentationQTable::interpolate_augmentation_function(), q3::AugmentationQTable::maximum_angular_momentum_q(), q3::GauntCoefficientTable::nonzero_product_indices, q3::Vector3D::norm(), q3::GauntCoefficientTable::number_of_nonzero_product_indices, q3::NormConservingPlanewaveBasis::number_of_projectors(), q3::NormConservingPlanewaveBasis::projector_mapping_table_, q3::NormConservingPlanewaveBasis::pseudopotential_set_, q3::ProjectorMappingTable::to_beta_index(), and q3::math::tolerances::zero_vector_magnitude_threshold.

Referenced by q3::PlanewaveBlochOverlapCalculator::add_augmentation_contribution().

◆ gaunt_table()

const GauntCoefficientTable & q3::UltrasoftPlanewaveBasis::gaunt_table ( ) const

inline

Access the Gaunt coefficient table.

Definition at line 76 of file init.hpp.

◆ has_q_functions_g()

bool q3::UltrasoftPlanewaveBasis::has_q_functions_g ( int species_type_index ) const

inline

Check whether Q(G) table is available for a species.

Definition at line 163 of file init.hpp.

Referenced by add_augmentation_charges().

◆ init()

void q3::UltrasoftPlanewaveBasis::init	(	const Crystal &	crystal_structure,
		const PlanewaveBasis &	planewave_basis,
		const PseudopotentialSet &	pseudopotential_set
	)

inlineoverridevirtual

Initialize all pseudopotential data from crystal, basis, and PP files. Must be called after PlanewaveBasis::generate().

Step 2: Tabulate beta functions on a uniform q-grid, with the table implementation chosen by what the datasets carry (see the factory).

Reimplemented from q3::NormConservingPlanewaveBasis.

Definition at line 46 of file init.hpp.

References q3::NormConservingPlanewaveBasis::init(), init_ultrasoft(), q3::NormConservingPlanewaveBasis::planewave_basis(), and q3::NormConservingPlanewaveBasis::pseudopotential_set().

Referenced by q3::UltrasoftPlanewaveKohnShamQuantumSystem::init().

◆ init_ultrasoft()

void q3::UltrasoftPlanewaveBasis::init_ultrasoft	(	const Crystal &	crystal_structure,
		const PlanewaveBasis &	planewave_basis
	)

inlineprotected

Definition at line 391 of file init.hpp.

References q3::NormConservingPlanewaveBasis::beta_function_table_, q3::Crystal::cell_volume(), q3::NormConservingPlanewaveBasis::cell_volume_, q3::generate_gaunt_coefficients(), q3::log::info(), q3::SphericalHarmonicsTable::init(), q3::AugmentationQTable::init(), q3::AugmentationQTable::maximum_angular_momentum_q(), q3::NormConservingPlanewaveBasis::mesh_sizes_, q3::NormConservingPlanewaveBasis::planewave_basis(), precompute_q_functions_g(), projector_degree_for_gaunt_table(), and q3::NormConservingPlanewaveBasis::pseudopotential_set_.

Referenced by init().

◆ interpolate_qfunc()

Real q3::UltrasoftPlanewaveBasis::interpolate_qfunc	(	int	species_type_index,
		int	projector_pair_index,
		Real	q_value
	)		const

inline

Interpolate augmentation function Q_{ij}(|q|) from the radial table.

Parameters

species_type_index	Species index
projector_pair_index	Linearized (i,j) pair index for the Q function
q_value	Wavevector magnitude in atomic units

Definition at line 89 of file init.hpp.

References q3::AugmentationQTable::interpolate_qfunc().

◆ is_ultrasoft()

bool q3::UltrasoftPlanewaveBasis::is_ultrasoft ( ) const

inlineoverridevirtual

Whether this basis uses ultrasoft pseudopotentials. Base class returns false; UltrasoftPlanewaveBasis returns true.

Reimplemented from q3::NormConservingPlanewaveBasis.

Definition at line 53 of file init.hpp.

◆ precompute_q_functions_g()

void q3::UltrasoftPlanewaveBasis::precompute_q_functions_g ( const PlanewaveBasis & planewave_basis )

inlineprotected

Pre-compute Q_ij(G) for all G-vectors and all channel pairs per species.

For each species s with N_ch channels and each pair (i,j): Q_ij(G) = Σ_{lm} gaunt(lm, lm_i, lm_j) * Q_ij^l(s, |G|)

Ylm(Ĝ) * (-i)^l

where Q_ij^l(s, |G|) is the radial Fourier transform of species s's augmentation function Q_ij in channel l (interpolate_augmentation_function).

Stored as q_functions_g_[species][pair_linear_index][g_index] where pair_linear_index = channel_i * n_channels + channel_j.

Memory: n_species × n_channels² × n_g_vectors × sizeof(Complex)

Capability, not type — a PAW species has Q_ij^L too (synthesized for datasets that ship only the augmentation charges) and needs the same Q(G) table.

Definition at line 477 of file init.hpp.

References q3::ProjectorMappingTable::angular_momentum(), q3::GauntCoefficientTable::coefficients, q3::ProjectorMappingTable::combined_lm_index(), q3::PlanewaveBasis::g_vector_magnitudes(), q3::log::info(), q3::AugmentationQTable::interpolate_augmentation_function(), q3::AugmentationQTable::maximum_angular_momentum_q(), q3::GauntCoefficientTable::nonzero_product_indices, q3::ProjectorMappingTable::number_of_channels(), q3::PlanewaveBasis::number_of_g_vectors(), q3::GauntCoefficientTable::number_of_nonzero_product_indices, q3::NormConservingPlanewaveBasis::number_of_projector_channels(), q3::NormConservingPlanewaveBasis::planewave_basis(), q3::NormConservingPlanewaveBasis::projector_mapping_table_, q3::NormConservingPlanewaveBasis::pseudopotential_set_, q3::PseudopotentialSet::size(), q3::ProjectorMappingTable::to_beta_index(), and q3::SphericalHarmonicsTable::value().

Referenced by init_ultrasoft().

◆ projector_degree_for_gaunt_table()

static int q3::UltrasoftPlanewaveBasis::projector_degree_for_gaunt_table	(	const PseudopotentialSet &	pseudopotential_set,
		int	maximum_angular_momentum_q
	)

inlinestaticprotected

The projector degree the Gaunt table is to be built for, checked against the product degree its consumers will index.

generate_gaunt_coefficients expands products of projector harmonics, so its argument is the projector degree and it evaluates harmonics up to twice that. What used to be passed instead was maximum_angular_momentum_q(), which is the maximum L of the Q expansion and is already twice the projector degree – so the table was sized for a product degree of 8 on any d-channel dataset, its collocation asked for harmonics at degrees 5 through 8, and compute_real_spherical_harmonic answered every one of them with 0.0.

Nothing observed it: the damaged entries are the ones no consumer indexes. Consumers reach product degree lmaxq and no further, and every entry at or below that sum was correct. The table was simply about nine times larger than anything reads, with the excess wrong.

Passing the projector degree makes the collocation exact over exactly the range consumers use, and – because it no longer doubles an already doubled degree – lets the harmonic refusal above degree 4 be added without rejecting every ultrasoft and PAW run.

The relation lmaxq <= 2 * projector_degree is what makes the smaller table still cover every index a consumer forms. It holds for every dataset Q3 reads, because lmaxq is derived from the projector channels; a dataset that declared otherwise would index past the end of the table, so it is checked rather than assumed.

Definition at line 442 of file init.hpp.

References q3::NormConservingPlanewaveBasis::pseudopotential_set(), and q3::PseudopotentialSet::size().

Referenced by init_ultrasoft().

◆ q_function_g_data()

const Complex * q3::UltrasoftPlanewaveBasis::q_function_g_data	(	int	species_type_index,
		int	channel_pair_index,
		int	number_of_g_vectors
	)		const

inline

Access the raw Q(G) data pointer for a given species and channel pair.

Returns a pointer to number_of_g_vectors contiguous Complex values. Used for BLAS-accelerated dot products in newd and addusdens.

Definition at line 154 of file init.hpp.

Referenced by add_augmentation_charges(), and q3::UltrasoftPlanewaveEffectivePotential::compute().

◆ q_function_g_value()

Complex q3::UltrasoftPlanewaveBasis::q_function_g_value	(	int	species_type_index,
		int	channel_pair_index,
		int	g_index,
		int	number_of_g_vectors
	)		const

inline

Access a pre-computed Q_ij(G) value.

Parameters

species_type_index	Species index
channel_pair_index	Linearized channel pair: channel_i * n_channels + channel_j
g_index	G-vector index
number_of_g_vectors	Total number of G-vectors (stride for flat indexing)

Returns: The pre-computed Q_ij(G) complex value

Definition at line 140 of file init.hpp.

◆ spherical_harmonics_table()

const SphericalHarmonicsTable & q3::UltrasoftPlanewaveBasis::spherical_harmonics_table ( ) const

inline

Access the precomputed spherical harmonics table.

Definition at line 69 of file init.hpp.

The documentation for this class was generated from the following file:

src/layer_3_representation/basis/planewave/pseudopotential/ultrasoft/init.hpp

Public Member Functions

Protected Member Functions

Static Protected Member Functions

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ UltrasoftPlanewaveBasis()

◆ ~UltrasoftPlanewaveBasis()

Member Function Documentation

◆ add_augmentation_charges()

◆ augmentation_overlap_integral()

◆ augmentation_q_table()

◆ compute_q_g_single()

◆ gaunt_table()

◆ has_q_functions_g()

◆ init()

◆ init_ultrasoft()

◆ interpolate_qfunc()

◆ is_ultrasoft()

◆ precompute_q_functions_g()

◆ projector_degree_for_gaunt_table()

◆ q_function_g_data()

◆ q_function_g_value()

◆ spherical_harmonics_table()