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.DFTCalculation object if not already done

  • update 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:

density

gpaw.new.density.Density

ibzwfs

gpaw.new.ibzwfs.IBZWaveFunctions

potential

gpaw.new.potential.Potential

scf_loop

gpaw.new.scf.SCFLoop

pot_calc

gpaw.new.pot_calc.PotentialCalculator

Overview:

  • atoms.calc:

    • dft:

      • density: nt_sR, D_asii, …

      • potential: vt_sR, dH_asii, …

      • pot_calc:

        • xc

        • poisson_solver

      • scf_loop:

        • eigensolver

        • hamiltonian

        • occ_calc

        • mixer

      • ibzwfs:

        • ibz:

          • symmetries

          • bz

        • _wfs_u[u]:

          • psit_nX

          • occ_n

See also: code.svg.

There are three ways to create a DFTCalculation object:

Naming convention for arrays

Commonly used indices:

index

description

a

Atom number

c

Unit cell axis-index (0, 1, 2)

v

xyz-index (0, 1, 2)

K

BZ k-point index

k

IBZ k-point index

q

IBZ k-point index (local, i.e. it starts at 0 on each processor)

s

Spin index (\(\sigma\))

s

Symmetry index

u

Combined spin and k-point index (local)

R

Three indices into the coarse 3D grid

r

Three indices into the fine 3D grid

G

Index of plane-wave coefficient (wave function expansion, ecut)

g

Index of plane-wave coefficient (densities, 2 * ecut)

h

Index of plane-wave coefficient (compensation charges, 8 * ecut)

X

R or G

x

r, g or h

x

Zero or more extra dimensions

M

LCAO orbital index (\(\mu\))

n

Band number

n

Principal quantum number

l

Angular momentum quantum number (s, p, d, …)

m

Magnetic quantum number (0, 1, …, 2*`ell` - 1)

L

l and m (L = l**2 + m)

j

Valence orbital number (n and l)

i

Valence orbital number (n, l and m)

q

j1 and j2 pair

p

i1 and i2 pair

r

CPU-rank

Examples:

density.D_asii

\(D_{\sigma,i_1,i_2}^a\)

AtomArrays

density.nt_sR

\(\tilde{n}_\sigma(\mathbf{r})\)

UGArray

ibzwfs._wfs_u[u].P_ani

\(P_{\sigma \mathbf{k} ni}^a\)

AtomArrays

ibzwfs._wfs_u[u].psit_nX

\(\tilde{\psi}_{\sigma \mathbf{k} n}(\mathbf{r})\)

UGArray | PWArray

ibzwfs._wfs_u[u].pt_aX

\(\tilde{p}_{\sigma \mathbf{k} i}^a(\mathbf{r}-\mathbf{R}^a)\)

AtomCenteredFunctions

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()

Size of uniform grid.

global_shape()

Actual size of uniform grid.

phase_factor_cd()

Phase factor for block-boundary conditions.

new()

Create new uniform grid description.

empty()

Create new UGArray object.

from_data()

blocks()

Yield views of blocks of data.

xyz()

Create array of (x, y, z) coordinates.

atom_centered_functions()

Create UGAtomCenteredFunctions object.

transformer()

Create transformer from one grid to another.

eikr()

Plane wave.

from_cell_and_grid_spacing()

Create UGDesc from grid-spacing.

fft_plans()

Create FFTW-plans.

ranks_from_fractional_positions()

ekin_max()

Maximum value of ekin so that all 0.5 * G^2 < ekin.

gradient_operator()

and the UGArray class:

new()

Create new UniforGridFunctions object of same kind.

xy()

Extract x, y values along line.

to_complex()

Return a copy with dtype=complex.

scatter_from()

Scatter data from rank-0 to all ranks.

gather()

Gather data from all ranks to rank-0.

fft()

Do FFT.

norm2()

Calculate integral over cell of absolute value squared.

integrate()

Integral of self or self times cc(other).

to_pbc_grid()

Convert to UniformGrid with pbc=(True, True, True).

multiply_by_eikr()

Multiply by \(exp(ik.r)\).

interpolate()

Interpolate to finer grid.

fft_restrict()

Restrict to coarser grid.

abs_square()

Add weighted absolute square of data to output array.

symmetrize()

Make data symmetric.

randomize()

Insert random numbers between -0.5 and 0.5 into data.

moment()

Calculate moment of data.

scaled()

Create new scaled UGArray object.

add_ked()

redist()

isosurface()

trace_inner_product()

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:

itemsize()

int([x]) -> integer

global_shape()

Tuple with one element: number of plane waves.

reciprocal_vectors()

