Skip to content

Unit Cell

UnitCell

Storage class for the lattice vectors of a crystal i.e. its unit cell.

Attributes:

Name Type Description
direct np.ndarray

the direct matrix of this unit cell i.e. the lattice vectors

reciprocal_lattice np.ndarray

the reciprocal matrix of this unit cell i.e. the reciprocal lattice vectors

inverse np.ndarray

the inverse matrix of this unit cell i.e. the transpose of reciprocal_lattice

lattice np.ndarray

an alias for direct

a: float property readonly

Length of lattice vector a

a_star: float property readonly

length of reciprocal lattice vector a*

abc_different: bool property readonly

are all of the lengths a, b, c different?

abc_equal: bool property readonly

are the lengths a, b, c all equal?

alpha: float property readonly

Angle between lattice vectors b and c

alpha_deg: float property readonly

Angle between lattice vectors b and c in degrees

alpha_star: float property readonly

Angle between reciprocal lattice vectors b and c

angles_different: bool property readonly

are all of the angles alpha, beta, gamma different?

b: float property readonly

Length of lattice vector b

b_star: float property readonly

length of reciprocal lattice vector b*

beta: float property readonly

Angle between lattice vectors a and c

beta_deg: float property readonly

Angle between lattice vectors a and c in degrees

beta_star: float property readonly

Angle between reciprocal lattice vectors a and c

c: float property readonly

Length of lattice vector c

c_star: float property readonly

length of reciprocal lattice vector c*

direct_homogeneous: ndarray property readonly

The direct matrix in homogeneous coordinates

gamma: float property readonly

Angle between lattice vectors a and b

gamma_deg: float property readonly

Angle between lattice vectors a and b in degrees

gamma_star: float property readonly

Angle between reciprocal lattice vectors a and c

is_cubic: bool property readonly

Returns true if all lengths are equal and all angles are 90 degrees

is_hexagonal: bool property readonly

Returns true if lengths a == b, a != c, alpha and beta == 90 and gamma == 120

is_monoclinic: bool property readonly

Returns true if angles alpha and gamma are equal

is_orthorhombic: bool property readonly

Returns true if all angles are 90 degrees

is_rhombohedral: bool property readonly

Returns true if all lengths are equal and all angles are equal

is_tetragonal: bool property readonly

Returns true if a, b are equal and all angles are 90 degrees

is_triclinic: bool property readonly

Returns true if and lengths are different

lattice: ndarray property readonly

The direct matrix of this unit cell i.e. vectors of the lattice

orthogonal: bool property readonly

returns true if the lattice vectors are orthogonal

parameters: ndarray property readonly

single vector of lattice side lengths and angles in degrees

reciprocal_lattice: ndarray property readonly

The reciprocal matrix of this unit cell i.e. vectors of the reciprocal lattice

v_a: ndarray property readonly

lattice vector a

v_a_star: ndarray property readonly

reciprocal lattice vector a*

v_b: ndarray property readonly

lattice vector a

v_b_star: ndarray property readonly

reciprocal lattice vector b*

v_c: ndarray property readonly

lattice vector a

v_c_star: ndarray property readonly

reciprocal lattice vector c*

__init__(self, vectors) special

Create a UnitCell object from a list of lattice vectors or a row major direct matrix. Unless otherwise specified, length units are Angstroms, and angular units are radians.

Parameters:

Name Type Description Default
vectors array_like

(3, 3) array of lattice vectors, row major i.e. vectors[0, :] is lattice vector A etc.

required
Source code in chmpy/crystal/unit_cell.py
22
23
24
25
26
27
28
29
30
31
32
def __init__(self, vectors):
    """
    Create a UnitCell object from a list of lattice vectors or
    a row major direct matrix. Unless otherwise specified, length
    units are Angstroms, and angular units are radians.

    Args:
        vectors (array_like): (3, 3) array of lattice vectors, row major i.e. vectors[0, :] is
            lattice vector A etc.
    """
    self.set_vectors(vectors)

cubic(length) classmethod

Construct a new cubic UnitCell from the provided side length.

Parameters:

Name Type Description Default
length float

Lattice side length a in Angstroms.

