quimb.tensor.tensor_3d#

Functions

calc_cell_map(cells)

calc_cell_sizes(coo_groups[, autogroup])

cell_to_sites(p)

Turn a cell ((i0, j0, k0), (di, dj, dk)) into the sites it contains.

gen_3d_bonds(Lx, Ly, Lz, steppers[, coo_filter])

Convenience function for tiling pairs of bond coordinates on a 3D lattice given a function like lambda i, j, k: (i + 1, j + 1, k + 1).

is_lone_coo(where)

Check if where has been specified as a single coordinate triplet.

sites_to_cell(sites)

Get the minimum covering cell for sites.

Classes

PEPS3D(arrays, *[, shape, tags, ...])

Projected Entangled Pair States object (3D).

Rotator3D(tn, xrange, yrange, zrange, from_which)

Object for rotating coordinates and various contraction functions so that the core algorithms only have to written once, but nor does the actual TN have to be modified.

TensorNetwork3D(ts, *[, virtual, ...])

TensorNetwork3DFlat(ts, *[, virtual, ...])

Mixin class for a 3D square lattice tensor network with a single tensor per site, for example, both PEPS and PEPOs.

TensorNetwork3DVector(ts, *[, virtual, ...])

Mixin class for a 3D square lattice vector TN, i.e. one with a single physical index per site.

class quimb.tensor.tensor_3d.PEPS3D(arrays, *, shape='urfdlbp', tags=None, site_ind_id='k{},{},{}', site_tag_id='I{},{},{}', x_tag_id='X{}', y_tag_id='Y{}', z_tag_id='Z{}', **tn_opts)[source]#

Projected Entangled Pair States object (3D).

Parameters
  • arrays (sequence of sequence of sequence of array) – The core tensor data arrays.

  • shape (str, optional) – Which order the dimensions of the arrays are stored in, the default 'urfdlbp' stands for (‘up’, ‘right’, ‘front’, down’, ‘left’, ‘behind’, ‘physical’) meaning (x+, y+, z+, x-, y-, z-, physical) respectively. Arrays on the edge of lattice are assumed to be missing the corresponding dimension.

  • tags (set[str], optional) – Extra global tags to add to the tensor network.

  • site_ind_id (str, optional) – String specifier for naming convention of site indices.

  • site_tag_id (str, optional) – String specifier for naming convention of site tags.

  • x_tag_id (str, optional) – String specifier for naming convention of x-slice tags.

  • y_tag_id (str, optional) – String specifier for naming convention of y-slice tags.

  • z_tag_id (str, optional) – String specifier for naming convention of z-slice tags.

classmethod empty(Lx, Ly, Lz, bond_dim, phys_dim=2, like='numpy', **peps3d_opts)[source]#

Create an empty 3D PEPS.

Parameters
  • Lx (int) – The number of x-slices.

  • Ly (int) – The number of y-slices.

  • Lz (int) – The number of z-slices.

  • bond_dim (int) – The bond dimension.

  • physical (int, optional) – The physical index dimension.

  • peps3d_opts – Supplied to PEPS3D.

Returns

psi

Return type

PEPS3D

classmethod from_fill_fn(fill_fn, Lx, Ly, Lz, bond_dim, phys_dim=2, **peps3d_opts)[source]#

Create a 3D PEPS from a filling function with signature fill_fn(shape).

Parameters
  • Lx (int) – The number of x-slices.

  • Ly (int) – The number of y-slices.

  • Lz (int) – The number of z-slices.

  • bond_dim (int) – The bond dimension.

  • physical (int, optional) – The physical index dimension.

  • peps_opts – Supplied to PEPS3D.

Returns

psi

Return type

PEPS3D

classmethod ones(Lx, Ly, Lz, bond_dim, phys_dim=2, like='numpy', **peps3d_opts)[source]#

Create a 3D PEPS whose tensors are filled with ones.

Parameters
  • Lx (int) – The number of x-slices.

  • Ly (int) – The number of y-slices.

  • Lz (int) – The number of z-slices.

  • bond_dim (int) – The bond dimension.

  • physical (int, optional) – The physical index dimension.

  • peps3d_opts – Supplied to PEPS3D.

Returns

psi

Return type

PEPS3D

permute_arrays(shape='urfdlbp')[source]#

Permute the indices of each tensor in this PEPS3D to match shape. This doesn’t change how the overall object interacts with other tensor networks but may be useful for extracting the underlying arrays consistently. This is an inplace operation.

Parameters

shape (str, optional) – A permutation of 'lrp' specifying the desired order of the left, right, and physical indices respectively.

classmethod rand(Lx, Ly, Lz, bond_dim, phys_dim=2, dtype='float64', seed=None, **peps3d_opts)[source]#

Create a random (un-normalized) 3D PEPS.