Returns reciprocal lattice vectors, G + k, in xyz coordinates.

kinetic_energies()

Kinetic energy of plane waves.

empty()

Create new PlaneWaveExpanions object.

from_data()

new()

Create new plane-wave expansion description.

indices()

Return indices into FFT-grid.

minimal_uniform_grid()

cut()

Cut out G-vectors with (G+k)^2/2<E_kin.

paste()

Paste G-vectors with (G+k)^2/2<E_kin into 3-D FFT grid and

map_indices()

Map from one (distributed) set of plane waves to smaller global set.

atom_centered_functions()

Create PlaneWaveAtomCenteredFunctions object.

and the PWArray class:

new()

Create new PWArray object of same kind.

copy()

Create a copy (surprise!).

sanity_check()

Sanity check for real-valued PW expansions.

matrix()

Matrix view of data.

ifft()

Do inverse FFT(s) to uniform grid(s).

interpolate()

gather()

Gather coefficients on master.

gather_all()

Gather coefficients from self[r] on rank r.

scatter_from()

Scatter plane-wave coefficients from rank-0 to all ranks.

scatter_everything_from()

Scatter everything from rank-0 to all ranks.

scatter_from_all()

Scatter all coefficients from rank r to self on other cores.

integrate()

Integral of self or self time cc(other).

norm2()

Calculate integral over cell.

abs_square()

Add weighted absolute square of self to output array.

to_pbc_grid()

randomize()

Insert random numbers between -0.5 and 0.5 into data.

moment()

boundary_value()

Calculate average value at boundary of box.

morph()

Transform expansion to new cell.

add_ked()

transform()

Symmetry-transform data.

trace_inner_product()

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')
../_images/acf_example.png

Matrix object

Here are the methods of the Matrix class:

new()

Create new matrix of same shape and dtype.

copy()

Create a copy.

is_distributed()

True if this matrix has nontrivial BLACS or GPU distribution.

multiply()

BLAS matrix-multiplication with other matrix.

redist()

Redistribute to other BLACS layout.

gather()

Gather the Matrix on rank zero.

scatter_from()

Scatter from rank-0.

scatter()

Construct a distributed Matrix object by scattering a raw 2D array

inv()

Inplace inversion.

invcholesky()

In-place inverse of Cholesky decomposition.

eigh()

Calculate eigenvectors and eigenvalues.

eighl()

Solve generalized eigenvalue problem.

complex_conjugate()

Inplace complex conjugation.

add_hermitian_conjugate()

Add hermitian conjugate to myself.

tril2full()

Fill in upper triangle from lower triangle.

add_to_diagonal()

Add list of numbers or single number to diagonal of matrix.

to_cpu()

Create new matrix object with values transferred from GPU to CPU.

to_xp()

Create new matrix object with data on GPU or CPU.

to_dtype()

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.

blocks(data)[source]

Yield views of blocks of data.

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}}\]
Parameters:

kpt_c (Sequence[float] | ndarray | None) – k-point in units of the reciprocal cell. Defaults to the UGDesc objects own k-point.

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:
  • dims (int | tuple[int, ...]) – Extra dimensions.

  • 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.

from_data(data)[source]
global_shape()[source]

Actual size of uniform grid.

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.

ranks_from_fractional_positions(relpos_ac)[source]
property size

Size of uniform grid.

transformer(other, stencil_range=3, xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>)[source]

Create transformer from one grid to another.

(for interpolation and restriction).

xyz()[source]

Create array of (x, y, z) coordinates.

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.

cut(array_Q)[source]

Cut out G-vectors with (G+k)^2/2<E_kin.

empty(dims=(), comm=SerialCommunicator(), xp=None)[source]

Create new PlaneWaveExpanions object.

Parameters:
  • dims (int | tuple[int, ...]) – Extra dimensions.

  • comm (_Communicator) – Distribute dimensions along this communicator.

from_data(data)[source]
global_shape()[source]

Tuple with one element: number of plane waves.

indices(shape)[source]

Return indices into FFT-grid.

itemsize: int = 16
kinetic_energies()[source]

Kinetic energy of plane waves.

\[|\mathbf{G}+\mathbf{k}|^{2} / 2\]
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]]
minimal_uniform_grid(n=1, factors=(2, 3, 5, 7))[source]
new(*, ecut=None, gcut=None, kpt=None, dtype=None, comm='inherit')[source]

Create new plane-wave expansion description.

paste(coef_G, array_Q)[source]

Paste G-vectors with (G+k)^2/2<E_kin into 3-D FFT grid and zero-pad.

reciprocal_vectors(xp=<module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/gpaw/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>)[source]