required

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
424
425
426
427
428
429
430
431
432
433
434
435
@classmethod
def cubic(cls, length):
    """
    Construct a new cubic UnitCell from the provided side length.

    Args:
        length (float): Lattice side length a in Angstroms.

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """
    return cls(np.eye(3) * length)

from_lengths_and_angles(lengths, angles, unit='radians') classmethod

Construct a new UnitCell from the provided lengths and angles.

Parameters:

Name Type Description Default
lengths array_like

Lattice side lengths (a, b, c) in Angstroms.

required
angles array_like

Lattice angles (alpha, beta, gamma) in provided units (default radians)

required
unit str

Unit for angles i.e. 'radians' or 'degrees' (default radians).

'radians'

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
@classmethod
def from_lengths_and_angles(cls, lengths, angles, unit="radians"):
    """
    Construct a new UnitCell from the provided lengths and angles.

    Args:
        lengths (array_like): Lattice side lengths (a, b, c) in Angstroms.
        angles (array_like): Lattice angles (alpha, beta, gamma) in provided units (default radians)
        unit (str, optional): Unit for angles i.e. 'radians' or 'degrees' (default radians).

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """
    uc = cls(np.eye(3))
    if unit == "radians":
        if np.any(np.abs(angles) > np.pi):
            LOG.warn(
                "Large angle in UnitCell.from_lengths_and_angles, "
                "are you sure your angles are not in degrees?"
            )
        uc.set_lengths_and_angles(lengths, angles)
    else:
        uc.set_lengths_and_angles(lengths, np.radians(angles))
    return uc

from_unique_parameters(params, cell_type='triclinic', **kwargs) classmethod

Construct a new unit cell from the unique parameters and the specified cell type.

Parameters:

Name Type Description Default
params Tuple

tuple of floats of unique parameters

required
cell_type str

the desired cell type

'triclinic'
Source code in chmpy/crystal/unit_cell.py
437
438
439
440
441
442
443
444
445
446
447
@classmethod
def from_unique_parameters(cls, params, cell_type="triclinic", **kwargs):
    """
    Construct a new unit cell from the unique parameters and
    the specified cell type.

    Args:
        params (Tuple): tuple of floats of unique parameters
        cell_type (str, optional): the desired cell type
    """
    return getattr(cls, cell_type)(*params)

hexagonal(*params, **kwargs) classmethod

Construct a new UnitCell from the provided side lengths and angles.

Parameters:

Name Type Description Default
params array_like

Lattice side lengths (a, c)

()

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
@classmethod
def hexagonal(cls, *params, **kwargs):
    """
    Construct a new UnitCell from the provided side lengths and angles.

    Args:
        params (array_like): Lattice side lengths (a, c)

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """
    assert len(params) == 2, "Requre 2 lengths for Hexagonal cell"
    unit = kwargs.pop("unit", "radians")
    unit = "radians"
    angles = [np.pi / 2, np.pi / 2, 2 * np.pi / 3]
    return cls.from_lengths_and_angles(
        (params[0], params[0], params[1]), angles, unit=unit, **kwargs
    )

monoclinic(*params, **kwargs) classmethod

Construct a new UnitCell from the provided side lengths and angle.

Parameters:

Name Type Description Default
params array_like

Lattice side lengths and angles (a, b, c, beta)

()

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
@classmethod
def monoclinic(cls, *params, **kwargs):
    """
    Construct a new UnitCell from the provided side lengths and angle.

    Args:
        params (array_like): Lattice side lengths and angles (a, b, c, beta)

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """

    assert (
        len(params) == 4
    ), "Requre three lengths and one angle for Monoclinic cell"
    unit = kwargs.get("unit", "radians")
    if unit != "radians":
        alpha, gamma = 90, 90
    else:
        alpha, gamma = np.pi / 2, np.pi / 2
    return cls.from_lengths_and_angles(
        params[:3], (alpha, params[3], gamma), **kwargs
    )

orthorhombic(*lengths, **kwargs) classmethod

Construct a new orthorhombic UnitCell from the provided side lengths.

Parameters:

Name Type Description Default
lengths array_like

Lattice side lengths (a, b, c) in Angstroms.

