symmray.flat.flat_fermionic_array ================================= .. py:module:: symmray.flat.flat_fermionic_array .. autoapi-nested-parse:: Fermionic symmetric arrays with flat backend. Classes ------- .. autoapisummary:: symmray.flat.flat_fermionic_array.FermionicArrayFlat symmray.flat.flat_fermionic_array.Z2FermionicArrayFlat Functions --------- .. autoapisummary:: symmray.flat.flat_fermionic_array.perm_to_swaps Module Contents --------------- .. py:function:: perm_to_swaps(perm) Convert a permutation to a list of swaps. :param perm: The permutation to convert. Should be a permutation of the elements in range(len(perm)). :type perm: tuple[int, ...] :returns: A tuple of swaps that will perform the permutation. :rtype: tuple[tuple[int, int], ...] .. py:class:: FermionicArrayFlat(sectors, blocks, indices, phases=None, label=None, dummy_modes=None, symmetry=None) Bases: :py:obj:`symmray.fermionic_common.FermionicCommon`, :py:obj:`symmray.flat.flat_array_common.FlatArrayCommon`, :py:obj:`symmray.flat.flat_data_common.FlatCommon`, :py:obj:`symmray.array_common.ArrayCommon`, :py:obj:`symmray.common.SymmrayCommon` Fermionic abelian symmetric array with flat backend. :param sectors: The stack of sector keys, with shape (num_blocks, ndim). Each row represents a sector of a corresponding block, and each column represents a charge in a given axis. :type sectors: array_like :param blocks: The stack of array blocks, with shape (num_blocks, *shape_block), i.e. `ndim + 1` dimensions, where the first dimension is the block index, which should match the first dimension of `sectors`, and the rest are the dimensions of individual blocks. :type blocks: array_like :param indices: Indices describing the dualness and any subindex information for each dimension of the array. If bools are supplied, they will be converted to a FlatIndex with the corresponding dualness, and no subindex information. :type indices: sequence[FlatIndex] :param phases: An array of +/- 1 phases, with shape (num_blocks,), giving the phase of each block. If not supplied, all phases are assumed to be +1. :type phases: array_like, optional :param label: An optional label for the array, potentially needed for ordering dummy odd fermionic modes. :type label: hashable, optional :param dummy_modes: A sequence of dummy fermionic modes representing to effectively prepend to the array. :type dummy_modes: sequence[FermionicOperator], optional :param symmetry: The symmetry of the array, if not using a specific symmetry class. :type symmetry: str or Symmetry, optional .. py:attribute:: __slots__ :value: ('_blocks', '_dummy_modes', '_indices', '_phases', '_sectors', '_symmetry', 'backend') .. py:attribute:: fermionic :value: True .. py:attribute:: static_symmetry :value: None .. py:attribute:: _dummy_modes :value: () .. py:property:: phases The phases for each block. .. py:property:: parity The total parity of the array, 0 for even, 1 for odd. .. py:method:: check() Check the internal consistency of the array. .. py:method:: new_with(sectors, blocks, indices, label=None) -> FermionicArrayFlat Create a new flat fermionic array of the same class as this one. Unlike `copy`, this does not copy over any existing data and drops by default `label`, `phases`, and `dummy_modes`. .. py:method:: copy(deep=False) -> FermionicArrayFlat Create a copy of the array. .. py:method:: copy_with(sectors=None, blocks=None, indices=None, phases=None) -> FermionicArrayFlat A copy of this fermionic flat array with some attributes replaced. Note that checks are not performed on the new properties, this is intended for internal use. .. py:method:: modify(sectors=None, blocks=None, indices=None, phases=None, dummy_modes=None) -> FermionicArrayFlat Modify this fermionic flat array in place with some attributes replaced. Note that checks are not performed on the new properties, this is intended for internal use. .. py:method:: set_params(params) Set the underlying array blocks. .. py:method:: to_pytree() Convert this flat fermionic array to a pytree purely of non-symmray containers and objects. .. py:method:: from_pytree(data) :classmethod: Create a flat fermionic array from a pytree purely of non-symmray containers and objects. .. py:method:: _map_blocks(fn_sector=None, fn_block=None) .. py:method:: from_blocks(blocks, indices, phases=None, label=None, dummy_modes=None, symmetry=None) -> FermionicArrayFlat :classmethod: Create a fermionic flat array from an explicit dictionary of blocks, and sequence of indices or duals. :param blocks: A dictionary mapping sector keys (tuples of charges) to blocks (arrays). :type blocks: dict[tuple[int, ...], array_like] :param indices: A sequence of indices describing the dualness and any subindex information for each dimension of the array. If bools are supplied, they will be converted to a FlatIndex with the corresponding dualness, and no subindex information. :type indices: sequence[FlatIndex] | sequence[bool] :param phases: A dictionary mapping sector keys to +/- 1 phases. If not supplied, all phases are assumed to be +1. :type phases: dict[tuple[int, ...], int], optional :param label: An optional label for the array, potentially needed for ordering dummy odd fermionic modes. :type label: hashable, optional :param dummy_modes: A sequence of dummy fermionic modes representing to effectively prepend to the array. :type dummy_modes: sequence[FermionicOperator], optional :param symmetry: The symmetry of the array, if not using a specific symmetry class. :type symmetry: str or Symmetry, optional .. py:method:: from_blocksparse(x: symmray.sparse.sparse_fermionic_array.FermionicArray, symmetry=None) :classmethod: Create a fermionic flat array from a fermionic blocksparse array. :param x: The fermionic blocksparse array to convert. :type x: FermionicArray :param symmetry: The symmetry to use. If not supplied, the symmetry of `x` is used. :type symmetry: str or Symmetry, optional .. py:method:: to_blocksparse() -> symmray.sparse.sparse_fermionic_array.FermionicArray Convert this fermionic flat array to a fermionic blocksparse array. :returns: The equivalent fermionic blocksparse array. :rtype: FermionicArray .. py:method:: sort_stack(axes=None, all_axes=None, inplace=False) -> FermionicArrayFlat Lexicgraphic sort the stack of sectors and blocks according to the values of charges in the specified axes, optionally filling in the rest of the axes with the remaining axes in the order they appear. :param axes: The axes to sort by. If a single integer is given, it will be interpreted as the axis to sort by. If a tuple of integers is given, it will be interpreted as the axes to sort by in order. Default is None, if all_axes is also None or True, this will sort all axes in their current order. :type axes: int | tuple[int, ...], optional :param all_axes: Whether to include all non-specified axes as tie-breakers, after the specified axes. If ``None``, the default, this will be True if `axes` is not supplied explicitly, and False otherwise. :type all_axes: bool, optional :param inplace: Whether to perform the operation inplace or return a new array. Default is False, which returns a new array. :type inplace: bool, optional .. py:method:: phase_sync(inplace=False) -> FermionicArrayFlat Multiply all lazy phases into the block arrays. :param inplace: Whether to perform the operation in place. :type inplace: bool, optional :returns: The resolved array, which now has all trivial (+1) phases. :rtype: FermionicArrayFlat .. py:method:: phase_flip(*axs, inplace=False) -> FermionicArrayFlat Flip the phase of all sectors with odd parity at the given axis. :param ax: The axis along which to flip the phase. :type ax: int :param inplace: Whether to perform the operation in place. :type inplace: bool, optional :returns: The phase-flipped array. :rtype: FermionicArrayFlat .. py:method:: phase_transpose(axes=None, inplace=False) -> FermionicArrayFlat Phase this fermionic array as if it were transposed virtually, i.e. the actual arrays are not transposed. Useful when one wants the actual data layout to differ from the required fermionic mode layout. :param axes: The permutation of axes, by default None. :type axes: tuple of int, optional :param inplace: Whether to perform the operation in place. :type inplace: bool, optional :rtype: FermionicArrayFlat .. py:method:: phase_sector(sector, inplace=False) :abstractmethod: .. py:method:: phase_global(inplace=False) Flip the global phase of the array. :param inplace: Whether to perform the operation in place. :type inplace: bool, optional :rtype: FermionicArray .. py:method:: _resolve_dummy_modes_conj(phase_permutation=True) Assuming we have effectively taken the conjugate of a fermionic array with dummy modes, get their new order and compute any phase changes coming from moving back to the beginning of the index order. .. py:method:: _resolve_dummy_modes_combine(a, b) Calculate the new combined dummy modes and any associated global phases combing from contracting two fermionic arrays `a` and `b`. This modifies this array in place. .. py:method:: _resolve_dummy_modes_squeeze(axes_squeeze) Assuming we are about to squeeze away `axes_squeeze`, compute the phases associated with moving them to the beginning of the array, and then turn them into dummy modes. .. py:method:: transpose(axes=None, phase=True, inplace=False) -> FermionicArrayFlat Transpose this flat fermionic array, by default accounting for the phases accumulated from swapping odd charges. :param axes: The permutation of axes, by default None. :type axes: tuple of int, optional :param phase: Whether to flip the phase of sectors whose odd charges undergo a odd permutation. By default True. :type phase: bool, optional :param inplace: Whether to perform the operation in place. :type inplace: bool, optional :returns: The transposed array. :rtype: FermionicArrayFlat .. py:method:: conj(phase_permutation=True, phase_dual=False, inplace=False) Conjugate this flat fermionic array. By default this include phases from both the virtual flipping of all axes, but *not* the conjugation of dual indices, such that:: ( tensordot_fermionic(x.conj(), x, ndim) == tensordot_fermionic(x, x.conj(), ndim) ) If all indices have matching dualness (i.e. all bra or all ket), *or* you set `phase_dual=True` then the above contractions will also be equal to ``x.norm() ** 2``. :param phase_permutation: Whether to flip the phase of sectors whose odd charges undergo a odd permutation due to *virtually* flipping the order of axes, by default True. :type phase_permutation: bool, optional :param phase_dual: Whether to flip the phase of dual indices, by default False. If a FermionicArrayFlat has a mix of dual and non-dual indices, and you are explicitly forming the norm, you may want to set this to True. But if it is part of a large tensor network you only need to flip the phase of true 'outer' dual indices. :type phase_dual: bool, optional :rtype: FermionicArrayFlat .. py:method:: dagger(phase_dual=False, inplace=False) Fermionic conjugate transpose. .. py:method:: eigh(phase_eigenvalues=True) -> tuple[symmray.flat.flat_vector.FlatVector, FermionicArrayFlat] Hermitian eigen-decomposition of this flat fermionic array. :param phase_eigenvalues: If True, any local phase will be absorbed into the eigenvalues, such that `U @ diag(w) @ U.H == a` always holds. By default True. If `False` then one of `U` or `U.H` should be phased individually to account for local phasesin the above expression. :type phase_eigenvalues: bool, optional :returns: * **eigenvalues** (*FlatVector*) -- The eigenvalues. * **eigenvectors** (*FermionicArrayFlat*) -- The abelian array of right eigenvectors. .. py:class:: Z2FermionicArrayFlat(sectors, blocks, indices, phases=None, label=None, dummy_modes=None, symmetry=None) Bases: :py:obj:`FermionicArrayFlat` Fermionic abelian symmetric array with flat backend. :param sectors: The stack of sector keys, with shape (num_blocks, ndim). Each row represents a sector of a corresponding block, and each column represents a charge in a given axis. :type sectors: array_like :param blocks: The stack of array blocks, with shape (num_blocks, *shape_block), i.e. `ndim + 1` dimensions, where the first dimension is the block index, which should match the first dimension of `sectors`, and the rest are the dimensions of individual blocks. :type blocks: array_like :param indices: Indices describing the dualness and any subindex information for each dimension of the array. If bools are supplied, they will be converted to a FlatIndex with the corresponding dualness, and no subindex information. :type indices: sequence[FlatIndex] :param phases: An array of +/- 1 phases, with shape (num_blocks,), giving the phase of each block. If not supplied, all phases are assumed to be +1. :type phases: array_like, optional :param label: An optional label for the array, potentially needed for ordering dummy odd fermionic modes. :type label: hashable, optional :param dummy_modes: A sequence of dummy fermionic modes representing to effectively prepend to the array. :type dummy_modes: sequence[FermionicOperator], optional :param symmetry: The symmetry of the array, if not using a specific symmetry class. :type symmetry: str or Symmetry, optional .. py:attribute:: static_symmetry