Returns reciprocal lattice vectors, G + k, in xyz coordinates.

class gpaw.core.atom_centered_functions.AtomCenteredFunctions(functions, relpos_ac, atomdist=None, xp=None)[source]
add_to(functions, coefs=1.0)[source]

Add atom-centered functions multiplied by coefs to functions.

property atomdist
derivative(functions, out=None)[source]

Calculate derivatives of integrals with respect to atom positions.

empty(dims=(), comm=SerialCommunicator())[source]

Create AtomsArray for coefficients.

integrate(functions, out=None, add_to=False)[source]

Calculate integrals of atom-centered functions multiplied by functions.

property layout
move(relpos_ac, atomdist)[source]

Move atoms to new positions.

new(desc, atomdist)[source]
stress_contribution(a, c=1.0)[source]
class gpaw.core.UGArray(grid, dims=(), comm=SerialCommunicator(), data=None, xp=None)[source]

Object for storing function(s) on a uniform grid.

Parameters:
  • grid (UGDesc) – Description of uniform grid.

  • dims (int | tuple[int, ...]) – Extra dimensions.

  • comm (MPIComm) – Distribute dimensions along this communicator.

  • data (np.ndarray | None) – Data array for storage.

abs_square(weights, out=None)[source]

Add weighted absolute square of data to output array.

add_ked(occ_n, taut_R)[source]

Add weighted absolute square of gradient of data to output array.

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.

fft_restrict(plan1=None, plan2=None, grid=None, out=None)[source]

Restrict to coarser grid.

Parameters:
  • plan1 (FFTPlans | None) – Plan for FFT.

  • plan2 (FFTPlans | None) – Plan for inverse FFT.

  • grid (UGDesc | None) – Target grid.

  • out (UGArray | None) – Target UGArray object.

gather(out=None, broadcast=False)[source]

Gather data from all ranks to rank-0.

integrate(other=None, skip_sum=False)[source]

Integral of self or self times cc(other).

interpolate(*, plan1=None, plan2=None, grid=None, out=None)[source]

Interpolate to finer grid.

Parameters:
  • plan1 (FFTPlans | None) – Plan for FFT (course grid).

  • plan2 (FFTPlans | None) – Plan for inverse FFT (fine grid).

  • grid (UGDesc | None) – Target grid.

  • out (UGArray | None) – Target UGArray object.

isosurface(show=True, **kwargs)[source]
moment()[source]

Calculate moment of data.

multiply_by_eikr(kpt_c=None)[source]

Multiply by \(exp(ik.r)\).

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}\]
randomize(seed=None)[source]

Insert random numbers between -0.5 and 0.5 into data.

redist(domain, comm1, comm2)[source]

Redistribute to new domain.

The “world” is spanned by:

(self.desc.comm, comm1)

and:

(domain.comm, comm2).
scaled(cell, values=1.0)[source]

Create new scaled UGArray object.

Unit cell axes are multiplied by \(cell\) and data by \(values\).

scatter_from(data=None)[source]

Scatter data from rank-0 to all ranks.

symmetrize(rotation_scc, translation_sc)[source]

Make data symmetric.

to_complex()[source]

Return a copy with dtype=complex.

to_pbc_grid()[source]

Convert to UniformGrid with pbc=(True, True, True).

trace_inner_product(other)[source]
xy(*axes)[source]

Extract x, y values along line.

Useful for plotting:

x, y = grid.xy(0, ..., 0)
plt.plot(x, y)
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.

add_ked(weights, out)[source]

Add weighted absolute square of gradient of data to output array.

copy()[source]
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
flat()[source]
gather(out=None, broadcast=False)[source]
gathergather()[source]
integrate(other=None)[source]
interpolate(*, plan1=None, plan2=None, grid=None, out=None)[source]
property matrix: Matrix
matrix_elements(other, *, out=None, symmetric='_default', function=None, domain_sum=True, cc=False)[source]
new(data=None, dims=None)[source]
norm2(kind='normal', weights=None, skip_sum=False)[source]
redist(domain, comm1, comm2)[source]

Redistribute to new domain.

The “world” is spanned by:

(self.desc.comm, comm1)

and:

(domain.comm, comm2).
sanity_check()[source]

Sanity check.

scatter_from(data=None)[source]
to_xp(xp)[source]
trace_inner_product(other)[source]
class gpaw.core.atom_arrays.AtomArrays(layout, dims=(), comm=SerialCommunicator(), data=None)[source]

AtomArrays object.

Parameters:
  • layout (AtomArraysLayout) – Layout-description.

  • dims (int | Sequence[int]) – Extra dimensions.

  • 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_axii and out_ani:

\[\sum^{}_{i} A^{a}_{ni} B^{a}_{ij} \rightarrow C^{a}_{nj}\]

If index is not None, block_diag_matrix_axii must have an extra dimension: \(B_{ij}^{ax}\) and x=index is used.

copy()[source]
gather(broadcast: Literal[False] = False, copy: bool = False) AtomArrays | None[source]
gather(broadcast: Literal[True], copy: bool = False) AtomArrays

Gather all atoms on master.

gathergather()[source]

Gather all atoms and extra dimensions on master.

get(a)[source]
items()[source]
keys()[source]
property matrix: Matrix
moved(atomdist)[source]
my_slice()[source]
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.

redist(atomdist, comm1, comm2)[source]
scatter_everything_from(a_axi=None)[source]

Scatter atoms and extra dimension.

scatter_from(data=None)[source]

Scatter atoms.

to_cpu()[source]
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.]])
to_lower_triangle()[source]

Convert \(N*N\) matrices to \(N*(N+1)/2\) vectors.

>>> a = AtomArraysLayout([(3, 3)]).empty()
>>> a[0][:] = [[11, 12, 13],
...            [12, 22, 23],
...            [13, 23, 33]]
>>> a.to_lower_triangle()[0]
array([11., 12., 22., 13., 23., 33.])
to_xp(xp)[source]
values()[source]
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:
  • dims (int | tuple[int, ...]) – Extra dimensions.

  • comm (_Communicator) – Distribute dimensions along this communicator.

new(atomdist=None, dtype=None, xp=None)[source]

Create new AtomsArrayLayout object with new atomdist.

sizes()[source]

Compute array sizes for all ranks.

>>> AtomArraysLayout([3, 4]).sizes()
([{0: 3, 1: 4}], array([7]))
zeros(dims=(), comm=SerialCommunicator())[source]
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])
classmethod from_number_of_atoms(natoms, comm=SerialCommunicator())[source]

Distribute atoms evenly.

>>> AtomDistribution.from_number_of_atoms(3).rank_a
array([0, 0, 0])
gather()[source]
class gpaw.core.PWArray(pw, dims=(), comm=SerialCommunicator(), data=None, xp=None)[source]

Object for storing function(s) as a plane-wave expansions.

Parameters:
  • pw (PWDesc) – Description of plane-waves.

  • dims (int | tuple[int, ...]) – Extra dimensions.

  • comm (MPIComm) – Distribute extra dimensions along this communicator.

  • data (np.ndarray | None) – Data array for storage.

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}\]
add_ked(occ_n, taut_R)[source]

Add weighted absolute square of gradient of data to output array.

boundary_value(axis)[source]

Calculate average value at boundary of box.

copy()[source]

Create a copy (surprise!).

gather(out=None, broadcast=False)[source]

Gather coefficients on master.

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.

integrate(other=None)[source]

Integral of self or self time cc(other).

interpolate(plan1=None, plan2=None, grid=None, out=None)[source]
property matrix: Matrix

Matrix view of data.

moment()[source]
morph(pw)[source]

Transform expansion to new cell.

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,\]
randomize(seed=None)[source]

Insert random numbers between -0.5 and 0.5 into data.

sanity_check()[source]

Sanity check for real-valued PW expansions.

Make sure the G=(0,0,0) coefficient doesn’t have an imaginary part.

scatter_everything_from(array, comm)[source]

Scatter everything from rank-0 to all ranks.

scatter_from(data=None)[source]

Scatter plane-wave coefficients from rank-0 to all ranks.

scatter_from_all(a_G)[source]

Scatter all coefficients from rank r to self on other cores.

to_pbc_grid()[source]
trace_inner_product(other)[source]
transform(U_cc, complex_conjugate=False, pw=None)[source]

Symmetry-transform data.

class gpaw.core.plane_waves.Empty(dims)[source]
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.

add_hermitian_conjugate(scale=1.0)[source]

Add hermitian conjugate to myself.

add_to_diagonal(d)[source]

Add list of numbers or single number to diagonal of matrix.

complex_conjugate()[source]

Inplace complex conjugation.

copy()[source]

Create a copy.

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 S is given, solve a generalized eigenvalue problem.

Parameters:
  • cc (bool) – Complex conjugate matrix before finding eigenvalues.

  • scalapack (tuple) – BLACS distribution for ScaLapack to use. Default is to do serial diagonalization.

  • limit (int | None) – Number of eigenvector and values to find. Defaults to all.

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.

inv(uplo='L')[source]

Inplace inversion.

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]])
is_distributed()[source]

