Generic iPEPS¶
Class describes iPEPS with rectangular unit cell, together with convenience functions for reading and writing such iPEPS from and into (JSON) files.
- class ipeps.ipeps.IPEPS(sites=None, vertexToSite=None, lX=None, lY=None, peps_args=<config.PEPSARGS object>, global_args=<config.GLOBALARGS object>)[source]¶
- Parameters:
sites (dict[tuple(int,int) : torch.tensor]) – map from elementary unit cell to on-site tensors
vertexToSite (function(tuple(int,int))->tuple(int,int)) – function mapping arbitrary vertex of a square lattice into a vertex within elementary unit cell
lX (int) – length of the elementary unit cell in X direction
lY (int) – length of the elementary unit cell in Y direction
peps_args (PEPSARGS) – ipeps configuration
global_args (GLOBALARGS) – global configuration
Member
sites
is a dictionary of non-equivalent on-site tensors indexed by tuple of coordinates (x,y) within the elementary unit cell. The index-position convetion for on-site tensors is defined as follows:u s |/ l--a--r <=> a[s,u,l,d,r] | d
where s denotes physical index, and u,l,d,r label four principal directions up, left, down, right in anti-clockwise order starting from up. Member
vertexToSite
is a mapping function from any vertex (x,y) on a square lattice passed in as tuple(int,int) to a corresponding vertex within elementary unit cell.On-site tensor of an IPEPS object
wfc
at vertex (x,y) is conveniently accessed through the member functionsite
, which internally usesvertexToSite
mapping:coord= (0,0) a_00= wfc.site(coord)
By combining the appropriate
vertexToSite
mapping function with elementary unit cell specified throughsites
, various tilings of a square lattice can be achieved:# Example 1: 1-site translational iPEPS sites={(0,0): a} def vertexToSite(coord): return (0,0) wfc= IPEPS(sites,vertexToSite) # resulting tiling: # y\x -2 -1 0 1 2 # -2 a a a a a # -1 a a a a a # 0 a a a a a # 1 a a a a a # Example 2: 2-site bipartite iPEPS sites={(0,0): a, (1,0): b} def vertexToSite(coord): x = (coord[0] + abs(coord[0]) * 2) % 2 y = abs(coord[1]) return ((x + y) % 2, 0) wfc= IPEPS(sites,vertexToSite) # resulting tiling: # y\x -2 -1 0 1 2 # -2 A b a b a # -1 B a b a b # 0 A b a b a # 1 B a b a b # Example 3: iPEPS with 3x2 unit cell with PBC sites={(0,0): a, (1,0): b, (2,0): c, (0,1): d, (1,1): e, (2,1): f} wfc= IPEPS(sites,lX=3,lY=2) # resulting tiling: # y\x -2 -1 0 1 2 # -2 b c a b c # -1 e f d e f # 0 b c a b c # 1 e f d e f
where in the last example a default setting for
vertexToSite
is used, which maps square lattice into elementary unit cell of sizelX
xlY
assuming periodic boundary conditions (PBC) along both X and Y directions.- add_noise(noise, noise_f=None)[source]¶
- Parameters:
noise (float) – magnitude of the noise
Take IPEPS and add random uniform noise with magnitude
noise
to all on-site tensors
- get_checkpoint()[source]¶
- Returns:
all data necessary to reconstruct the state. In this case member
sites
- Return type:
dict[tuple(int,int): torch.tensor]
This function is called by optimizer to create checkpoints during the optimization process.
- get_parameters()[source]¶
- Returns:
variational parameters of iPEPS
- Return type:
iterable
This function is called by optimizer to access variational parameters of the state.
- load_checkpoint(checkpoint_file)[source]¶
- Parameters:
checkpoint_file (str) – path to checkpoint file
Initializes the state according to the supplied checkpoint file.
Note
The vertexToSite mapping function is not a part of checkpoint and must be provided either when instantiating IPEPS_ABELIAN or afterwards.
- site(coord)[source]¶
- Parameters:
coord (tuple(int,int)) – tuple (x,y) specifying vertex on a square lattice
- Returns:
on-site tensor corresponding to the vertex (x,y)
- Return type:
torch.tensor
- write_to_file(outputfile, aux_seq=[0, 1, 2, 3], tol=1e-14, normalize=False)[source]¶
Writes state to file. See
write_ipeps()
.
- ipeps.ipeps.extend_bond_dim(state, new_d)[source]¶
- Parameters:
state (IPEPS) – wavefunction to modify
new_d (int) – new enlarged auxiliary bond dimension
- Returns:
wavefunction with enlarged auxiliary bond dimensions
- Return type:
Take IPEPS and enlarge all auxiliary bond dimensions of all on-site tensors up to size
new_d
- ipeps.ipeps.read_ipepo(jsonfile, vertexToSite=None, aux_seq=[0, 1, 2, 3], peps_args=<config.PEPSARGS object>, global_args=<config.GLOBALARGS object>)[source]¶
- Parameters:
jsonfile (str or Path object) – input file describing iPEPO in json format
vertexToSite (function(tuple(int,int))->tuple(int,int)) – function mapping arbitrary vertex of a square lattice into a vertex within elementary unit cell
aux_seq (list[int]) – array specifying order of auxiliary indices of on-site tensors stored in jsonfile
peps_args (PEPSARGS) – ipeps configuration
global_args (GLOBALARGS) – global configuration
- Returns:
wavefunction
- Return type:
IPEPO
A simple PBC
vertexToSite
function is used by defaultParameter
aux_seq
defines the expected order of auxiliary indices in input file relative to the convention fixed in tn-torch:0 1A3 <=> [up, left, down, right]: aux_seq=[0,1,2,3] 2 for alternative order, eg. 1 0A2 <=> [left, up, right, down]: aux_seq=[1,0,3,2] 3
- ipeps.ipeps.read_ipeps(jsonfile, vertexToSite=None, aux_seq=[0, 1, 2, 3], peps_args=<config.PEPSARGS object>, global_args=<config.GLOBALARGS object>)[source]¶
- Parameters:
jsonfile (str or Path object) – input file describing iPEPS in json format
vertexToSite (function(tuple(int,int))->tuple(int,int)) – function mapping arbitrary vertex of a square lattice into a vertex within elementary unit cell
aux_seq (list[int]) – array specifying order of auxiliary indices of on-site tensors stored in jsonfile
peps_args (PEPSARGS) – ipeps configuration
global_args (GLOBALARGS) – global configuration
- Returns:
wavefunction
- Return type:
A simple PBC
vertexToSite
function is used by defaultParameter
aux_seq
defines the expected order of auxiliary indices in input file relative to the convention fixed in tn-torch:0 1A3 <=> [up, left, down, right]: aux_seq=[0,1,2,3] 2 for alternative order, eg. 1 0A2 <=> [left, up, right, down]: aux_seq=[1,0,3,2] 3
- ipeps.ipeps.write_ipepo(state, outputfile, tol=1e-14, normalize=False, peps_args=<config.PEPSARGS object>, global_args=<config.GLOBALARGS object>)[source]¶
- Parameters:
state (IPEPO) – operator to write out in json format
outputfile – target file
aux_seq (list[int]) – array specifying order in which the auxiliary indices of on-site tensors will be stored in the outputfile
tol (float) – minimum magnitude of tensor elements which are written out
normalize (bool) – if True, on-site tensors are normalized before writing
Parameter
aux_seq
defines the order of auxiliary indices relative to the convention fixed in tn-torch in which the tensor elements are written out:0 1A3 <=> [up, left, down, right]: aux_seq=[0,1,2,3] 2 for alternative order, eg. 1 0A2 <=> [left, up, right, down]: aux_seq=[1,0,3,2] 3
- ipeps.ipeps.write_ipeps(state, outputfile, aux_seq=[0, 1, 2, 3], tol=1e-14, normalize=False, peps_args=<config.PEPSARGS object>, global_args=<config.GLOBALARGS object>)[source]¶
- Parameters:
state (IPEPS) – wavefunction to write out in json format
outputfile – target file
aux_seq (list[int]) – array specifying order in which the auxiliary indices of on-site tensors will be stored in the outputfile
tol (float) – minimum magnitude of tensor elements which are written out
normalize (bool) – if True, on-site tensors are normalized before writing
Parameter
aux_seq
defines the order of auxiliary indices relative to the convention fixed in tn-torch in which the tensor elements are written out:0 1A3 <=> [up, left, down, right]: aux_seq=[0,1,2,3] 2 for alternative order, eg. 1 0A2 <=> [left, up, right, down]: aux_seq=[1,0,3,2] 3