()

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
542
543
544
545
546
547
548
549
550
551
552
553
554
555
@classmethod
def orthorhombic(cls, *lengths, **kwargs):
    """
    Construct a new orthorhombic UnitCell from the provided side lengths.

    Args:
        lengths (array_like): Lattice side lengths (a, b, c) in Angstroms.

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """

    assert len(lengths) == 3, "Requre three lengths for Orthorhombic cell"
    return cls(np.diag(lengths))

rhombohedral(*params, **kwargs) classmethod

Construct a new UnitCell from the provided side lengths and angles.

Parameters:

Name Type Description Default
params array_like

Lattice side length a and angle alpha c

()

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
528
529
530
531
532
533
534
535
536
537
538
539
540
@classmethod
def rhombohedral(cls, *params, **kwargs):
    """
    Construct a new UnitCell from the provided side lengths and angles.

    Args:
        params (array_like): Lattice side length a and angle alpha c

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """
    assert len(params) == 2, "Requre 1 length and 1 angle for Rhombohedral cell"
    return cls.from_lengths_and_angles([params[0]] * 3, [params[1]] * 3, **kwargs)

set_lengths_and_angles(self, lengths, angles)

Modify this unit cell by setting the lattice vectors according to lengths a, b, c and angles alpha, beta, gamma of a parallelipiped.

Parameters:

Name Type Description Default
lengths array_like

array of (a, b, c), the unit cell side lengths in Angstroms.

required
angles array_like

array of (alpha, beta, gamma), the unit cell angles lengths in radians.

required
Source code in chmpy/crystal/unit_cell.py
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
def set_lengths_and_angles(self, lengths, angles):
    """
    Modify this unit cell by setting the lattice vectors
    according to lengths a, b, c and angles alpha, beta, gamma of
    a parallelipiped.

    Args:
        lengths (array_like): array of (a, b, c), the unit cell side lengths in Angstroms.
        angles (array_like): array of (alpha, beta, gamma), the unit cell angles lengths
            in radians.
    """
    self.lengths = lengths
    self.angles = angles
    a, b, c = self.lengths
    ca, cb, cg = np.cos(self.angles)
    sg = np.sin(self.angles[2])
    v = self.volume()
    self.direct = np.transpose(
        [
            [a, b * cg, c * cb],
            [0, b * sg, c * (ca - cb * cg) / sg],
            [0, 0, v / (a * b * sg)],
        ]
    )
    r = [
        [1 / a, 0.0, 0.0],
        [-cg / (a * sg), 1 / (b * sg), 0],
        [
            b * c * (ca * cg - cb) / v / sg,
            a * c * (cb * cg - ca) / v / sg,
            a * b * sg / v,
        ],
    ]
    self.inverse = np.array(r)
    self._set_cell_type()

set_vectors(self, vectors)

Modify this unit cell by setting the lattice vectors according to those provided. This is performed by setting the lattice parameters (lengths and angles) based on the provided vectors, such that it results in a consistent basis without directly matrix inverse (and typically losing precision), and as the SHELX file/CIF output will be relying on these lengths/angles anyway, it is important to have these consistent.

Parameters:

Name Type Description Default
vectors array_like

(3, 3) array of lattice vectors, row major i.e. vectors[0, :] is lattice vector A etc.

required
Source code in chmpy/crystal/unit_cell.py
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
def set_vectors(self, vectors):
    """
    Modify this unit cell by setting the lattice vectors
    according to those provided. This is performed by setting the
    lattice parameters (lengths and angles) based on the provided vectors,
    such that it results in a consistent basis without directly
    matrix inverse (and typically losing precision), and
    as the SHELX file/CIF output will be relying on these
    lengths/angles anyway, it is important to have these consistent.


    Args:
        vectors (array_like): (3, 3) array of lattice vectors, row major i.e. vectors[0, :] is
            lattice vector A etc.
    """
    self.direct = vectors
    params = zeros(6)
    a, b, c = np.linalg.norm(self.direct, axis=1)
    u_a = vectors[0, :] / a
    u_b = vectors[1, :] / b
    u_c = vectors[2, :] / c
    alpha = np.arccos(np.clip(np.vdot(u_b, u_c), -1, 1))
    beta = np.arccos(np.clip(np.vdot(u_c, u_a), -1, 1))
    gamma = np.arccos(np.clip(np.vdot(u_a, u_b), -1, 1))
    params[3:] = np.degrees([alpha, beta, gamma])
    self.lengths = [a, b, c]
    self.angles = [alpha, beta, gamma]
    self.inverse = np.linalg.inv(self.direct)
    self._set_cell_type()