True if this matrix has nontrivial BLACS or GPU distribution.

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.

static scatter(data, *, dist, root=0)[source]

Construct a distributed Matrix object by scattering a raw 2D array from ‘root’ rank. The ‘dist’ argument must specify the communicator and wanted distribution in same way as in the Matrix constructor Empty ‘dist’ argument is not allowed!

scatter_from(other)[source]

Scatter from rank-0.

to_cpu()[source]

Create new matrix object with values transferred from GPU to CPU.

to_dtype(dtype)[source]

Convert to new data type.

to_xp(xp)[source]

Create new matrix object with data on GPU or CPU.

tril2full()[source]

Fill in upper triangle from lower triangle.

For a real matrix:

a ? ?    a b d
b c ? -> b c e
d e f    d e f

For a complex matrix, the complex conjugate of the lower part will be inserted into the upper part.

class gpaw.core.matrix.MatrixDistribution[source]
add_hermitian_conjugate(a, scale)[source]
all_data_on_rank_zero = True
bc: int
br: int
columns: int
comm: _Communicator
desc: ndarray
eigh_serial(H, S, eigs, cc=False, limit=None)[source]
eighl(H, L)[source]
full_shape: tuple[int, int]
matrix(dtype=None, data=None)[source]
multiply(alpha, a, opa, b, opb, beta, c, symmetric)[source]
my_row_range()[source]

Return indices for range of my rows.

>>> Matrix(2, 2).dist.my_row_range()
(0, 2)
new(M, N)[source]
redist(other, m1data, m2data)[source]
rows: int
shape: tuple[int, int]
simple = True
to_xp(xp)[source]

Input-parameter objects

class gpaw.dft.PW(ecut=340, *, qspiral=None, dedecut=None, dtype='double', force_complex_dtype=False, **kwargs)[source]

PW-mode.

Parameters:

ecut (float) – Plane-wave cutoff energy in eV.

class gpaw.dft.LCAO(*, dtype='double', force_complex_dtype=False, interpolation=117)[source]
class gpaw.dft.FD(*, nn=3, dtype='double', force_complex_dtype=False)[source]
class gpaw.dft.Mode(*, dtype='double', force_complex_dtype=False, interpolation=117)[source]
class gpaw.dft.ExtensionInput[source]
class gpaw.dft.MonkhorstPack(size=None, density=None, gamma=None, even=None)[source]
class gpaw.dft.Mixer(params)[source]
class gpaw.dft.Occupations(params)[source]
class gpaw.dft.PoissonSolver(params)[source]
class gpaw.dft.XC(name, **kwargs)[source]

DFT-components

class gpaw.new.calculation.DFTCalculation(atoms, ibzwfs, density, potential, setups, scf_loop, pot_calc, log, params, energies=None)[source]
ase_calculator()[source]

Create ASE-compatible GPAW calculator.

calculate_dipole()[source]

Calculate dipole moment in Angstrom.

calculate_energy()[source]

Calculate total energy in eV.

calculate_forces()[source]

Calculate atomic forces in eV / Angstrom.

calculate_stress()[source]

Calculate stress in eV / Angstrom^3.

change(*, xc=None, eigensolver=None, mixer=None, occupations=None, convergence=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.

densities()[source]
electrostatic_potential()[source]
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).

get_state()[source]
iconverge(maxiter=None, calculate_forces=None)[source]
magmoms()[source]
move_atoms(atoms)[source]
new(atoms, params, log=None)[source]

Create new DFTCalculation object.

property state
vacuum_level()[source]
wave_function(band, kpt=0, spin=None, periodic=False, broadcast=True)[source]
wave_functions(n1=0, n2=None, kpt=0, spin=None, periodic=False, broadcast=True, _pad=True)[source]
workfunctions(*, _vacuum_level=None)[source]
write_converged()[source]
write_gpw_file(filename, *, include_wfs=False, precision='double', include_projections=True)[source]

Write calculator object to a file.

Parameters:
  • filename – File to be written.

  • include_wfs (bool) – Use include_wfs=True to include wave functions in the file.

  • precision (str) – ‘double’ (the default) or ‘single’.

  • include_projections (bool) – Use include_projections=False to not include the PAW-projections.

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]
calculate_compensation_charge_coefficients()[source]
calculate_dipole_moment(relpos_ac)[source]
calculate_magnetic_moments()[source]
calculate_orbital_magnetic_moments()[source]
classmethod from_data_and_setups(nt_sR, taut_sR, D_asii, charge, setups, 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]
gather()[source]
move(relpos_ac, atomdist)[source]
property nct_R
new(new_grid, pw, relpos_ac=None, atomdist=None)[source]
normalize(background_charge)[source]
redist(grid, xdesc, atomdist, comm1, comm2)[source]
symmetrize(symmetries)[source]
property tauct_R
update(ibzwfs, ked=False)[source]
update_ked(ibzwfs, symmetrize=True)[source]
write_to_gpw(writer, flags)[source]
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.

