Introduction to GPAW internals
Warning
This page describes “New GPAW”!
Full picture
The ase.Atoms object has an
gpaw.new.ase_interface.ASECalculator object attached
created with the gpaw.dft.GPAW() function:
>>> from ase import Atoms
>>> from gpaw.dft import GPAW, PW
>>> atoms = Atoms('H2',
... positions=[(0, 0, 0), (0, 0, 0.75)],
... cell=[2, 2, 3],
... pbc=True)
>>> atoms.calc = GPAW(mode=PW(ecut=400.0), txt='h2.txt')
>>> atoms.calc
ASECalculator(mode=PW(ecut=400.0))
The atoms.calc object manages a
gpaw.new.calculation.DFTCalculation object that does the actual work.
When we do this:
>>> e = atoms.get_potential_energy()
the gpaw.new.ase_interface.ASECalculator.get_potential_energy()
method gets called (atoms.calc.get_potential_energy(atoms))
and the following will happen:
create
gpaw.new.calculation.DFTCalculationobject if not already doneupdate positions/unit cell if they have changed
start SCF loop and converge if needed
calculate energy
store a copy of the atoms
DFT-calculation object
An instance of the gpaw.new.calculation.DFTCalculation class has
the following attributes:
|
|
|
|
|
|
|
|
|
Overview:
atoms.calc:dft:density:nt_sR,D_asii, …potential:vt_sR,dH_asii, …pot_calc:xcpoisson_solver
scf_loop:eigensolverhamiltonianocc_calcmixer
ibzwfs:ibz:symmetriesbz
_wfs_u[u]:psit_nXocc_n…
See also: code.svg.
There are three ways to create a DFTCalculation
object:
Via the the
gpaw.dft.GPAW()function which will create angpaw.new.ase_interface.ASECalculatorobject that has adftattribute:atoms.calc = GPAW(<parameters>) atoms.get_potential_energy() dft = atoms.calc.dft
Directly using the
DFTCalculationconstructor (not recommended):dft = DFTCalculation(...)
Using the
gpaw.dft.DFT()function:dft = DFT(atoms, <parameters>)
Via the
Parametersobject:dft = Parameters(<parameters>).dft_calculation(atoms)
The
Parametersis used by both thegpaw.dft.DFT()andgpaw.dft.GPAW()functions to handle:error checking
normalization
backwards compatibility and deprecation warnings
Naming convention for arrays
Commonly used indices:
index
description
aAtom number
cUnit cell axis-index (0, 1, 2)
vxyz-index (0, 1, 2)
KBZ k-point index
kIBZ k-point index
qIBZ k-point index (local, i.e. it starts at 0 on each processor)
sSpin index (\(\sigma\))
sSymmetry index
uCombined spin and k-point index (local)
RThree indices into the coarse 3D grid
rThree indices into the fine 3D grid
GIndex of plane-wave coefficient (wave function expansion,
ecut)
gIndex of plane-wave coefficient (densities,
2 * ecut)
hIndex of plane-wave coefficient (compensation charges,
8 * ecut)
X
RorG
x
r,gorh
xZero or more extra dimensions
MLCAO orbital index (\(\mu\))
nBand number
nPrincipal quantum number
lAngular momentum quantum number (s, p, d, …)
mMagnetic quantum number (0, 1, …, 2*`ell` - 1)
L
landm(L = l**2 + m)
jValence orbital number (
nandl)
iValence orbital number (
n,landm)
q
j1andj2pair
p
i1andi2pair
rCPU-rank
Examples:
|
\(D_{\sigma,i_1,i_2}^a\) |
|
|
\(\tilde{n}_\sigma(\mathbf{r})\) |
|
|
\(P_{\sigma \mathbf{k} ni}^a\) |
|
|
\(\tilde{\psi}_{\sigma \mathbf{k} n}(\mathbf{r})\) |
|
|
\(\tilde{p}_{\sigma \mathbf{k} i}^a(\mathbf{r}-\mathbf{R}^a)\) |
Domain descriptors
GPAW has two different container types for storing one or more functions in a unit cell (wave functions, electron densities, …):
Uniform grids
A uniform grid can be created with the UGDesc class:
>>> import numpy as np
>>> from gpaw.core import UGDesc
>>> a = 4.0
>>> n = 20
>>> grid = UGDesc(cell=a * np.eye(3),
... size=(n, n, n))
Given a UGDesc object, one can create
UGArray objects like this
>>> func_R = grid.empty()
>>> func_R.data.shape
(20, 20, 20)
>>> func_R.data[:] = 1.0
>>> grid.zeros((3, 2)).data.shape
(3, 2, 20, 20, 20)
Here are the methods of the UGDesc class:
Size of uniform grid. |
|
Actual size of uniform grid. |
|
Phase factor for block-boundary conditions. |
|
Create new uniform grid description. |
|
Create new UGArray object. |
|
… |
|
Yield views of blocks of data. |
|
Create array of (x, y, z) coordinates. |
|
Create UGAtomCenteredFunctions object. |
|
Create transformer from one grid to another. |
|
Plane wave. |
|
Create UGDesc from grid-spacing. |
|
Create FFTW-plans. |
|
… |
|
Maximum value of ekin so that all 0.5 * G^2 < ekin. |
|
… |
and the UGArray class:
Create new UniforGridFunctions object of same kind. |
|
Extract x, y values along line. |
|
Return a copy with dtype=complex. |
|
Scatter data from rank-0 to all ranks. |
|
Gather data from all ranks to rank-0. |
|
Do FFT. |
|
Calculate integral over cell of absolute value squared. |
|
Integral of self or self times cc(other). |
|
Convert to UniformGrid with |
|
Multiply by \(exp(ik.r)\). |
|
Interpolate to finer grid. |
|
Restrict to coarser grid. |
|
Add weighted absolute square of data to output array. |
|
Make data symmetric. |
|
Insert random numbers between -0.5 and 0.5 into data. |
|
Calculate moment of data. |
|
Create new scaled UGArray object. |
|
… |
|
… |
|
… |
|
… |
Plane waves
A set of plane-waves are characterized by a cutoff energy and a uniform grid:
>>> from gpaw.core import PWDesc
>>> pw = PWDesc(ecut=100, cell=grid.cell)
>>> func_G = pw.empty()
>>> func_R.fft(out=func_G)
PWArray(pw=PWDesc(ecut=100 <coefs=1536/1536>, cell=[4.0, 4.0, 4.0], pbc=[1, 1, 1], comm=0/1, dtype=float64), dims=())
>>> G = pw.reciprocal_vectors()
>>> G.shape
(1536, 3)
>>> G[0]
array([0., 0., 0.])
>>> func_G.data[0]
np.complex128(1+0j)
>>> func_G.ifft(out=func_R)
UGArray(grid=UGDesc(size=[20, 20, 20], cell=[4.0, 4.0, 4.0], pbc=[1, 1, 1], comm=0/1, dtype=float64), dims=())
>>> round(func_R.data[0, 0, 0], 15)
np.float64(1.0)
Here are the methods of the PWDesc class:
int([x]) -> integer |
|
Tuple with one element: number of plane waves. |
|
Returns reciprocal lattice vectors, G + k, in xyz coordinates. |
|
Kinetic energy of plane waves. |
|
Create new PlaneWaveExpanions object. |
|
… |
|
Create new plane-wave expansion description. |
|
Return indices into FFT-grid. |
|
… |
|
Cut out G-vectors with (G+k)^2/2<E_kin. |
|
Paste G-vectors with (G+k)^2/2<E_kin into 3-D FFT grid and |
|
Map from one (distributed) set of plane waves to smaller global set. |
|
Create PlaneWaveAtomCenteredFunctions object. |
and the PWArray class:
Create new PWArray object of same kind. |
|
Create a copy (surprise!). |
|
Sanity check for real-valued PW expansions. |
|
Matrix view of data. |
|
Do inverse FFT(s) to uniform grid(s). |
|
… |
|
Gather coefficients on master. |
|
Gather coefficients from self[r] on rank r. |
|
Scatter plane-wave coefficients from rank-0 to all ranks. |
|
Scatter everything from rank-0 to all ranks. |
|
Scatter all coefficients from rank r to self on other cores. |
|
Integral of self or self time cc(other). |
|
Calculate integral over cell. |
|
Add weighted absolute square of self to output array. |
|
… |
|
Insert random numbers between -0.5 and 0.5 into data. |
|
… |
|
Calculate average value at boundary of box. |
|
Transform expansion to new cell. |
|
… |
|
Symmetry-transform data. |
|
… |
Atoms-arrays
As an example, here is how to store the PAW atomic density-matrices for a water molecule (\(D_{\sigma,i_1,i_2}^a\)):
>>> nspins = 2
>>> D_asii = AtomArraysLayout(
... [(5, 5), (5, 5), (13, 13)],
... # dtype=float,
... # xp=np,
... # atomdist=...
... ).zeros(nspins)
>>> D_asii.data.shape
(2, 219)
Matrix elements
>>> psit_nG = pw.zeros(5)
>>> def T(psit_nG):
... """Kinetic energy operator."""
... out = psit_nG.new()
... out.data[:] = psit_nG.desc.ekin_G * psit_nG.data
... return out
>>> H_nn = psit_nG.matrix_elements(psit_nG, function=T)
Same as:
>>> Tpsit_nG = T(psit_nG)
>>> psit_nG.matrix_elements(Tpsit_nG, symmetric=True)
Matrix(float64: 5x5)
but faster.
Atom-centered functions
# creates: acf_example.png
import numpy as np
import matplotlib.pyplot as plt
from gpaw.core import UGDesc
alpha = 4.0
rcut = 2.0
l = 0
gauss = (l, rcut, lambda r: (4 * np.pi)**0.5 * np.exp(-alpha * r**2))
grid = UGDesc(cell=[4.0, 2.0, 2.0], size=[40, 20, 20])
pos = [[0.25, 0.5, 0.5], [0.75, 0.5, 0.5]]
acf_aR = grid.atom_centered_functions([[gauss], [gauss]], pos)
coef_as = acf_aR.empty(dims=(2,))
coef_as[0] = [[1], [-1]]
coef_as[1] = [[2], [1]]
print(coef_as.data, coef_as[0])
f_sR = grid.zeros(2)
acf_aR.add_to(f_sR, coef_as)
x = grid.xyz()[:, 10, 10, 0]
y1, y2 = f_sR.data[:, :, 10, 10]
ax = plt.subplot(1, 1, 1)
ax.plot(x, y1, 'o-')
ax.plot(x, y2, 'x-')
# plt.show()
plt.savefig('acf_example.png')
Matrix object
Here are the methods of the Matrix class:
Create new matrix of same shape and dtype. |
|
Create a copy. |
|
True if this matrix has nontrivial BLACS or GPU distribution. |
|
BLAS matrix-multiplication with other matrix. |
|
Redistribute to other BLACS layout. |
|
Gather the Matrix on rank zero. |
|
Scatter from rank-0. |
|
Construct a distributed Matrix object by scattering a raw 2D array |
|
Inplace inversion. |
|
In-place inverse of Cholesky decomposition. |
|
Calculate eigenvectors and eigenvalues. |
|
Solve generalized eigenvalue problem. |
|
Inplace complex conjugation. |
|
Add hermitian conjugate to myself. |
|
Fill in upper triangle from lower triangle. |
|
Add list of numbers or single number to diagonal of matrix. |
|
Create new matrix object with values transferred from GPU to CPU. |
|
Create new matrix object with data on GPU or CPU. |
|
Convert to new data type. |
A simple example that we can run with MPI on 4 cores:
from gpaw.core.matrix import Matrix
from gpaw.mpi import world
a = Matrix(5, 5, dist=(world, 2, 2, 2))
a.data[:] = world.rank
print(world.rank, a.data.shape)
Here, we have created a 5x5 Matrix of floats distributed on a 2x2
BLACS grid with a block size of 2 and we then print the shapes of the ndarrays,
which looks like this (in random order):
1 (2, 3)
2 (3, 2)
3 (2, 2)
0 (3, 3)
Let’s create a new matrix b and redistribute
from
a to b:
b = a.new(dist=(None, 1, 1, None))
a.redist(b)
if world.rank == 0:
print(b.array)
This will output:
[[ 0. 0. 2. 2. 0.]
[ 0. 0. 2. 2. 0.]
[ 1. 1. 3. 3. 1.]
[ 1. 1. 3. 3. 1.]
[ 0. 0. 2. 2. 0.]]
Matrix-matrix multiplication works like this:
c = a.multiply(a, opb='T')
API
Core
- class gpaw.core.UGDesc(*, cell, size, pbc=(True, True, True), zerobc=(False, False, False), kpt=None, comm=SerialCommunicator(), decomp=None, dtype=None)[source]
Description of 3D uniform grid.
- Parameters:
cell (ArrayLike1D | ArrayLike2D) – Unit cell given as three floats (orthorhombic grid), six floats (three lengths and the angles in degrees) or a 3x3 matrix (units: bohr).
size (ArrayLike1D) – Number of grid points along axes.
pbc – Periodic boundary conditions flag(s).
zerobc – Zero-boundary conditions flag(s). Skip first grid-point (assumed to be zero).
comm (MPIComm) – Communicator for domain-decomposition.
kpt (Vector | None) – K-point for Block-boundary conditions specified in units of the reciprocal cell.
decomp (Sequence[Sequence[int]] | None) – Decomposition of the domain.
dtype – Data-type (float or complex).
- atom_centered_functions(functions, positions, *, qspiral_v=None, atomdist=None, integrals=None, save_memory=True, cut=False, xp=None)[source]
Create UGAtomCenteredFunctions object.
- eikr(kpt_c=None, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>)[source]
Plane wave.
\[e^{i\mathbf{k}\cdot \mathbf{r}}\]
- ekin_max()[source]
Maximum value of ekin so that all 0.5 * G^2 < ekin.
In 1D, this will be 0.5*(pi/h)^2 where h is the grid-spacing.
- empty(dims=(), comm=SerialCommunicator(), xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>)[source]
Create new UGArray object.
- Parameters:
comm (_Communicator) – Distribute dimensions along this communicator.
- fft_plans(flags=0, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>, dtype=None)[source]
Create FFTW-plans.
- classmethod from_cell_and_grid_spacing(cell, grid_spacing, pbc=(True, True, True), kpt=None, comm=SerialCommunicator(), dtype=None)[source]
Create UGDesc from grid-spacing.
- gradient_operator(v, *, scale=1.0, n=1, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>)[source]
- new(*, kpt=None, dtype=None, comm='inherit', size=None, pbc=None, zerobc=None, decomp=None)[source]
Create new uniform grid description.
- property phase_factor_cd
Phase factor for block-boundary conditions.
- property size
Size of uniform grid.
- class gpaw.core.PWDesc(*, ecut=None, gcut=None, cell, kpt=None, comm=SerialCommunicator(), dtype=None)[source]
Description of plane-wave basis.
- Parameters:
ecut (float | None) – Cutoff energy for kinetic energy of plane waves (units: hartree).
cell (ArrayLike1D | ArrayLike2D) – Unit cell given as three floats (orthorhombic grid), six floats (three lengths and the angles in degrees) or a 3x3 matrix (units: bohr).
comm (MPIComm) – Communicator for distribution of plane-waves.
kpt (Vector | None) – K-point for Block-boundary conditions specified in units of the reciprocal cell.
dtype – Data-type (float or complex).
- atom_centered_functions(functions, positions, *, qspiral_v=None, atomdist=None, integrals=None, save_memory=True, cut=False, xp=None)[source]
Create PlaneWaveAtomCenteredFunctions object.
- empty(dims=(), comm=SerialCommunicator(), xp=None)[source]
Create new PlaneWaveExpanions object.
- Parameters:
comm (_Communicator) – Distribute dimensions along this communicator.
- map_indices(other)[source]
Map from one (distributed) set of plane waves to smaller global set.
Say we have 9 G-vector on two cores:
5 3 4 . 3 4 0 . . 2 0 1 -> rank=0: 2 0 1 rank=1: . . . 8 6 7 . . . 3 1 2
and we want a mapping to these 5 G-vectors:
3 2 0 1 4
On rank=0: the return values are:
[0, 1, 2, 3], [[0, 1, 2, 3], [4]]
and for rank=1:
[1], [[0, 1, 2, 3], [4]]
- new(*, ecut=None, gcut=None, kpt=None, dtype=None, comm='inherit')[source]
Create new plane-wave expansion description.
- class gpaw.core.atom_centered_functions.AtomCenteredFunctions(functions, relpos_ac, atomdist=None, xp=None)[source]
- property atomdist
- derivative(functions, out=None)[source]
Calculate derivatives of integrals with respect to atom positions.
- integrate(functions, out=None, add_to=False)[source]
Calculate integrals of atom-centered functions multiplied by functions.
- property layout
- class gpaw.core.UGArray(grid, dims=(), comm=SerialCommunicator(), data=None, xp=None)[source]
Object for storing function(s) on a uniform grid.
- Parameters:
- fft(plan=None, pw=None, out=None)[source]
Do FFT.
- Returns:
PWArray with values
\[C(\mathbf{G}) = \frac{1}{Ω} \int d\mathbf{r} e^{-i\mathbf{G}\cdot \mathbf{r}} f(\mathbf{r}),\]where \(C(\bG)\) are the plane wave coefficients and Ω is the cell volume.
- Parameters:
plan – Plan for FFT.
pw – Target PW description.
out – Target PWArray object.
- new(data=None, zeroed=False, dims=None)[source]
Create new UniforGridFunctions object of same kind.
- Parameters:
data – Array to use for storage.
zeroed – If True, set data to zero.
dims – Extra dimensions (bands, spin, etc.), required if data does not fit the full array.
- norm2(kind='normal', weights=None, skip_sum=False)[source]
Calculate integral over cell of absolute value squared.
\[\int |a(\mathbf{r})|^{2} d\mathbf{r}\]
- redist(domain, comm1, comm2)[source]
Redistribute to new domain.
The “world” is spanned by:
(self.desc.comm, comm1)
and:
(domain.comm, comm2).
- class gpaw.core.arrays.XArray(dims, myshape, comm, domain_comm, data, dv, dtype, xp=None)[source]
- abs_square(weights, out)[source]
Add weighted absolute square of data to output array.
See also XKCD: 849.
- create_work_buffer(data_buffer)[source]
Create new Distributed array object of same kind, to be used as a buffer array when doing sliced operations.
- Parameters:
data_buffer (ndarray) – Array to use for storage.
- desc: DomainType
- matrix_elements(other, *, out=None, symmetric='_default', function=None, domain_sum=True, cc=False)[source]
- class gpaw.core.atom_arrays.AtomArrays(layout, dims=(), comm=SerialCommunicator(), data=None)[source]
AtomArrays object.
- Parameters:
layout (AtomArraysLayout) – Layout-description.
comm (MPIComm) – Distribute dimensions along this communicator.
data (np.ndarray | None) – Data array for storage.
- block_diag_multiply(block_diag_matrix_axii, out_ani, index=None)[source]
Multiply by block diagonal matrix.
with A, B and C referring to
self,block_diag_matrix_axiiandout_ani:\[\sum^{}_{i} A^{a}_{ni} B^{a}_{ij} \rightarrow C^{a}_{nj}\]If index is not None,
block_diag_matrix_axiimust have an extra dimension: \(B_{ij}^{ax}\) and x=index is used.
- gather(broadcast: Literal[False] = False, copy: bool = False) AtomArrays | None[source]
- gather(broadcast: Literal[True], copy: bool = False) AtomArrays
Gather all atoms on master.
- new(*, comm='inherit', layout=None, data=None, xp=None)[source]
Create new AtomArrays object of same kind.
- Parameters:
layout – Layout-description.
data – Array to use for storage.
- to_full()[source]
Convert \(N(N+1)/2\) vectors to \(N\times N\) matrices.
>>> a = AtomArraysLayout([6]).empty() >>> a[0][:] = [1, 2, 3, 4, 5, 6] >>> a.to_full()[0] array([[1., 2., 4.], [2., 3., 5.], [4., 5., 6.]])
- class gpaw.core.atom_arrays.AtomArraysLayout(shapes, atomdist=SerialCommunicator(), dtype=<class 'float'>, xp=None)[source]
Description of layout of atom arrays.
- Parameters:
shapes (Sequence[int | tuple[int, ...]]) – Shapse of arrays - one for each atom.
atomdist (AtomDistribution | MPIComm) – Distribution of atoms.
dtype – Data-type (float or complex).
- empty(dims=(), comm=SerialCommunicator())[source]
Create new AtomArrays object.
- Parameters:
comm (_Communicator) – Distribute dimensions along this communicator.
- new(atomdist=None, dtype=None, xp=None)[source]
Create new AtomsArrayLayout object with new atomdist.
- class gpaw.core.atom_arrays.AtomDistribution(ranks, comm=SerialCommunicator())[source]
Atom-distribution.
- Parameters:
ranks (ArrayLike1D) – List of ranks, one rank per atom.
comm (MPIComm) – MPI-communicator.
- classmethod from_atom_indices(atom_indices, comm=SerialCommunicator(), *, natoms=None)[source]
Create distribution from atom indices.
>>> AtomDistribution.from_atom_indices([0, 1, 2]).rank_a array([0, 0, 0])
- class gpaw.core.PWArray(pw, dims=(), comm=SerialCommunicator(), data=None, xp=None)[source]
Object for storing function(s) as a plane-wave expansions.
- Parameters:
- abs_square(weights, out, _slow=False)[source]
Add weighted absolute square of self to output array.
With \(a_n(G)\) being self and \(w_n\) the weights:
\[out(\mathbf{r}) \leftarrow out(\mathbf{r}) + \sum^{}_{n} |FFT^{-1} [a_{n} (\mathbf{G})]|^{2} w_{n}\]
- gather_all(out)[source]
Gather coefficients from self[r] on rank r.
On rank r, an array of all G-vector coefficients will be returned. These will be gathered from self[r] on all the cores.
- ifft(*, plan=None, grid=None, grid_spacing=None, out=None, periodic=False)[source]
Do inverse FFT(s) to uniform grid(s).
- Returns:
UGArray with values
\[f(\mathbf{r}) = \sum^{}_{G} c(G) e ^{i\mathbf{G}\cdot \mathbf{R}}\]- Parameters:
plan – Plan for inverse FFT.
grid – Target grid.
out – Target UGArray object.
- new(data=None, dims=None)[source]
Create new PWArray object of same kind.
- Parameters:
data – Array to use for storage.
dims – Extra dimensions (bands, spin, etc.), required if data does not fit the full array.
- norm2(kind='normal', weights=None, skip_sum=False)[source]
Calculate integral over cell.
For kind=’normal’ we calculate:
\[\int |a(\mathbf{r})|^{2} d\mathbf{r} = \sum^{}_{G} |c_{G} |^{2} V,\]where V is the volume of the unit cell.
And for kind=’kinetic’:
\[\frac{1}{2} \sum^{}_{G} |c_{G} |^{2} G^{2} V,\]
- class gpaw.core.matrix.Matrix(M, N, *, dtype=None, data=None, dist=None, xp=None)[source]
Matrix object.
- Parameters:
M (int) – Rows.
N (int) – Columns.
dtype – Data type (float or complex).
dist (MatrixDistribution | MPIComm | tuple | None) – BLACS distribution given as (communicator, rows, columns, blocksize) tuple. Default is None meaning no distribution.
data (ArrayLike2D | None) – Numpy ndarray to use for storage. By default, a new ndarray will be allocated.
- eigh(S=None, *, cc=False, scalapack=(SerialCommunicator(), 1, 1, 0), limit=None)[source]
Calculate eigenvectors and eigenvalues.
Matrix must be symmetric/hermitian and stored in lower half. If
Sis given, solve a generalized eigenvalue problem.
- eighl(L, comm2=SerialCommunicator())[source]
Solve generalized eigenvalue problem.
With \(H\) being self, we solve for the eigenvectors \(C\) and the eigenvalues \(Λ\) (a diagonal matrix):
\[HC = SCΛ,\]where \(L\) is a lower triangle matrix such that:
\[LSL^{†} = 1\cdot\]The solution has these three steps:
\[\tilde{H} = LHL^{†} , \tilde{H}\tilde{C} = \tilde{C}Λ, C = L^{†} \tilde{C}\cdot\]Note that \(H\) must be the full matrix not just half of it!
- gather(broadcast: Literal[True] = True) Matrix[source]
- gather(broadcast: Literal[False] = False) Matrix | None
Gather the Matrix on rank zero.
- invcholesky()[source]
In-place inverse of Cholesky decomposition.
Calculate a lower triangle matrix \(L\) where:
\[LSL^{†} = 1,\]and \(S\) is self. Only the lower part of \(S\) is used.
>>> S = Matrix(2, 2, data=[[1.0, np.nan], ... [0.1, 1.0]]) >>> S.invcholesky() >>> S.data array([[ 1. , 0. ], [-0.10050378, 1.00503782]])
- multiply(other, alpha=1.0, opa='N', opb='N', out=None, data_buffer=None, beta=0.0, symmetric=False)[source]
BLAS matrix-multiplication with other matrix.
- new(dist='inherit', data=None)[source]
Create new matrix of same shape and dtype.
Default is to use same BLACS distribution. Use dist to use another distribution.
- redist(other)[source]
Redistribute to other BLACS layout. \(other\) is the output, newly distributed matrix.
Input-parameter objects
DFT-components
- class gpaw.new.calculation.DFTCalculation(atoms, ibzwfs, density, potential, setups, scf_loop, pot_calc, log, params, energies=None)[source]
- change_mode(mode, *, nbands=None)[source]
In-place convertion from one mode to another.
Only LCAO to PW or FD mode implemented!
- converge(maxiter=None, steps=99999999999999999)[source]
Converge to self-consistent solution of Kohn-Sham equation.
- classmethod from_gpw_file(filename, log=None, comm=None, parallel=None, dtype=None, force_complex_dtype=False, object_hooks=None)[source]
- classmethod from_parameters(atoms, params, comm, log=None)[source]
Create DFTCalculation object from parameters and atoms.
- gather(txt='-')[source]
Gather calculation data from DFTCalculation object on master and return new DFTCalculation (only on master, None everywhere else).
- property state
- class gpaw.new.calculation.DFTState(ibzwfs, density, potential, energies)[source]
State of a Kohn-Sham calculation.
- class gpaw.new.density.Density(nt_sR, taut_sR, D_asii, charge, nvalence, delta_aiiL, delta0_a, N0_aii, n_aj, l_aj, nct_aX, tauct_aX)[source]
- classmethod from_superposition(*, grid, nct_aX, tauct_aX, atomdist, setups, basis_set, magmom_av, ncomponents, charge=0.0, hund=False, mgga=False)[source]
- property nct_R
- property tauct_R
- class gpaw.new.ibzwfs.IBZWaveFunctions(ibz, *, ncomponents, wfs_u, kpt_comm=SerialCommunicator(), kpt_band_comm=SerialCommunicator(), comm=SerialCommunicator())[source]
Collection of wave function objects for k-points in the IBZ.
- classmethod create(*, ibz, ncomponents, create_wfs_func, kpt_comm=SerialCommunicator(), kpt_band_comm=SerialCommunicator(), comm=SerialCommunicator())[source]
- get_all_electron_wave_function(band, kpt=0, spin=0, grid_spacing=0.05, skip_paw_correction=False)[source]
- get_max_shape(global_shape=False)[source]
Find the largest wave function array shape.
For a PW-calculation, this shape could depend on k-point.
- property mode
- class gpaw.new.potential.Potential(vt_sR, dH_asii, dedtaut_sR, vHt_x=None, e_stress=nan)[source]
- class gpaw.new.pot_calc.PotentialCalculator(xc, poisson_solver, setups, *, relpos_ac, extensions=None, soc=False)[source]
- property charge
- class gpaw.new.scf.SCFLoop(hamiltonian, occ_calc, eigensolver, mixer, comm, convergence, maxiter)[source]
- class gpaw.dft.Parameters(*, mode, basis=None, charge=None, convergence=None, eigensolver=None, experimental=None, extensions=None, gpts=None, h=None, hund=None, interpolation=None, kpts=None, magmoms=None, maxiter=None, mixer=None, nbands=None, occupations=None, parallel=None, poissonsolver=None, random=None, setups=None, soc=None, spinpol=None, symmetry=None, xc=None, external=None)[source]
DFT-parameters object.
>>> p = Parameters(mode=PW(400)) >>> p mode=PW(ecut=400) >>> p.charge 0.0 >>> p.xc XC(name='LDA') >>> from ase.build import molecule >>> atoms = molecule('H2', vacuum=3.0) >>> dft = p.dft_calculation(atoms, txt=None) >>> atoms.calc = dft.ase_calculator()
- Parameters:
basis (str | dict[str | int | None, str] | None) – Basis-set. Used for LCAO calculations and wave-function initial guess for PW and FD calculations. Default is to use the PAW pseudo partial-waves.
charge (float | None) – Total charge of the system in units of \(|e|\).
convergence (dict | None) – SCF-convergence criteria.
eigensolver (str | dict | Eigensolver | None) – Eigensolver. Default for PW and FD mode is
'davidson'.gpts (Sequence[int] | None) – Number of real-space grid-points for wave-functions (three integers).
h (float | None) –
Grid-spacing for wave-function grid (Å). Default value is 0.2 Å for LCAO and FD mode calculations. For a PW-mode calculation, we use the formula \(h=γh_0\) with \(γ \simeq 1.4\) and:
\[h_0 = \frac{\pi}{\sqrt{8E_c}}.\]Ideally, we would use \(\gamma=1\), but in practice, 1.4 is a good compromise between accuracy and efficiency. In eV and Å units we have \(h_0=3.07/\sqrt{E_c}\).
hund (bool | None) – Use Hund’s rule for initial magnetic moments.
experimental (dict | None) – Experimental stuff.
extensions (Sequence[ExtensionInput] | None) – Extensions (D3, …).
interpolation (int | Literal['fft'] | None) –
…
kpts (KptsType | MonkhorstPack | None) – Brillouin-zone sampling. Default is Γ-point only.
magmoms (Sequence[float] | Sequence[Sequence[float]] | None) – Initial magnetic moments for non-collinear calculations.
maxiter (int | None) – Maximum number of allowed SCF-iterations. Default is 333.
occupations (dict | Occupations | None) –
…
parallel (dict | None) – Parallelization strategy. Example: Force parallelization over
'kptwith{'band': 1, 'domain': 1}.poissonsolver (dict | PoissonSolver | None) –
…
random (bool | None) – Use random numbers for initial wave functions.
…
soc (bool | None) – Enable spin-orbit coupling.
spinpol (bool | None) – Force spin-polarized calculation.
symmetry (str | dict | Symmetry | None) – Use of symmetry. Default is to use …
xc (str | dict | XC | None) – XC-functional. Default is PZ-LDA.
- class gpaw.new.pwfd.wave_functions.PWFDWaveFunctions(psit_nX, *, spin, q, k, setups, relpos_ac, atomdist, domain_band_comm, weight=1.0, ncomponents=1, qspiral_v=None)[source]
- property P_ani: AtomArrays
PAW projections.
\[\langle \tilde{p}^{a}_{i}|\tilde{ψ}_{n} \rangle\]
- build_hamiltonian(Ht, dH, psit2_nX, calculate_energy=True)[source]
psit2_nX will be used as a buffer for the wave functions.
Ht(in, out):
\[\tilde{H} = \hat{T} + \tilde{v}\]dH:
\[\langle \tilde{𝜓}_{m}|\tilde{p} \rangle ΔH^{a}_{ij} \langle \tilde{p}_{j}|\tilde{𝜓} \rangle\]
- collect_bands_and_domain(n1=0, n2=0)[source]
Collect range of bands to master of band and domain comms.
- copy()[source]
Make a copy of the wave functions.
The buffer for the wave functions themselves is copied, so that it can be modified without changing the original.
Useful, e.g. in RTTDDFT, where the wave functions are copied in order to be restored later.
- dipole_matrix_elements()[source]
Calculate dipole matrix-elements.
\[\mathbf{μ}_{mn} = \int \tilde{𝜓}_{m} \tilde{𝜓}_{n} \mathbf{r}d\mathbf{r} + \sum^{}_{aij} P^{a}_{im} P^{a}_{jn} Δ\mathbf{μ}^{a}_{ij}\]- Returns:
matrix elements in atomic units.
- Return type:
Array3D
- orthonormalize(psit2_nX=None)[source]
Orthonormalize wave functions.
Computes the overlap matrix:
\[S_{mn} = \int \tilde{ψ}_{m}(\mathbf{r})^{*} \tilde{ψ}_{n}(\mathbf{r}) d\mathbf{r} + \sum^{}_{aij} (P^{a}_{im} )^{*} P^{a}_{jn} ΔS^{a}_{ij}\]With \(LSL^\dagger=1\), we update the wave functions and projections inplace like this:
\[Ψ_{m} \leftarrow \sum^{}_{n} L^{*}_{mn} Ψ_{n} ,\]and:
\[P^{a}_{mi} \leftarrow \sum^{}_{n} L^{*}_{mn} P^{a}_{ni} \cdot\]
- property pt_aiX: AtomCenteredFunctions
PAW projector functions.
\[\tilde{p}^{a}_{i} (\mathbf{r})\]
- subspace_diagonalize(Ht, dH, psit2_nX, data_buffer=None, scalapack_parameters=(None, 1, 1, 0), nocc=None, eigenvalues_only=False, calculate_energy=True)[source]
If data_buffer is None, psit2_nX will be used as a buffer for the wave functions.
Ht(in, out):
\[\tilde{H} = \hat{T} + \tilde{v}\]dH:
\[\langle \tilde{𝜓}_{m}|\tilde{p} \rangle ΔH^{a}_{ij} \langle \tilde{p}_{j}|\tilde{𝜓} \rangle\]
- class gpaw.new.ase_interface.ASECalculator(params, *, log, dft=None, atoms=None)[source]
This is the ASE-calculator frontend for doing a GPAW calculation.
- calculate_property(atoms, prop)[source]
Calculate (if not already calculated) a property.
The
propstring must be one ofenergy
forces
stress
magmom
magmoms
dipole
- property density
- property dft: DFTCalculation
- dos(soc=False, theta=0.0, phi=0.0, shift_fermi_level=True)[source]
Create DOS-calculator.
Default is to
shift_fermi_levelto 0.0 eV. Forsoc=True, angles can be given in degrees.
- get_wannier_localization_matrix(nbands, dirG, kpoint, nextkpoint, G_I, spin)[source]
Calculate integrals for maximally localized Wannier functions.
- property hamiltonian
- iconverge(atoms, *, need_wfs=False)[source]
Iterate to self-consistent solution.
Will also calculate “cheap” properties: energy, magnetic moments and dipole moment.
- implemented_properties = ['energy', 'free_energy', 'forces', 'stress', 'dipole', 'magmom', 'magmoms']
- property initialized
- name = 'gpaw'
- old = False
- property parameters
- property results
- property setups
- property spos_ac
- property symmetry
- property wfs
- property world
- gpaw.dft.DFT(atoms, *, mode, basis=None, charge=None, convergence=None, eigensolver=None, experimental=None, extensions=None, gpts=None, h=None, hund=None, interpolation=None, kpts=None, magmoms=None, maxiter=None, mixer=None, nbands=None, occupations=None, parallel=None, poissonsolver=None, random=None, setups=None, soc=None, spinpol=None, symmetry=None, xc=None, txt='-', communicator=None)[source]
Create a DFTCalculation object.
See
gpaw.dft.Parametersfor the complete list of parameters.
- gpaw.dft.GPAW(filename=None, *, basis=None, charge=None, convergence=None, eigensolver=None, experimental=None, extensions=None, gpts=None, h=None, hund=None, interpolation=None, kpts=None, magmoms=None, maxiter=None, mixer=None, mode=None, nbands=None, occupations=None, parallel=None, poissonsolver=None, random=None, setups=None, soc=None, spinpol=None, symmetry=None, xc=None, txt='?', communicator=None, object_hooks=None, legacy_gpaw=None, external=None, background_charge=None)[source]
Create ASE-compatible GPAW calculator.
See
gpaw.dft.Parametersfor the complete list of parameters.- Parameters:
filename (str | Path | IO[str] | None) – Name of gpw-file to restart from.
txt (str | Path | IO[str] | None) – Text log-file. Use
Nonefor no logging and'-'for using standard out.communicator (MPIComm | None) – MPI-communicator. Default is to use
gpaw.mpi.world.object_hooks – Dictionary of hook-functions to create custom parameter-objects.
- gpaw.new.pwfd.move_wfs.move_wave_functions(oldrelpos_ac, newrelpos_ac, P_ani, psit_nX, setups)[source]
Move wavefunctions with atoms according to PAW basis
Wavefunctions are approximated as:
\[\tilde{ψ}_{n}(\mathbf{r}) = \sum^{}_{ai} \tilde{φ}^{a}_{i} (\mathbf{r}) \langle \tilde{p}^{a}_{i}|\tilde{ψ}_{n} \rangle ,\]where i runs over the bound partial-waves only. This quantity is then subtracted and re-added at the new positions.
- class gpaw.new.symmetry.Symmetries(*, cell, rotations=None, translations=None, atommaps=None, tolerance=None, _backwards_compatible=False)[source]
Symmetries object.
“Rotations” here means rotations, mirror and inversion operations.
Units of “cell” and “tolerance” should match.
>>> sym = Symmetries.from_cell([1, 2, 3]) >>> sym.has_inversion True >>> len(sym) 8 >>> sym2 = sym.analyze_positions([[0, 0, 0], [0, 0, 0.4]], ids=[1, 2]) >>> sym2.has_inversion False >>> len(sym2) 4
- classmethod from_cell(cell, *, pbc=(True, True, True), tolerance=None, guarantee_group=True, _backwards_compatible=False)[source]
- classmethod from_cell_and_atoms(cell, *, pbc=(True, True, True), tolerance=None, _backwards_compatible=False, relative_positions, ids, symmorphic=True, guarantee_group=True)[source]
- classmethod from_cell_and_atoms_spglib(cell, *, pbc=(True, True, True), tolerance=None, _backwards_compatible=False, relative_positions, ids, symmorphic=True)[source]
- property gcd_c
- property has_inversion
- property symmorphic
- class gpaw.new.brillouin.IBZ(symmetries, bz, ibz2bz, bz2ibz, weights, bz2bz_Ks=None, s_K=None, time_reversal_K=None)[source]
FFTW
Python wrapper for FFTW3 library
- class gpaw.fftw.NumpyFFTPlans(size_c, dtype, empty=<function empty>)[source]
Numpy fallback.
- gpaw.fftw.check_fft_size(n, factors=[2, 3, 5, 7])[source]
Check if n is an efficient fft size.
Efficient means that n can be factored into small primes (2, 3, 5, 7).
>>> check_fft_size(17) False >>> check_fft_size(18) True
- gpaw.fftw.create_plans(size_c, dtype, flags=0, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>)[source]
Create plan-objects for FFT and inverse FFT.
- gpaw.fftw.empty(shape, dtype=<class 'float'>)[source]
numpy.empty() equivalent with 16 byte alignment.