quimb.tensor.tensor_gen#
Generate specific tensor network states and operators.
Functions
Hyper tensor network representation of the 2D classical ising model partition function. |
|
Hyper tensor network representation of the 3D classical ising model partition function. |
|
Build a hyper tensor network representation of a classical ising model partition function by specifying graph edges. |
|
|
Create a hyper tensor network from a '.cnf' or '.wcnf' file - i.e. a model counting or weighted model counting instance specification. |
|
XXZ-Hamiltonian in MPO form. |
|
XY-Hamiltonian in MPO form. |
|
Hamiltonian of one-dimensional bilinear biquadratic chain in MPO form, see PhysRevB.93.184428. |
|
Heisenberg Hamiltonian in MPO form. |
|
Ising Hamiltonian in MPO form. |
|
The many-body-localized spin hamiltonian in MPO form. |
|
Generate an identity MPO of size |
|
Return an identity matrix operator with the same physical index and inds/tags as |
|
Return an MPO of bond dimension 1 representing the product of raw operators given in |
|
Generate a random matrix product state. |
|
Generate a random hermitian matrix product operator. |
|
Generate a zeros MPO of size |
|
Return a zeros matrix product operator with the same physical index and inds/tags as |
|
A computational basis state in Matrix Product State form. |
|
Build the chi=2 OBC MPS representation of the GHZ state. |
|
Generate the neel state in Matrix Product State form. |
|
Generate a product state in MatrixProductState form, i,e, with bond dimension 1, from single site vectors described by |
|
Generate a random computation basis state, like '01101001010'. |
|
Generate a random matrix product state. |
|
A product state for sampling tensor network traces. |
|
Build the chi=2 OBC MPS representation of the W state. |
|
The all-zeros MPS state, of given bond-dimension. |
The tensor network representation of the 2D classical ising model partition function. |
|
|
A scalar 2D lattice tensor network initialized with empty tensors. |
|
A scalar 2D lattice tensor network with tensors filled by a function. |
|
A random scalar 2D lattice tensor network. |
|
A scalar 2D lattice tensor network with every element set to |
Tensor network representation of the 3D classical ising model partition function. |
|
|
A scalar 3D lattice tensor network initialized with empty tensors. |
|
A scalar 3D lattice tensor network with tensors filled by a function. |
|
A random scalar 3D lattice tensor network. |
|
A scalar 2D lattice tensor network with every element set to |
Build a regular tensor network representation of a classical ising model partition function by specifying graph edges. |
|
|
Make a tensor network from sequence of graph edges that counts the number of ways to cover the graph exactly with dimers. |
|
Create a tensor network from a sequence of edges defining a graph, and a 'fill' function that maps shapes to data. |
|
Create a tensor network from a sequence of edges defining a graph, initialized with empty tensors. |
|
Create a random tensor network with geometry defined from a sequence of edges defining a graph. |
|
Create a tensor network from a sequence of edges defining a graph, initialized with a constant value. |
|
Create a random regular tensor network. |
|
The magnetic field term for the classical ising model. |
|
The interaction term for the classical ising model. |
|
The single effective TN site for the classical ising model. |
|
The sqrt factorized interaction term for the classical ising model. |
|
|
|
Get the set of PARAFAC tensors representing satisfiability of |
|
|
|
|
|
XXZ-Hamiltonian in |
|
XY-Hamiltonian in |
|
Hamiltonian of one-dimensional bilinear biquadratic chain in LocalHam1D form, see PhysRevB.93.184428. |
|
Heisenberg Hamiltonian in |
|
Ising Hamiltonian in |
|
The many-body-localized spin hamiltonian in |
|
Heisenberg Hamiltonian in |
|
Ising Hamiltonian in |
|
Heisenberg Hamiltonian in |
|
Heisenberg Hamiltonian in |
Check if |
|
|
Get the array representing satisfiability of |
|
Get the set of MPS tensors representing satisfiability of |
|
Get the set of PARAFAC arrays representing satisfiability of |
|
Get the tensor representing satisfiability of |
|
Generate a random tensor with specified shape and inds, and randomly 'phased' (distributed on the unit circle) data, such that |
|
Generate a random tensor with specified shape and inds. |
|
Generate tensor(s) for a spin hamiltonian MPO. |
Classes
|
Class for easily building custom spin hamiltonians in MPO or LocalHam1D form. |
- quimb.tensor.tensor_gen.HTN2D_classical_ising_partition_function(Lx, Ly, beta, h=0.0, j=1.0, ind_id='s{},{}', cyclic=False)[source]#
Hyper tensor network representation of the 2D classical ising model partition function. The indices will be shared by 4 or 5 tensors depending on whether
h
is non-zero. As opposed to the ‘normal’ tensor network, here each classical spin is still a single index, which is easier to contract exactly.- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
beta (float) – The inverse temperature.
h (float, optional) – The magnetic field strength.
j (float, optional) – The interaction strength, positive being ferromagnetic.
cyclic (bool or (bool, bool), optional) – Whether to use periodic boundary conditions. X and Y can be specified separately using a tuple.
ind_id (str, optional) – How to label the indices i.e.
ind_id.format(i, j)
, each of which corresponds to a single classical spin.
- Return type
- quimb.tensor.tensor_gen.HTN3D_classical_ising_partition_function(Lx, Ly, Lz, beta, j=1.0, h=0.0, cyclic=False, ind_id='s{},{},{}')[source]#
Hyper tensor network representation of the 3D classical ising model partition function. The indices will be shared by 6 or 7 tensors depending on whether
h
is non-zero. As opposed to the ‘normal’ tensor network, here each classical spin is still a single index, which is easier to contract exactly.- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
Lz (int) – Length of side z.
beta (float) – The inverse temperature.
j (float, optional) – The interaction strength, positive being ferromagnetic.
h (float, optional) – The magnetic field strength.
cyclic (bool or (bool, bool, bool), optional) – Whether to use periodic boundary conditions. X, Y and Z can be specified separately using a tuple.
ind_id (str, optional) – How to label the indices i.e.
ind_id.format(i, j, k)
, each of which corresponds to a single classical spin.
- Return type
- quimb.tensor.tensor_gen.HTN_classical_partition_function_from_edges(edges, beta, j=1.0, h=0.0, site_ind_id='s{}', site_tag_id='I{}', bond_tag_id='B{},{}')[source]#
Build a hyper tensor network representation of a classical ising model partition function by specifying graph edges. There will be a single tensor per interaction rather than per site, as well as a single tensor for each site, if
h != 0.0
.- Parameters
edges (sequence of tuple[hashable, hashable]) – The graph edges, as a sequence of pairs of hashable objects, for example integers, representing the nodes. You can redundantly specify
(u, v)
and(v, u)
and only one edge will be added.beta (float, optional) – The inverse temperature.
j (float, or callable, optional) – The interaction strength, positive being ferromagnetic. If a callable should have the signature
j(node_a, node_b)
and return a float.h (float, or callable, optional) – The magnetic field strength. If a callable should have the signature
h(node)
and return a float.site_ind_id (str, optional) – A string formatter for naming tensor indices like
site_ind_id.format(node)
.site_tag_id (str, optional) – A string formatter for naming tensor tags like
site_tag_id.format(node)
.bond_tag_id (str, optional) – A string formatter for naming tensor tags like
bond_tag_id.format(node_a, node_b)
.
- Return type
- quimb.tensor.tensor_gen.HTN_from_cnf(fname, mode='parafac')[source]#
Create a hyper tensor network from a ‘.cnf’ or ‘.wcnf’ file - i.e. a model counting or weighted model counting instance specification.
- Parameters
fname (str) – Path to a ‘.cnf’ or ‘.wcnf’ file.
mode ({'parafac', 'mps', 'dense', int}, optional) –
How to represent the clauses:
’parafac’ - N rank-2 tensors connected by a single hyper index. You could further call
hyperinds_resolve()
for more options to convert the hyper index into a (decomposed) COPY-tensor.’mps’ - N rank-3 tensors connected along a 1D line.
’dense’ - contract the hyper index.
int - use the ‘parafac’ mode, but only if the length of a clause is larger than this threshold.
- Returns
htn
- Return type
- quimb.tensor.tensor_gen.MPO_ham_XXZ(L, delta, jxy=1.0, *, S=0.5, cyclic=False, **mpo_opts)[source]#
XXZ-Hamiltonian in MPO form.
\[H_\mathrm{XXZ} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} + \Delta \sigma^Z_i \sigma^Z_{i + 1} ) - B_Z \sum_{i} \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
delta (float) – The ZZ-interaction strength. Positive is antiferromagnetic.
jxy (float, optional) – The X- and Y- interaction strength, defaults to 1. Positive is antiferromagnetic.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
mpo_opts – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_ham_XY(L, j=1.0, bz=0.0, *, S=0.5, cyclic=False, **mpo_opts)[source]#
XY-Hamiltonian in MPO form.
\[H_\mathrm{XY} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} ) - B_x \sum_{i} \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
j (float or (float, float), optional) – The XX and YY interaction strength. Positive is antiferromagnetic.
bz (float, optional) – The Z-magnetic field strength.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts (mpo_opts or) – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_ham_bilinear_biquadratic(L=None, theta=0, *, S=0.5, cyclic=False, compress=True, **mpo_opts)[source]#
Hamiltonian of one-dimensional bilinear biquadratic chain in MPO form, see PhysRevB.93.184428.
- Parameters
L (int) – The number of sites.
theta (float or (float, float), optional) – The parameter for linear and non-linear term of interaction strength, defaults to 0.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a hamiltonian with periodic boundary conditions or not, default is open boundary conditions.
mpo_opts – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_ham_heis(L, j=1.0, bz=0.0, *, S=0.5, cyclic=False, **mpo_opts)[source]#
Heisenberg Hamiltonian in MPO form.
\[H_\mathrm{Heis} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} + J_Z \sigma^Z_i \sigma^Z_{i + 1} ) - B_Z \sum_{i} \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
j (float or (float, float, float), optional) – The XX, YY and ZZ interaction strength. Positive is antiferromagnetic.
bz (float, optional) – The Z-magnetic field strength.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
mpo_opts – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_ham_ising(L, j=1.0, bx=0.0, *, S=0.5, cyclic=False, **mpo_opts)[source]#
Ising Hamiltonian in MPO form.
\[H_\mathrm{Ising} = J \sum_{i} \sigma^Z_i \sigma^Z_{i + 1} - B_x \sum_{i} \sigma^X_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
j (float, optional) – The ZZ interaction strength. Positive is antiferromagnetic.
bx (float, optional) – The X-magnetic field strength.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts (mpo_opts or) – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_ham_mbl(L, dh, j=1.0, seed=None, S=0.5, *, cyclic=False, dh_dist='s', dh_dim=1, beta=None, **mpo_opts)[source]#
The many-body-localized spin hamiltonian in MPO form.
\[H_\mathrm{MBL} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} + J_Z \sigma^Z_i \sigma^Z_{i + 1} ) - \sum_{i} h_i \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – Number of spins.
dh (float) – Random noise strength.
j (float, or (float, float, float), optional) – Interaction strength(s) e.g. 1 or (1., 1., 0.5). Positive is antiferromagnetic.
seed (int, optional) – Random number to seed the noise with.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Whether to use periodic boundary conditions - default is False.
dh_dist ({'s', 'g', 'qp'}, optional) – Whether to use sqaure, guassian or quasiperiodic noise.
beta (float, optional) – Frequency of the quasirandom noise, only if
dh_dist='qr'
.mpo_opts – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_identity(L, phys_dim=2, dtype='float64', cyclic=False, **mpo_opts)[source]#
Generate an identity MPO of size
L
.- Parameters
L (int) – The number of sites.
phys_dim (int, optional) – The physical (site) dimensions, defaults to 2.
dtype ({float, complex} or numpy dtype, optional) – Data type of the tensor network.
cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
mpo_opts – Supplied to
MatrixProductOperator
.
- quimb.tensor.tensor_gen.MPO_identity_like(mpo, **mpo_opts)[source]#
Return an identity matrix operator with the same physical index and inds/tags as
mpo
.
- quimb.tensor.tensor_gen.MPO_product_operator(arrays, cyclic=False, **mpo_opts)[source]#
Return an MPO of bond dimension 1 representing the product of raw operators given in
arrays
.- Parameters
arrays (sequence of 2D array_like) – The operators to form a tensor product of.
cyclic (bool, optional) – Whether to generate a cyclic MPO or not.
mpo_opts – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_rand(L, bond_dim, phys_dim=2, normalize=True, cyclic=False, herm=False, dtype='float64', **mpo_opts)[source]#
Generate a random matrix product state.
- Parameters
L (int) – The number of sites.
bond_dim (int) – The bond dimension.
phys_dim (int, optional) – The physical (site) dimensions, defaults to 2.
normalize (bool, optional) – Whether to normalize the operator such that
trace(A.H @ A) == 1
.cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
dtype ({float, complex} or numpy dtype, optional) – Data type of the tensor network.
herm (bool, optional) – Whether to make the matrix hermitian (or symmetric if real) or not.
mpo_opts – Supplied to
MatrixProductOperator
.
- quimb.tensor.tensor_gen.MPO_rand_herm(L, bond_dim, phys_dim=2, normalize=True, dtype='float64', **mpo_opts)[source]#
Generate a random hermitian matrix product operator. See
MPO_rand
.
- quimb.tensor.tensor_gen.MPO_zeros(L, phys_dim=2, dtype='float64', cyclic=False, **mpo_opts)[source]#
Generate a zeros MPO of size
L
.- Parameters
L (int) – The number of sites.
phys_dim (int, optional) – The physical (site) dimensions, defaults to 2.
dtype ({float, complex} or numpy dtype, optional) – Data type of the tensor network.
cyclic (bool, optional) – Generate a MPO with periodic boundary conditions or not, default is open boundary conditions.
mpo_opts – Supplied to
MatrixProductOperator
.
- Return type
- quimb.tensor.tensor_gen.MPO_zeros_like(mpo, **mpo_opts)[source]#
Return a zeros matrix product operator with the same physical index and inds/tags as
mpo
.- Parameters
mpo (MatrixProductOperator) – The MPO to copy the shape of.
- Return type
- quimb.tensor.tensor_gen.MPS_computational_state(binary, dtype='float64', cyclic=False, **mps_opts)[source]#
A computational basis state in Matrix Product State form.
- quimb.tensor.tensor_gen.MPS_ghz_state(L, dtype='float64', **mps_opts)[source]#
Build the chi=2 OBC MPS representation of the GHZ state.
- Parameters
L (int) – Number of qubits.
dtype ({'float64', 'complex128', 'float32', 'complex64'}, optional) – The underlying data type.
mps_opts – Supplied to
MatrixProductState
.
- quimb.tensor.tensor_gen.MPS_neel_state(L, down_first=False, dtype='float64', **mps_opts)[source]#
Generate the neel state in Matrix Product State form.
- quimb.tensor.tensor_gen.MPS_product_state(arrays, cyclic=False, **mps_opts)[source]#
Generate a product state in MatrixProductState form, i,e, with bond dimension 1, from single site vectors described by
arrays
.
- quimb.tensor.tensor_gen.MPS_rand_computational_state(L, dtype='float64', **mps_opts)[source]#
Generate a random computation basis state, like ‘01101001010’.
- Parameters
L (int) – The number of qubits.
seed (int, optional) – The seed to use.
dtype ({float, complex} or numpy dtype, optional) – Data type of the tensor network.
mps_opts – Supplied to
MatrixProductState
.
- quimb.tensor.tensor_gen.MPS_rand_state(L, bond_dim, phys_dim=2, normalize=True, cyclic=False, dtype='float64', trans_invar=False, **mps_opts)[source]#
Generate a random matrix product state.
- Parameters
L (int) – The number of sites.
bond_dim (int) – The bond dimension.
phys_dim (int, optional) – The physical (site) dimensions, defaults to 2.
normalize (bool, optional) – Whether to normalize the state.
cyclic (bool, optional) – Generate a MPS with periodic boundary conditions or not, default is open boundary conditions.
dtype ({float, complex} or numpy dtype, optional) – Data type of the tensor network.
trans_invar (bool (optional)) – Whether to generate a translationally invariant state, requires cyclic=True.
mps_opts – Supplied to
MatrixProductState
.
- quimb.tensor.tensor_gen.MPS_sampler(L, dtype=<class 'complex'>, squeeze=True, **mps_opts)[source]#
A product state for sampling tensor network traces. Seen as a vector it has the required property that
psi.H @ psi == d
always for hilbert space sized
.
- quimb.tensor.tensor_gen.MPS_w_state(L, dtype='float64', **mps_opts)[source]#
Build the chi=2 OBC MPS representation of the W state.
- Parameters
L (int) – Number of qubits.
dtype ({'float64', 'complex128', 'float32', 'complex64'}, optional) – The underlying data type.
mps_opts – Supplied to
MatrixProductState
.
- quimb.tensor.tensor_gen.MPS_zero_state(L, bond_dim=1, phys_dim=2, cyclic=False, dtype='float64', **mps_opts)[source]#
The all-zeros MPS state, of given bond-dimension.
- Parameters
L (int) – The number of sites.
bond_dim (int, optional) – The bond dimension, defaults to 1.
phys_dim (int, optional) – The physical (site) dimensions, defaults to 2.
cyclic (bool, optional) – Generate a MPS with periodic boundary conditions or not, default is open boundary conditions.
dtype ({float, complex} or numpy dtype, optional) – Data type of the tensor network.
mps_opts – Supplied to
MatrixProductState
.
- class quimb.tensor.tensor_gen.SpinHam1D(S=0.5, cyclic=False)[source]#
Class for easily building custom spin hamiltonians in MPO or LocalHam1D form. Currently limited to nearest neighbour interactions (and single site terms). It is possible to set ‘default’ translationally invariant terms, but also terms acting on specific sites only (which take precedence). It is also possible to build a sparse matrix version of the hamiltonian (obviously for small sizes only).
- Parameters
Examples
Initialize the spin hamiltonian builder:
>>> builder = SpinHam1D(S=3 / 2)
Add some two-site terms:
>>> builder += 0.5, '+', '-' >>> builder += 0.5, '-', '+' >>> builder += 1.0, 'Z', 'Z'
Add a single site term:
>>> builder -= 0.3, 'Z'
Build a MPO version of the hamiltonian for use with DMRG:
>>> mpo_ham = builder.build_mpo(100) >>> mpo_ham <MatrixProductOperator(tensors=100, L=100, max_bond=5)>
Build a LocalHam1D version of the hamiltonian for use with TEBD:
>>> builder.build_local_ham(100) <LocalHam1D(L=100, cyclic=False)>
You can also set terms for specific sites (this overides any of the ‘default’, translationally invariant terms set as above):
>>> builder[10, 11] += 0.75, '+', '-' >>> builder[10, 11] += 0.75, '-', '+' >>> builder[10, 11] += 1.5, 'Z', 'Z'
Or specific one-site terms (which again overides any default single site terms set above):
>>> builder[10] += 3.7, 'Z' >>> builder[11] += 0.0, 'I' # '0' term turns off field
- add_term(factor, *operators)[source]#
Add another term to the expression to be built.
- Parameters
factor (scalar) – Scalar factor to multiply this term by.
*operators (str or array) – The operators to use. Can specify one or two for single or two site terms respectively. Can use strings, which are supplied to
spin_operator()
, or actual arrays as long as they have the correct dimension.
- build_local_ham(L=None, **local_ham_1d_opts)[source]#
Build a nearest neighbour interactor instance of this spin hamiltonian of size
L
. See alsoLocalHam1D
.- Parameters
L (int, optional) – The number of spins, if the hamiltonian only has two-site terms this is optional.
- Return type
- build_mpo(L, upper_ind_id='k{}', lower_ind_id='b{}', site_tag_id='I{}', tags=None, bond_name='')[source]#
Build an MPO instance of this spin hamiltonian of size
L
. See alsoMatrixProductOperator
.
- quimb.tensor.tensor_gen.TN2D_classical_ising_partition_function(Lx, Ly, beta, j=1.0, h=0.0, cyclic=False, site_tag_id='I{},{}', row_tag_id='ROW{}', col_tag_id='COL{}')[source]#
The tensor network representation of the 2D classical ising model partition function.
- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
beta (float) – The inverse temperature.
j (float, optional) – The interaction strength, positive being ferromagnetic.
h (float, optional) – The magnetic field strength.
cyclic (bool or (bool, bool), optional) – Whether to use periodic boundary conditions. X and Y can be specified separately using a tuple.
site_tag_id (str, optional) – String specifier for naming convention of site tags.
row_tag_id (str, optional) – String specifier for naming convention of row tags.
col_tag_id (str, optional) – String specifier for naming convention of column tags.
- Return type
- quimb.tensor.tensor_gen.TN2D_empty(Lx, Ly, D, cyclic=False, site_tag_id='I{},{}', row_tag_id='ROW{}', col_tag_id='COL{}', dtype='float64')[source]#
A scalar 2D lattice tensor network initialized with empty tensors.
- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool), optional) – Whether to use periodic boundary conditions. X and Y can be specified separately using a tuple.
site_tag_id (str, optional) – String specifier for naming convention of site tags.
row_tag_id (str, optional) – String specifier for naming convention of row tags.
col_tag_id (str, optional) – String specifier for naming convention of column tags.
dtype (str, optional) – The data type of the tensors.
- Return type
- quimb.tensor.tensor_gen.TN2D_from_fill_fn(fill_fn, Lx, Ly, D, cyclic=False, site_tag_id='I{},{}', row_tag_id='ROW{}', col_tag_id='COL{}')[source]#
A scalar 2D lattice tensor network with tensors filled by a function.
- Parameters
fill_fn (callable) – A function with signature
fill_fn(shape) -> array
, used to fill each tensor.Lx (int) – Length of side x.
Ly (int) – Length of side y.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool), optional) – Whether to use periodic boundary conditions. X and Y can be specified separately using a tuple.
site_tag_id (str, optional) – String specifier for naming convention of site tags.
row_tag_id (str, optional) – String specifier for naming convention of row tags.
col_tag_id (str, optional) – String specifier for naming convention of column tags.
- Return type
- quimb.tensor.tensor_gen.TN2D_rand(Lx, Ly, D, cyclic=False, site_tag_id='I{},{}', row_tag_id='ROW{}', col_tag_id='COL{}', dtype='float64')[source]#
A random scalar 2D lattice tensor network.
- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool), optional) – Whether to use periodic boundary conditions. X and Y can be specified separately using a tuple.
site_tag_id (str, optional) – String specifier for naming convention of site tags.
row_tag_id (str, optional) – String specifier for naming convention of row tags.
col_tag_id (str, optional) – String specifier for naming convention of column tags.
dtype (dtype, optional) – Data type of the random arrays.
seed (int, optional) – A random seed.
- Return type
- quimb.tensor.tensor_gen.TN2D_with_value(value, Lx, Ly, D, cyclic=False, site_tag_id='I{},{}', row_tag_id='ROW{}', col_tag_id='COL{}', dtype=None)[source]#
A scalar 2D lattice tensor network with every element set to
value
. This usesnumpy.broadcast_to
and therefore essentially no memory.- Parameters
value (scalar) – The value to fill the tensors with.
Lx (int) – Length of side x.
Ly (int) – Length of side y.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool), optional) – Whether to use periodic boundary conditions. X and Y can be specified separately using a tuple.
site_tag_id (str, optional) – String specifier for naming convention of site tags.
row_tag_id (str, optional) – String specifier for naming convention of row tags.
col_tag_id (str, optional) – String specifier for naming convention of column tags.
dtype (str, optional) – The data type of the tensors.
- Return type
- quimb.tensor.tensor_gen.TN3D_classical_ising_partition_function(Lx, Ly, Lz, beta, j=1.0, h=0.0, cyclic=False, site_tag_id='I{},{},{}', x_tag_id='X{}', y_tag_id='Y{}', z_tag_id='Z{}')[source]#
Tensor network representation of the 3D classical ising model partition function.
- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
Lz (int) – Length of side z.
beta (float) – The inverse temperature.
j (float, optional) – The interaction strength, positive being ferromagnetic.
h (float, optional) – The magnetic field strength.
cyclic (bool or (bool, bool, bool), optional) – Whether to use periodic boundary conditions. X, Y and Z can be specified separately using a tuple.
site_tag_id (str, optional) – String formatter specifying how to label each site.
- Return type
- quimb.tensor.tensor_gen.TN3D_empty(Lx, Ly, Lz, D, cyclic=False, site_tag_id='I{},{},{}', x_tag_id='X{}', y_tag_id='Y{}', z_tag_id='Z{}', dtype='float64')[source]#
A scalar 3D lattice tensor network initialized with empty tensors.
- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
Lz (int) – Length of side z.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool, bool), optional) – Whether to use periodic boundary conditions. X, Y and Z can be specified separately using a tuple.
site_tag_id (str, optional) – String formatter specifying how to label each site.
dtype (dtype, optional) – Data type of the random arrays.
seed (int, optional) – Random seed.
- Return type
- quimb.tensor.tensor_gen.TN3D_from_fill_fn(fill_fn, Lx, Ly, Lz, D, cyclic=False, site_tag_id='I{},{},{}', x_tag_id='X{}', y_tag_id='Y{}', z_tag_id='Z{}')[source]#
A scalar 3D lattice tensor network with tensors filled by a function.
- Parameters
fill_fn (callable) – A function with signature
fill_fn(shape) -> array
, used to fill each tensor.Lx (int) – Length of side x.
Ly (int) – Length of side y.
Lz (int) – Length of side z.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool, bool), optional) – Whether to use periodic boundary conditions. X, Y and Z can be specified separately using a tuple.
site_tag_id (str, optional) – String formatter specifying how to label each site.
dtype (dtype, optional) – Data type of the random arrays.
- Return type
- quimb.tensor.tensor_gen.TN3D_rand(Lx, Ly, Lz, D, cyclic=False, site_tag_id='I{},{},{}', x_tag_id='X{}', y_tag_id='Y{}', z_tag_id='Z{}', dtype='float64')[source]#
A random scalar 3D lattice tensor network.
- Parameters
Lx (int) – Length of side x.
Ly (int) – Length of side y.
Lz (int) – Length of side z.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool, bool), optional) – Whether to use periodic boundary conditions. X, Y and Z can be specified separately using a tuple.
site_tag_id (str, optional) – String formatter specifying how to label each site.
dtype (dtype, optional) – Data type of the random arrays.
seed (int, optional) – Random seed.
- Return type
- quimb.tensor.tensor_gen.TN3D_with_value(value, Lx, Ly, Lz, D, cyclic=False, site_tag_id='I{},{},{}', x_tag_id='X{}', y_tag_id='Y{}', z_tag_id='Z{}', dtype=None)[source]#
A scalar 2D lattice tensor network with every element set to
value
. This usesnumpy.broadcast_to
and therefore essentially no memory.- Parameters
value (scalar) – The value to fill the tensors with.
Lx (int) – Length of side x.
Ly (int) – Length of side y.
Lz (int) – Length of side z.
D (int) – The bond dimension connecting sites.
cyclic (bool or (bool, bool, bool), optional) – Whether to use periodic boundary conditions. X, Y and Z can be specified separately using a tuple.
site_tag_id (str, optional) – String formatter specifying how to label each site.
dtype (dtype, optional) – Data type of the random arrays.
seed (int, optional) – Random seed.
- Return type
- quimb.tensor.tensor_gen.TN_classical_partition_function_from_edges(edges, beta, j=1.0, h=0.0, site_tag_id='I{}', bond_ind_id='b{},{}')[source]#
Build a regular tensor network representation of a classical ising model partition function by specifying graph edges. There will be a single tensor per site.
- Parameters
edges (sequence of tuple[hashable, hashable]) – The graph edges, as a sequence of pairs of hashable objects, for example integers, representing the nodes. You can redundantly specify
(u, v)
and(v, u)
and only one edge will be added.beta (float, optional) – The inverse temperature.
j (float, or callable, optional) – The interaction strength, positive being ferromagnetic. If a callable should have the signature
j(node_a, node_b)
and return a float.h (float, or callable, optional) – The magnetic field strength. If a callable should have the signature
h(node)
and return a float.site_tag_id (str, optional) – A string formatter for naming tensor tags like
site_ind_id.format(node)
.bond_ind_id (str, optional) – A string formatter for naming the indices bewteen tensors like
bond_ind_id.format(node_a, node_b)
.
- Return type
- quimb.tensor.tensor_gen.TN_dimer_covering_from_edges(edges, cover_count=1, site_tag_id='I{}', bond_ind_id='b{}, {}', dtype=<class 'float'>)[source]#
Make a tensor network from sequence of graph edges that counts the number of ways to cover the graph exactly with dimers. See https://arxiv.org/abs/1805.10598 for the construction.
- Parameters
edges (sequence of tuple) – The edges, each item should be a pair of hashable objects describing nodes linked.
cover_count (int, optional) – The exact number of times each node must be ‘covered’. For example 1 for a standard dimer covering or 2 for ‘ice rules’.
site_tag_id (str, optional) – A string formatter for naming tensor tags like
site_ind_id.format(node)
.bond_ind_id (str, optional) – A string formatter for naming the indices bewteen tensors like
bond_ind_id.format(node_a, node_b)
.
- Return type
- quimb.tensor.tensor_gen.TN_from_edges_and_fill_fn(fill_fn, edges, D, phys_dim=None, site_tag_id='I{}', site_ind_id='k{}')[source]#
Create a tensor network from a sequence of edges defining a graph, and a ‘fill’ function that maps shapes to data.
- Parameters
fill_fn (callable) – A function with signature
fill_fn(shape) -> array
, used to fill each tensor.edges (sequence of tuple[hashable, hashable]) – The graph edges, as a sequence of pairs of hashable objects, for example integers, representing the nodes. You can redundantly specify
(u, v)
and(v, u)
and only one edge will be added.D (int) – The bond dimension connecting tensors.
phys_dim (int, optional) – If not
None
, give each tensor a ‘physical’, free index of this size at each node.site_tag_id (str, optional) – String with formatter to tag sites.
site_ind_id (str, optional) – String with formatter to tag indices (if
phys_dim
specified).
- Return type
- quimb.tensor.tensor_gen.TN_from_edges_empty(edges, D, phys_dim=None, site_tag_id='I{}', site_ind_id='k{}', dtype='float64')[source]#
Create a tensor network from a sequence of edges defining a graph, initialized with empty tensors.
- Parameters
edges (sequence of tuple[hashable, hashable]) – The graph edges, as a sequence of pairs of hashable objects, for example integers, representing the nodes. You can redundantly specify
(u, v)
and(v, u)
and only one edge will be added.D (int) – The bond dimension connecting tensors.
phys_dim (int, optional) – If not
None
, give each tensor a ‘physical’, free index of this size at each node.site_tag_id (str, optional) – String with formatter to tag sites.
site_ind_id (str, optional) – String with formatter to tag indices (if
phys_dim
specified).dtype (str, optional) – The data type of the tensors.
- Return type
- quimb.tensor.tensor_gen.TN_from_edges_rand(edges, D, phys_dim=None, seed=None, dtype='float64', site_tag_id='I{}', site_ind_id='k{}')[source]#
Create a random tensor network with geometry defined from a sequence of edges defining a graph.
- Parameters
G (sequence of tuple[node, node]) – The edges defining a graph, each element should be a pair of nodes described by hashable objects.
D (int) – The bond dimension connecting tensors.
phys_dim (int, optional) – If not
None
, give each tensor a ‘physical’, free index of this size to mimic a wavefunction oflen(G)
sites.seed (int, optional) – A random seed.
site_tag_id (str, optional) – String with formatter to tag sites.
site_ind_id (str, optional) – String with formatter to tag indices (if
phys_dim
specified).
- Return type
- quimb.tensor.tensor_gen.TN_from_edges_with_value(value, edges, D, phys_dim=None, site_tag_id='I{}', site_ind_id='k{}', dtype=None)[source]#
Create a tensor network from a sequence of edges defining a graph, initialized with a constant value. This uses
numpy.broadcast_to
and therefore essentially no memory.- Parameters
value (scalar) – The value to fill the tensors with.
edges (sequence of tuple[hashable, hashable]) – The graph edges, as a sequence of pairs of hashable objects, for example integers, representing the nodes. You can redundantly specify
(u, v)
and(v, u)
and only one edge will be added.D (int) – The bond dimension connecting tensors.
phys_dim (int, optional) – If not
None
, give each tensor a ‘physical’, free index of this size at each node.site_tag_id (str, optional) – String with formatter to tag sites.
site_ind_id (str, optional) – String with formatter to tag indices (if
phys_dim
specified).dtype (str, optional) – The data type of the tensors.
- Return type
- quimb.tensor.tensor_gen.TN_rand_reg(n, reg, D, phys_dim=None, seed=None, dtype='float64', site_tag_id='I{}', site_ind_id='k{}')[source]#
Create a random regular tensor network.
- Parameters
n (int) – The number of tensors.
reg (int) – The degree of the tensor network (how many tensors each tensor connects to).
D (int) – The bond dimension connecting tensors.
phys_dim (int, optional) – If not
None
, give each tensor a ‘physical’, free index of this size to mimic a wavefunction ofn
sites.seed (int, optional) – A random seed.
site_tag_id (str, optional) – String with formatter to tag sites.
site_ind_id (str, optional) – String with formatter to tag indices (if
phys_dim
specified).
- Return type
- quimb.tensor.tensor_gen.classical_ising_H_matrix(beta, h=0.0)[source]#
The magnetic field term for the classical ising model.
- quimb.tensor.tensor_gen.classical_ising_S_matrix(beta, j=1.0)[source]#
The interaction term for the classical ising model.
- quimb.tensor.tensor_gen.classical_ising_T_matrix(beta, j=1.0, h=0.0, directions='lrud', asymm=None)[source]#
The single effective TN site for the classical ising model.
- quimb.tensor.tensor_gen.classical_ising_sqrtS_matrix(beta, j=1.0, asymm=None)[source]#
The sqrt factorized interaction term for the classical ising model. If
j
is negative you can supplyasymm='l'
or'r'
to keep the matrix real, but it must be paired with the opposite in a tensor network.
- quimb.tensor.tensor_gen.clause_parafac_tensors(ndim, m, inds, tags=None)[source]#
Get the set of PARAFAC tensors representing satisfiability of
ndim
clauses with unsatisfied condition encoded inm
labelled byinds
andtags
.
- quimb.tensor.tensor_gen.ham_1d_XXZ(L=None, delta=None, jxy=1.0, *, S=0.5, cyclic=False, **local_ham_1d_opts)[source]#
XXZ-Hamiltonian in
LocalHam1D
form.\[H_\mathrm{XXZ} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} + \Delta \sigma^Z_i \sigma^Z_{i + 1} ) - B_Z \sum_{i} \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
delta (float) – The ZZ-interaction strength. Positive is antiferromagnetic.
jxy (float, optional) – The X- and Y- interaction strength, defaults to 1. Positive is antiferromagnetic.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a hamiltonian with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts – Supplied to
LocalHam1D
.
- Return type
- quimb.tensor.tensor_gen.ham_1d_XY(L=None, j=1.0, bz=0.0, *, S=0.5, cyclic=False, **local_ham_1d_opts)[source]#
XY-Hamiltonian in
LocalHam1D
form.\[H_\mathrm{XY} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} ) - B_Z \sum_{i} \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
j (float or (float, float), optional) – The XX and YY interaction strength. Positive is antiferromagnetic.
bz (float, optional) – The Z-magnetic field strength.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a hamiltonian with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts – Supplied to
LocalHam1D
.
- Return type
- quimb.tensor.tensor_gen.ham_1d_bilinear_biquadratic(L=None, theta=0, *, S=0.5, cyclic=False, **local_ham_1d_opts)[source]#
Hamiltonian of one-dimensional bilinear biquadratic chain in LocalHam1D form, see PhysRevB.93.184428.
- Parameters
L (int) – The number of sites.
theta (float or (float, float), optional) – The parameter for linear and non-linear term of interaction strength, defaults to 0.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a hamiltonian with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts – Supplied to
LocalHam1D
.
- Return type
- quimb.tensor.tensor_gen.ham_1d_heis(L=None, j=1.0, bz=0.0, *, S=0.5, cyclic=False, **local_ham_1d_opts)[source]#
Heisenberg Hamiltonian in
LocalHam1D
form.\[H_\mathrm{Heis} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} + J_Z \sigma^Z_i \sigma^Z_{i + 1} ) - B_Z \sum_{i} \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
j (float or (float, float, float), optional) – The XX, YY and ZZ interaction strength. Positive is antiferromagnetic.
bz (float, optional) – The Z-magnetic field strength.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a hamiltonian with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts – Supplied to
LocalHam1D
.
- Return type
- quimb.tensor.tensor_gen.ham_1d_ising(L=None, j=1.0, bx=0.0, *, S=0.5, cyclic=False, **local_ham_1d_opts)[source]#
Ising Hamiltonian in
LocalHam1D
form.\[H_\mathrm{Ising} = J \sum_{i} \sigma^Z_i \sigma^Z_{i + 1} - B_x \sum_{i} \sigma^X_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – The number of sites.
j (float, optional) – The ZZ interaction strength. Positive is antiferromagnetic.
bx (float, optional) – The X-magnetic field strength.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Generate a hamiltonian with periodic boundary conditions or not, default is open boundary conditions.
local_ham_1d_opts (mpo_opts or) – Supplied to
LocalHam1D
.
- Return type
- quimb.tensor.tensor_gen.ham_1d_mbl(L, dh, j=1.0, seed=None, S=0.5, *, cyclic=False, dh_dist='s', dh_dim=1, beta=None, **local_ham_1d_opts)[source]#
The many-body-localized spin hamiltonian in
LocalHam1D
form.\[H_\mathrm{MBL} = \sum_{i} ( J_X \sigma^X_i \sigma^X_{i + 1} + J_Y \sigma^Y_i \sigma^Y_{i + 1} + J_Z \sigma^Z_i \sigma^Z_{i + 1} ) - \sum_{i} h_i \sigma^Z_i\]Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
L (int) – Number of spins.
dh (float) – Random noise strength.
j (float, or (float, float, float), optional) – Interaction strength(s) e.g. 1 or (1., 1., 0.5). Positive is antiferromagnetic.
seed (int, optional) – Random number to seed the noise with.
S ({1/2, 1, 3/2, ...}, optional) – The underlying spin of the system, defaults to 1/2.
cyclic (bool, optional) – Whether to use periodic boundary conditions - default is False.
dh_dist ({'s', 'g', 'qp'}, optional) – Whether to use sqaure, guassian or quasiperiodic noise.
beta (float, optional) – Frequency of the quasirandom noise, only if
dh_dist='qr'
.local_ham_1d_opts – Supplied to
LocalHam1D
.
- Return type
- quimb.tensor.tensor_gen.ham_2d_heis(Lx, Ly, j=1.0, bz=0.0, **local_ham_2d_opts)[source]#
Heisenberg Hamiltonian in
LocalHam2D
. form.\[H_\mathrm{Heis} = \sum_{<ij>} ( J_X \sigma^X_i \sigma^X_{j} + J_Y \sigma^Y_i \sigma^Y_{j} + J_Z \sigma^Z_i \sigma^Z_{j} ) - B_Z \sum_{i} \sigma^Z_{i}\]for nearest neighbors \(<ij>\). Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
- Return type
- quimb.tensor.tensor_gen.ham_2d_ising(Lx, Ly, j=1.0, bx=0.0, **local_ham_2d_opts)[source]#
Ising Hamiltonian in
LocalHam2D
form.\[H_\mathrm{Ising} = J \sum_{<ij>} \sigma^Z_i \sigma^Z_{j} - B_x \sum_{i} \sigma^X_i\]for nearest neighbors \(<ij>\). Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
- Return type
- quimb.tensor.tensor_gen.ham_2d_j1j2(Lx, Ly, j1=1.0, j2=0.5, bz=0.0, **local_ham_2d_opts)[source]#
Heisenberg Hamiltonian in
LocalHam2D
. form.\[H_\mathrm{Heis} = \sum_{<ij>} ( J_{1,X} \sigma^X_i \sigma^X_{j} + J_{1,Y} \sigma^Y_i \sigma^Y_{j} + J_{1,Z} \sigma^Z_i \sigma^Z_{j} ) + \sum_{<<ij>>} ( J_{2,X} \sigma^X_i \sigma^X_{j} + J_{2,Y} \sigma^Y_i \sigma^Y_{j} + J_{2,Z} \sigma^Z_i \sigma^Z_{j} ) - B_Z \sum_{i} \sigma^Z_{i}\]for nearest neighbors \(<ij>\) and diagonal next nearest neighbors \(<<ij>>\). Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
Lx (int) – The number of rows.
Ly (int) – The number of columns.
j2 (float or (float, float, float), optional) – The nearest neighbor XX, YY and ZZ interaction strength. Positive is antiferromagnetic.
j2 – The diagonal next nearest nearest XX, YY and ZZ interaction strength. Positive is antiferromagnetic.
bz (float, optional) – The Z-magnetic field strength.
local_ham_2d_opts – Supplied to
LocalHam2D
.
- Return type
- quimb.tensor.tensor_gen.ham_3d_heis(Lx, Ly, Lz, j=1.0, bz=0.0, **local_ham_3d_opts)[source]#
Heisenberg Hamiltonian in
LocalHam3D
. form.\[H_\mathrm{Heis} = \sum_{<ij>} ( J_X \sigma^X_i \sigma^X_{j} + J_Y \sigma^Y_i \sigma^Y_{j} + J_Z \sigma^Z_i \sigma^Z_{j} ) - B_Z \sum_{i} \sigma^Z_{i}\]for nearest neighbors \(<ij>\). Note the default convention of antiferromagnetic interactions and spin operators not Pauli matrices.
- Parameters
Lx (int) – The number of x-planes.
Ly (int) – The number of y-planes.
Ly – The number of z-planes.
j (float or (float, float, float), optional) – The XX, YY and ZZ interaction strength. Positive is antiferromagnetic.
bz (float, optional) – The Z-magnetic field strength.
local_ham_3d_opts – Supplied to
LocalHam3D
.
- Return type
LocalHam3D
- quimb.tensor.tensor_gen.maybe_make_real(X)[source]#
Check if
X
is real, if so, convert to contiguous array.
- quimb.tensor.tensor_gen.or_clause_data(ndim, m=0, dtype=<class 'float'>, q=2)[source]#
Get the array representing satisfiability of
ndim
clauses with unsatisfied condition encoded inm
.
- quimb.tensor.tensor_gen.or_clause_mps_tensors(ndim, m, inds, tags=None)[source]#
Get the set of MPS tensors representing satisfiability of
ndim
clauses with unsatisfied condition encoded inm
labelled byinds
andtags
.
- quimb.tensor.tensor_gen.or_clause_parafac_data(ndim, m)[source]#
Get the set of PARAFAC arrays representing satisfiability of
ndim
clauses with unsatisfied condition encoded inm
.
- quimb.tensor.tensor_gen.or_clause_tensor(ndim, m, inds, tags=None)[source]#
Get the tensor representing satisfiability of
ndim
clauses with unsatisfied condition encoded inm
labelled byinds
andtags
.
- quimb.tensor.tensor_gen.rand_phased(shape, inds, tags=None, dtype=<class 'complex'>)[source]#
Generate a random tensor with specified shape and inds, and randomly ‘phased’ (distributed on the unit circle) data, such that
T.H @ T == T.norm()**2 == T.size
.- Parameters
shape (sequence of int) – Size of each dimension.
inds (sequence of str) – Names of each dimension.
tags (sequence of str) – Labels to tag this tensor with.
dtype ({'complex128', 'complex64'}, optional) – The underlying data type - can only be complex.
- Return type
- quimb.tensor.tensor_gen.rand_tensor(shape, inds, tags=None, dtype='float64', left_inds=None)[source]#
Generate a random tensor with specified shape and inds.
- Parameters
shape (sequence of int) – Size of each dimension.
inds (sequence of str) – Names of each dimension.
tags (sequence of str) – Labels to tag this tensor with.
dtype ({'float64', 'complex128', 'float32', 'complex64'}, optional) – The underlying data type.
left_inds (sequence of str, optional) – Which, if any, indices to group as ‘left’ indices of an effective matrix. This can be useful, for example, when automatically applying unitary constraints to impose a certain flow on a tensor network but at the atomistic (Tensor) level.
- Return type
- quimb.tensor.tensor_gen.spin_ham_mpo_tensor(one_site_terms, two_site_terms, S=0.5, left_two_site_terms=None, which=None, cyclic=False)[source]#
Generate tensor(s) for a spin hamiltonian MPO.
- Parameters
one_site_terms (sequence of (scalar, operator)) – The terms that act on a single site, each
operator
can be a string suitable to be sent tospin_operator()
or an actual 2d-array.two_site_terms (sequence of (scalar, operator operator)) – The terms that act on two neighbouring sites, each
operator
can be a string suitable to be sent tospin_operator()
or an actual 2d-array.S (fraction, optional) – What size spin to use, defaults to spin-1/2.
left_two_site_terms (sequence of (scalar, operator operator), optional) – If the interaction to the left of this site has different spin terms then the equivalent list of terms for that site.
which ({None, 'L', 'R', 'A'}, optional) – If
None
, generate the middle tensor, if ‘L’ a left-end tensor, if ‘R’ a right-end tensor and if ‘A’ all three.cyclic (bool, optional) – Whether to use periodic boundary conditions - default False.
- Returns
The middle, left, right or all three MPO tensors.
- Return type