add_to_density(nt_sR, D_asii)[source]

Compute density and add to nt_sR and D_asii.

add_to_ked(taut_sR)[source]
calculate_kinetic_energy(hamiltonian, density)[source]
calculate_occs(occ_calc, nelectrons, fix_fermi_level=False)[source]
classmethod create(*, ibz, ncomponents, create_wfs_func, kpt_comm=SerialCommunicator(), kpt_band_comm=SerialCommunicator(), comm=SerialCommunicator())[source]
property fermi_level: float
forces(potential, hamiltonian, D_asii)[source]
gather()[source]
get_all_eigs_and_occs(broadcast=False)[source]
get_all_electron_wave_function(band, kpt=0, spin=0, grid_spacing=0.05, skip_paw_correction=False)[source]
get_eigs_and_occs(kpt=0, spin=0)[source]
get_homo_lumo(spin=None)[source]

Return HOMO and LUMO eigenvalues.

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.

get_wfs_on_master(*, kpt=0, spin=0, n1=0, n2=0)[source]
has_wave_functions()[source]
make_sure_wfs_are_read_from_gpw_file()[source]
property mode
move(relpos_ac, atomdist)[source]
normalize_density(density)[source]
number_of_occupied_bands(tolerance=1e-05)[source]
orthonormalize(work_array_nX=None)[source]
summary(log)[source]
write(writer, flags)[source]

Write fermi-level(s), eigenvalues, occupation numbers, …

… k-points, symmetry information, projections and possibly also the wave functions.

write_summary(log)[source]
zero_padded_iter()[source]

Load-balancing for kpt_comm when doing hybrid functionals.

class gpaw.new.potential.Potential(vt_sR, dH_asii, dedtaut_sR, vHt_x=None, e_stress=nan)[source]
copy()[source]

Make a deep copy of the potential, preserving the parallel distribution of data.

Useful, e.g. in RTTDDFT, where the potential is copied in order to be restored later.

deltaH(P_ani, out_ani, spin)[source]
gather()[source]
get_vacuum_level()[source]
move(atomdist)[source]

Move atoms inplace.

redist(grid, desc, atomdist, comm1, comm2)[source]
write_to_gpw(writer, flags)[source]
class gpaw.new.pot_calc.PotentialCalculator(xc, poisson_solver, setups, *, relpos_ac, extensions=None, soc=False)[source]
calculate(density, ibzwfs=None, vHt_x=None, kpt_band_comm=None)[source]
calculate_charges(vHt_x)[source]
calculate_pseudo_potential(density, ibzwfs, vHt_x)[source]
calculate_without_orbitals(density, ibzwfs=None, vHt_x=None, kpt_band_comm=None)[source]
property charge
move(relpos_ac, atomdist)[source]
restrict(a_r, a_R=None)[source]
class gpaw.new.scf.SCFLoop(hamiltonian, occ_calc, eigensolver, mixer, comm, convergence, maxiter)[source]
iterate(ibzwfs, density, potential, energies, pot_calc, *, maxiter=None, calculate_forces=None, log=None)[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:
  • mode (str | dict | Mode) – PW, LCAO or FD mode.

  • 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.

  • mixer (dict | Mixer | None) – Density-mixing scheme.

  • nbands (int | str | None) – Number of bands.

  • occupations (dict | Occupations | None) –

  • parallel (dict | None) – Parallelization strategy. Example: Force parallelization over 'kpt with {'band': 1, 'domain': 1}.

  • poissonsolver (dict | PoissonSolver | None) –

  • random (bool | None) – Use random numbers for initial wave functions.

  • setups (str | dict | None) –

  • 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.

dft_calculation(atoms, txt='-', communicator=None)[source]
dft_component_builder(atoms, *, comm=None, log=None)[source]
dft_info(atoms)[source]
todict()[source]
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\]
add_to_density(nt_sR, D_asii)[source]
add_to_ked(taut_sR)[source]
array_shape(global_shape=False)[source]
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\]
canonical_transformation(H_nm, psit2_nX, data_buffer)[source]
collect_bands(n1=0, n2=0)[source]

Collect range of bands to master of band-comm.

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