Parameters
  • Lx (int) – The number of x-slices.

  • Ly (int) – The number of y-slices.

  • Lz (int) – The number of z-slices.

  • bond_dim (int) – The bond dimension.

  • physical (int, optional) – The physical index dimension.

  • dtype (dtype, optional) – The dtype to create the arrays with, default is real double.

  • seed (int, optional) – A random seed.

  • peps_opts – Supplied to PEPS3D.

Returns

psi

Return type

PEPS3D

class quimb.tensor.tensor_3d.Rotator3D(tn, xrange, yrange, zrange, from_which)[source]#

Object for rotating coordinates and various contraction functions so that the core algorithms only have to written once, but nor does the actual TN have to be modified.

class quimb.tensor.tensor_3d.TensorNetwork3D(ts, *, virtual=False, check_collisions=True)[source]#
property Lx#

The number of x-slices.

property Ly#

The number of y-slices.

property Lz#

The number of z-slices.

canonize_plane(xrange, yrange, zrange, equalize_norms=False, canonize_opts=None, **gen_pair_opts)[source]#

Canonize every pair of tensors within a subrange, optionally specifying a order to visit those pairs in.

compress_plane(xrange, yrange, zrange, max_bond=None, cutoff=1e-10, equalize_norms=False, compress_opts=None, **gen_pair_opts)[source]#

Compress every pair of tensors within a subrange, optionally specifying a order to visit those pairs in.

contract_boundary(max_bond=None, *, cutoff=1e-10, canonize=True, max_separation=1, sequence=None, xmin=None, xmax=None, ymin=None, ymax=None, zmin=None, zmax=None, optimize='auto-hq', equalize_norms=False, compress_opts=None, inplace=False, **contract_boundary_opts)[source]#
flatten(fuse_multibonds=True, inplace=False)[source]#

Contract all tensors corresponding to each site into one.

gen_bond_coos()[source]#

Generate pairs of coordinates for all the bonds in this 3D TN.

gen_pairs(xrange=None, yrange=None, zrange=None, xreverse=False, yreverse=False, zreverse=False, coordinate_order='xyz', xstep=None, ystep=None, zstep=None, stepping_order='xyz', step_only=None)[source]#

Helper function for generating pairs of cooordinates for all bonds within a certain range, optionally specifying an order.

Parameters
  • xrange ((int, int), optional) – The range of allowed values for the x, y, and z coordinates.

  • yrange ((int, int), optional) – The range of allowed values for the x, y, and z coordinates.

  • zrange ((int, int), optional) – The range of allowed values for the x, y, and z coordinates.

  • xreverse (bool, optional) – Whether to reverse the order of the x, y, and z sweeps.

  • yreverse (bool, optional) – Whether to reverse the order of the x, y, and z sweeps.

  • zreverse (bool, optional) – Whether to reverse the order of the x, y, and z sweeps.

  • coordinate_order (str, optional) – The order in which to sweep the x, y, and z coordinates. Earlier dimensions will change slower. If the corresponding range has size 1 then that dimension doesn’t need to be specified.

  • xstep (int, optional) – When generating a bond, step in this direction to yield the neighboring coordinate. By default, these follow xreverse, yreverse, and zreverse respectively.

  • ystep (int, optional) – When generating a bond, step in this direction to yield the neighboring coordinate. By default, these follow xreverse, yreverse, and zreverse respectively.

  • zstep (int, optional) – When generating a bond, step in this direction to yield the neighboring coordinate. By default, these follow xreverse, yreverse, and zreverse respectively.

  • stepping_order (str, optional) – The order in which to step the x, y, and z coordinates to generate bonds. Does not need to include all dimensions.

  • step_only (int, optional) – Only perform the ith steps in stepping_order, used to interleave canonizing and compressing for example.

Yields

coo_a, coo_b (((int, int, int), (int, int, int)))

gen_site_coos()[source]#

Generate coordinates for all the sites in this 3D TN.

gen_site_coos_present()[source]#

Return only the sites that are present in this TN.

Examples

>>> tn = qtn.TN3D_rand(4, 4, 4, 2)
>>> tn_sub = tn.select_local('I1,2,3', max_distance=1)
>>> list(tn_sub.gen_site_coos_present())
[(0, 2, 3), (1, 1, 3), (1, 2, 2), (1, 2, 3), (1, 3, 3), (2, 2, 3)]
get_ranges_present()[source]#

Return the range of site coordinates present in this TN.

Returns

xrange, yrange, zrange – The minimum and maximum site coordinates present in each direction.

Return type

tuple[tuple[int, int]]

Examples

>>> tn = qtn.TN3D_rand(4, 4, 4, 2)
>>> tn_sub = tn.select_local('I1,2,3', max_distance=1)
>>> tn_sub.get_ranges_present()
((0, 2), (1, 3), (2, 3))
maybe_convert_coo(coo)[source]#