tetragonal(*params, **kwargs) classmethod

Construct a new UnitCell from the provided side lengths and angles.

Parameters:

Name Type Description Default
params array_like

Lattice side lengths (a, c)

()

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
@classmethod
def tetragonal(cls, *params, **kwargs):
    """
    Construct a new UnitCell from the provided side lengths and angles.

    Args:
        params (array_like): Lattice side lengths (a, c)

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """
    assert len(params) == 2, "Requre 2 lengths for Tetragonal cell"
    unit = kwargs.get("unit", "radians")
    if unit != "radians":
        angles = [90] * 3
    else:
        angles = [np.pi / 2] * 3
    return cls.from_lengths_and_angles(
        (params[0], params[0], params[1]), angles, **kwargs
    )

to_cartesian(self, coords)

Transform coordinates from fractional space (a, b, c) to Cartesian space (x, y, z). The x-direction will be aligned along lattice vector A.

Parameters:

Name Type Description Default
coords ndarray

(N, 3) array of fractional coordinates

required

Returns:

Type Description
ndarray

np.ndarray: (N, 3) array of Cartesian coordinates

Source code in chmpy/crystal/unit_cell.py
51
52
53
54
55
56
57
58
59
60
61
62
63
def to_cartesian(self, coords: np.ndarray) -> np.ndarray:
    """
    Transform coordinates from fractional space (a, b, c)
    to Cartesian space (x, y, z). The x-direction will be aligned
    along lattice vector A.

    Args:
        coords (array_like): (N, 3) array of fractional coordinates

    Returns:
        np.ndarray: (N, 3) array of Cartesian coordinates
    """
    return np.dot(coords, self.direct)

to_fractional(self, coords)

Transform coordinates from Cartesian space (x, y, z) to fractional space (a, b, c). The x-direction will is assumed be aligned along lattice vector A.

Parameters:

Name Type Description Default
coords ndarray

an (N, 3) array of Cartesian coordinates

required

Returns:

Type Description
ndarray

np.ndarray: (N, 3) array of fractional coordinates

Source code in chmpy/crystal/unit_cell.py
65
66
67
68
69
70
71
72
73
74
75
76
77
def to_fractional(self, coords: np.ndarray) -> np.ndarray:
    """
    Transform coordinates from Cartesian space (x, y, z)
    to fractional space (a, b, c). The x-direction will is assumed
    be aligned along lattice vector A.

    Args:
        coords (array_like): an (N, 3) array of Cartesian coordinates

    Returns:
        np.ndarray: (N, 3) array of fractional coordinates
    """
    return np.dot(coords, self.inverse)

triclinic(*params, **kwargs) classmethod

Construct a new UnitCell from the provided side lengths and angles.

Parameters:

Name Type Description Default
params array_like

Lattice side lengths and angles (a, b, c, alpha, beta, gamma)

()

Returns:

Type Description
UnitCell

A new unit cell object representing the provided lattice.

Source code in chmpy/crystal/unit_cell.py
449
450
451
452
453
454
455
456
457
458
459
460
461
462
@classmethod
def triclinic(cls, *params, **kwargs):
    """
    Construct a new UnitCell from the provided side lengths and angles.

    Args:
        params (array_like): Lattice side lengths and angles (a, b, c, alpha, beta, gamma)

    Returns:
        UnitCell: A new unit cell object representing the provided lattice.
    """

    assert len(params) == 6, "Requre three lengths and angles for Triclinic cell"
    return cls.from_lengths_and_angles(params[:3], params[3:], **kwargs)

volume(self)

The volume of the unit cell, in cubic Angstroms

Source code in chmpy/crystal/unit_cell.py
196
197
198
199
200
def volume(self) -> float:
    """The volume of the unit cell, in cubic Angstroms"""
    a, b, c = self.lengths
    ca, cb, cg = np.cos(self.angles)
    return a * b * c * np.sqrt(1 - ca * ca - cb * cb - cg * cg + 2 * ca * cb * cg)