force_contribution(potential, F_av)[source]
classmethod from_wfs(wfs, psit_nX, relpos_ac=None, atomdist=None, domain_band_comm=None)[source]
gather_wave_function_coefficients()[source]
morph(desc, relpos_ac, atomdist)[source]
move(relpos_ac, atomdist, move_wave_functions)[source]
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})\]
receive(rank, comm)[source]
send(rank, comm)[source]
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\]
subspace_eigenvalues(H_nn, scalapack_params=(None, 1, 1, 0))[source]
to_uniform_grid_wave_functions(grid, basis, *, xp=None)[source]
class gpaw.new.ase_interface.ASECalculator(params, *, log, dft=None, atoms=None)[source]

This is the ASE-calculator frontend for doing a GPAW calculation.

property atoms: Atoms
band_structure()[source]

Create band-structure object for plotting.

calculate(atoms, properties=None, system_changes=None)[source]
calculate_property(atoms, prop)[source]

Calculate (if not already calculated) a property.

The prop string must be one of

  • energy

  • forces

  • stress

  • magmom

  • magmoms

  • dipole

calculation_required(atoms, properties)[source]
check_state(atoms, tol=1e-12)[source]
converge_wave_functions()[source]
create_new_calculation(atoms)[source]
create_new_calculation_from_old(atoms)[source]
property density
property dft: DFTCalculation
diagonalize_full_hamiltonian(nbands=None, scalapack=None, expert=None)[source]
dos(soc=False, theta=0.0, phi=0.0, shift_fermi_level=True)[source]

Create DOS-calculator.

Default is to shift_fermi_level to 0.0 eV. For soc=True, angles can be given in degrees.

eigenvalues()[source]
fixed_density(*, txt='-', update_fermi_level=False, **kwargs)[source]
get_all_electron_density(spin=None, gridrefinement=1, broadcast=True, skip_core=False)[source]
get_atomic_electrostatic_potentials()[source]
get_atoms()[source]
get_bz_k_points()[source]
get_bz_to_ibz_map()[source]

Return indices from BZ to IBZ.

get_dipole_moment(atoms=None)[source]
get_effective_potential(spin=0, broadcast=True)[source]
get_eigenvalues(kpt=0, spin=0, broadcast=True)[source]
get_electrostatic_corrections()[source]
get_electrostatic_potential()[source]
get_fermi_level()[source]
get_fermi_levels()[source]
get_forces(atoms=None)[source]
get_homo_lumo(spin=None)[source]
get_ibz_k_points()[source]
get_k_point_weights()[source]
get_magnetic_moment(atoms=None)[source]
get_magnetic_moments(atoms=None)[source]
get_non_collinear_magnetic_moment(atoms=None)[source]
get_non_collinear_magnetic_moments(atoms=None)[source]
get_nonselfconsistent_energies(type='beefvdw')[source]
get_number_of_bands()[source]
get_number_of_electrons()[source]
get_number_of_grid_points()[source]
get_number_of_iterations()[source]
get_number_of_spins()[source]
get_occupation_numbers(kpt=0, spin=0, broadcast=True, raw=False)[source]
get_orbital_magnetic_moments()[source]

Return the orbital magnetic moment vector for each atom.

get_potential_energy(atoms=None, force_consistent=False)[source]
get_property(name, atoms=None, allow_calculation=True)[source]
get_pseudo_density(spin=None, gridrefinement=1, broadcast=True)[source]
get_pseudo_wave_function(band, kpt=0, spin=None, periodic=False, broadcast=True, pad=True)[source]
get_reference_energy()[source]
get_stress(atoms=None)[source]
get_wannier_localization_matrix(nbands, dirG, kpoint, nextkpoint, G_I, spin)[source]

Calculate integrals for maximally localized Wannier functions.

get_xc_difference(xcparams)[source]

Calculate non-selfconsistent XC-energy difference.

get_xc_functional()[source]
gs_adapter()[source]
property hamiltonian
icalculate(atoms, system_changes=None)[source]
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']
initial_wannier(initialwannier, kpointgrid, fixedstates, edf, spin, nbands)[source]
initialize(atoms)[source]
initialize_positions(atoms=None)[source]
property initialized
move_atoms(atoms)[source]
name = 'gpaw'
new(**kwargs)[source]
occupations()[source]
old = False
property parameters
property results
set(eigensolver)[source]
property setups
property spos_ac
property symmetry
todict()[source]
property wfs
property world
write(filename, mode='', precision='double', include_projections=True)[source]

Write calculator object to a file.

Parameters:
  • filename (str | Path) – File to be written.

  • mode (str) – Write mode. Use mode='all' to include wave functions in the file.

  • precision (str) – ‘double’ (the default) or ‘single’.

  • include_projections (bool) – Use include_projections=False to not include the PAW-projections.

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.Parameters for the complete list of parameters.