Check if coo is a tuple of three ints and convert to the corresponding site tag if so.

property nsites#

The total number of sites.

site_tag(i, j, k)[source]#

The name of the tag specifiying the tensor at site (i, j, k).

property site_tag_id#

The string specifier for tagging each site of this 3D TN.

property site_tags#

All of the Lx * Ly site tags.

valid_coo(coo, xrange=None, yrange=None, zrange=None)[source]#

Check whether coo is in-bounds.

Parameters
  • coo ((int, int, int), optional) – The coordinates to check.

  • xrange ((int, int), optional) – The range of allowed values for the x, y, and z coordinates.

  • yrange ((int, int), optional) – The range of allowed values for the x, y, and z coordinates.

  • zrange ((int, int), optional) – The range of allowed values for the x, y, and z coordinates.

Return type

bool

property x_tag_id#

The string specifier for tagging each x-slice of this 3D TN.

property x_tags#

A tuple of all of the Lx different x-slice tags.

property y_tag_id#

The string specifier for tagging each y-slice of this 3D TN.

property y_tags#

A tuple of all of the Ly different y-slice tags.

property z_tag_id#

The string specifier for tagging each z-slice of this 3D TN.

property z_tags#

A tuple of all of the Lz different z-slice tags.

class quimb.tensor.tensor_3d.TensorNetwork3DFlat(ts, *, virtual=False, check_collisions=True)[source]#

Mixin class for a 3D square lattice tensor network with a single tensor per site, for example, both PEPS and PEPOs.

bond(coo1, coo2)[source]#

Get the name of the index defining the bond between sites at coo1 and coo2.

bond_size(coo1, coo2)[source]#

Return the size of the bond between sites at coo1 and coo2.

class quimb.tensor.tensor_3d.TensorNetwork3DVector(ts, *, virtual=False, check_collisions=True)[source]#

Mixin class for a 3D square lattice vector TN, i.e. one with a single physical index per site.

gate(G, where, contract=False, tags=None, info=None, inplace=False, **compress_opts)[source]#

Apply a gate G to sites where, preserving the outer site inds.

phys_dim(i=None, j=None, k=None)[source]#

Get the size of the physical indices / a specific physical index.

property site_inds#

All of the site inds.

to_dense(*inds_seq, **contract_opts)[source]#

Return the dense ket version of this 3D vector, i.e. a qarray with shape (-1, 1).

quimb.tensor.tensor_3d.cell_to_sites(p)[source]#

Turn a cell ((i0, j0, k0), (di, dj, dk)) into the sites it contains.

Examples

>>> cell_to_sites([(3, 4), (2, 2)])
((3, 4), (3, 5), (4, 4), (4, 5))
quimb.tensor.tensor_3d.gen_3d_bonds(Lx, Ly, Lz, steppers, coo_filter=None)[source]#

Convenience function for tiling pairs of bond coordinates on a 3D lattice given a function like lambda i, j, k: (i + 1, j + 1, k + 1).

Parameters
  • Lx (int) – The number of x-slices.

  • Ly (int) – The number of y-slices.

  • Lz (int) – The number of z-slices.

  • steppers (callable or sequence of callable) – Function(s) that take args (i, j, k) and generate another coordinate, thus defining a bond.

  • coo_filter (callable) – Function that takes args (i, j, k) and only returns True if this is to be a valid starting coordinate.

Yields

bond (tuple[tuple[int, int, int], tuple[int, int, int]]) – A pair of coordinates.

Examples

Generate nearest neighbor bonds:

>>> for bond in gen_3d_bonds(2, 2, 2, [lambda i, j, k: (i + 1, j, k),
...                                    lambda i, j, k: (i, j + 1, k),
...                                    lambda i, j, k: (i, j, k + 1)]):
...     print(bond)
((0, 0, 0), (1, 0, 0))
((0, 0, 0), (0, 1, 0))
((0, 0, 0), (0, 0, 1))
((0, 0, 1), (1, 0, 1))
((0, 0, 1), (0, 1, 1))
((0, 1, 0), (1, 1, 0))
((0, 1, 0), (0, 1, 1))
((0, 1, 1), (1, 1, 1))
((1, 0, 0), (1, 1, 0))
((1, 0, 0), (1, 0, 1))
((1, 0, 1), (1, 1, 1))
((1, 1, 0), (1, 1, 1))
quimb.tensor.tensor_3d.is_lone_coo(where)[source]#

Check if where has been specified as a single coordinate triplet.

quimb.tensor.tensor_3d.sites_to_cell(sites)[source]#

Get the minimum covering cell for sites.

Examples

>>> sites_to_cell([(1, 3, 3), (2, 2, 4)])
((1, 2, 3) , (2, 2, 2))