Parameters:
  • atoms (Atoms) – ASE-Atoms object.

  • txt (str | Path | IO[str] | None) – Text log-file. Use None for no logging and '-' for using standard out.

  • communicator (_Communicator | None) – MPI-communicator. Default is to use gpaw.mpi.world.

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.Parameters for 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 None for 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
analyze_positions(relative_positions, ids, *, symmorphic=True)[source]
check_grid(N_c)[source]

Check that symmetries are commensurate with grid.

check_positions(fracpos_ac)[source]
classmethod from_atoms(atoms, *, ids=None, symmorphic=True, tolerance=None)[source]
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
group_check()[source]

Sanity check.

property has_inversion
lcm()[source]

Find lowest common multiple compatible with translations.

summary(log, verbose=True)[source]
symmetrize_forces(F0_av)[source]

Symmetrize forces.

property symmorphic
with_atom_maps(relpos_ac, id_a)[source]
class gpaw.new.brillouin.IBZ(symmetries, bz, ibz2bz, bz2ibz, weights, bz2bz_Ks=None, s_K=None, time_reversal_K=None)[source]
ranks(comm, nspins=1)[source]

Distribute k-points over MPI-communicator.

summary(log, verbose=True)[source]
class gpaw.new.brillouin.BZPoints(points)[source]
reduce(symmetries, *, comm=None, strict=True, use_time_reversal=True, tolerance=1e-07)[source]

Find irreducible set of k-points.

class gpaw.new.brillouin.MonkhorstPackKPoints(size, shift=(0, 0, 0))[source]
class gpaw.new.extensions.Extension[source]

Extension object.

Used for jellium, solvation, solvated jellium model, D3, …

build(builder)[source]
charge = 0.0
create_poisson_solver(grid, pw, *, charge, xp)[source]
force_contribution(nt_r, vHt_r)[source]
get_energy_contributions()[source]
move_atoms(relpos_ac)[source]
name = 'unnamed extension'
post_scf_convergence(ibzwfs, nelectrons, occ_calc, mixer, log)[source]

Allow for environment to “converge”.

stress_contribution()[source]
update1(nt_r)[source]

Hook called right before solving the Poisson equation.

update1pw(nt_g)[source]

PW-mode hook called right before solving the Poisson equation.

update2(nt_r, vHt_r, vt_sr)[source]

Calculate environment energy.

update_non_local_hamiltonian(D_sii, setup, atom_index, dH_sii)[source]
update_potential(vt_sR, density)[source]

FFTW

Python wrapper for FFTW3 library

class gpaw.fftw.FFTPlans(size_c, dtype, empty=<function empty>)[source]
class gpaw.fftw.FFTPlan(in_R, out_R, sign, flags=0)[source]

FFT 3d transform.

class gpaw.fftw.FFTWPlan(in_R, out_R, sign, flags=0)[source]

FFTW3 3d transform.

class gpaw.fftw.NumpyFFTPlan(in_R, out_R, sign, flags=0)[source]

Numpy fallback.

class gpaw.fftw.NumpyFFTPlans(size_c, dtype, empty=<function empty>)[source]

Numpy fallback.

fft()[source]

Do FFT from tmp_R to tmp_Q.

>>> plans = create_plans([4, 1, 1], float)
>>> plans.tmp_R[:, 0, 0] = [1, 0, 1, 0]
>>> plans.fft()
>>> plans.tmp_Q[:, 0, 0]
array([2.+0.j, 0.+0.j, 2.+0.j, 0.+0.j])
ifft()[source]

Do inverse FFT from tmp_Q to tmp_R.

>>> plans = create_plans([4, 1, 1], complex)
>>> plans.tmp_Q[:, 0, 0] = [0, 1j, 0, 0]
>>> plans.ifft()
>>> plans.tmp_R[:, 0, 0]
array([ 0.+1.j, -1.+0.j,  0.-1.j,  1.+0.j])
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.

gpaw.fftw.get_efficient_fft_size(N, n=1, factors=[2, 3, 5, 7])[source]

Return smallest efficient fft size.

Must be greater than or equal to N and divisible by n.

>>> get_efficient_fft_size(17)
18
gpaw.fftw.have_fftw()[source]

Did we compile with FFTW?

BLAS

gpaw.utilities.blas.mmm(alpha, a, opa, b, opb, beta, c)[source]
gpaw.utilities.blas.rk(alpha, a, beta, c, trans='c')[source]

NOTE: Fills in the entire c matrix, unlike the C-blas version that only fills in the lower triangle.