PlaneWaveExcitationEngine
full name: tenpy.algorithms.plane_wave_excitation.PlaneWaveExcitationEngine
parent module:
tenpy.algorithms.plane_wave_excitation
type: class
Inheritance Diagram
Methods
|
|
Compute the energy of excited states |
|
Gives an approximate prediction for the required memory usage. |
|
Return necessary data to resume a |
|
Infinite sum to the left, see Eq. |
|
Infinite sum to the right, see Eq. |
|
Initial guess for the X tensors within a fixed charge sector. |
|
Resume a run that was interrupted. |
|
|
Run the plane-wave algorithm to find excited states of the given model. |
Initialize algorithm from another algorithm instance of a different class. |
- class tenpy.algorithms.plane_wave_excitation.PlaneWaveExcitationEngine(psi, model, options, **kwargs)[source]
Bases:
Algorithm
Base engine to compute quasiparticle excitations for uniform MPS.
Parameters are the same as for
Algorithm
.Options
- config PlaneWaveExcitationEngine
option summary E_boost in MultiSitePlaneWaveExcitationEngine.run
uniform strength of the energy boosts (instead of specifying a list), defau [...]
Dictionary as returned by ``self.env.get_initialization_data()`` from [...]
Energy shift from contracting the infinite environments. If `None`, compute [...]
Lanczos parameters as described in :cfg:config:`KrylovBased`.
max_N_sites_per_ring (from Algorithm) in Algorithm
Threshold for raising errors on too many sites per ring. Default ``18``. [...]
sum_iterations in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Maximum number of iterations for the explicit summation (default sum_iterat [...]
sum_method in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Whether to explicitly sum the environment by applying the unit cell tensors [...]
sum_tol in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Convergence criterion for the explicit summation.
trunc_params (from Algorithm) in Algorithm
Truncation parameters as described in :cfg:config:`truncation`.
- option lanczos_params: dict
Lanczos parameters as described in
KrylovBased
.
- option lambda_C1: float
Energy shift from contracting the infinite environments. If None, compute it again.
- option init_env_data: dict
Dictionary as returned by
self.env.get_initialization_data()
fromget_initialization_data()
.
- psi
The uniform MPS for which we compute (orthogonal) excitations.
- Type:
- run(p, qtotal_change=None, orthogonal_to=[], E_boosts=[], num_ev=1)[source]
Run the plane-wave algorithm to find excited states of the given model.
- Parameters:
p (float) – The momentum of the state; for unit cells larger than 1, we already include the factor of the smaller Brillouin zone: p*L.
qtotal_change (list of int) – Charge sectors for each of the defined charges.
orthogonal_to (list of list of
Array
) – Find excitations orthogonal to previously found tensors X.E_boosts (list of float) – energy boosts for orthogonal states
num_ev (int) – Number of eigenvalues and eigenvectors, that we extract from a single Arnoldi/ Lanczos run
- Returns:
Es (list of float) – Energies of the lowest-energy excitations. Number equal to num_ev.
psis (list of
MomentumMPS
) – MomentumMPS corresponding to the lowest-energy excitations. Number equal to num_ev.
Options
- config PlaneWaveExcitationEngine
option summary E_boost in MultiSitePlaneWaveExcitationEngine.run
uniform strength of the energy boosts (instead of specifying a list), defau [...]
init_env_data in PlaneWaveExcitationEngine
Dictionary as returned by ``self.env.get_initialization_data()`` from [...]
lambda_C1 in PlaneWaveExcitationEngine
Energy shift from contracting the infinite environments. If `None`, compute [...]
lanczos_params in PlaneWaveExcitationEngine
Lanczos parameters as described in :cfg:config:`KrylovBased`.
max_N_sites_per_ring (from Algorithm) in Algorithm
Threshold for raising errors on too many sites per ring. Default ``18``. [...]
sum_iterations in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Maximum number of iterations for the explicit summation (default sum_iterat [...]
sum_method in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Whether to explicitly sum the environment by applying the unit cell tensors [...]
sum_tol in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Convergence criterion for the explicit summation.
trunc_params (from Algorithm) in Algorithm
Truncation parameters as described in :cfg:config:`truncation`.
- option E_boost: float
uniform strength of the energy boosts (instead of specifying a list), default to E_boost=100
- resume_run()[source]
Resume a run that was interrupted.
In case we saved an intermediate result at a
checkpoint
, this function allows to resume therun()
of the algorithm (after re-initialization with the resume_data). Since most algorithms just have a while loop with break conditions, the default behavior implemented here is to just callrun()
.
- infinite_sum_right(p, X)[source]
Infinite sum to the right, see Eq. (194) in [vanderstraeten2019]
- pfloat
The momentum of the state; for unit cells larger than 1, we already include the factor of the smaller Brillouin zone: p*L.
- Xlist of
Array
Current excitation tensors for each site of the unit cell.
- Returns:
R_sum – Array representing the right environment including the B tensors at each site.
- Return type:
Options
- config PlaneWaveExcitationEngine
option summary E_boost in MultiSitePlaneWaveExcitationEngine.run
uniform strength of the energy boosts (instead of specifying a list), defau [...]
init_env_data in PlaneWaveExcitationEngine
Dictionary as returned by ``self.env.get_initialization_data()`` from [...]
lambda_C1 in PlaneWaveExcitationEngine
Energy shift from contracting the infinite environments. If `None`, compute [...]
lanczos_params in PlaneWaveExcitationEngine
Lanczos parameters as described in :cfg:config:`KrylovBased`.
max_N_sites_per_ring (from Algorithm) in Algorithm
Threshold for raising errors on too many sites per ring. Default ``18``. [...]
sum_iterations in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Maximum number of iterations for the explicit summation (default sum_iterat [...]
sum_method in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Whether to explicitly sum the environment by applying the unit cell tensors [...]
sum_tol in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Convergence criterion for the explicit summation.
trunc_params (from Algorithm) in Algorithm
Truncation parameters as described in :cfg:config:`truncation`.
-
option sum_method:
explicit
|GMRES
Whether to explicitly sum the environment by applying the unit cell tensors until convergence (default) or solving the geometric series with the GMRES method.
- option sum_tol: float
Convergence criterion for the explicit summation.
- option sum_iterations: int
Maximum number of iterations for the explicit summation (default sum_iterations=100).
-
option sum_method:
- infinite_sum_left(p, X)[source]
Infinite sum to the left, see Eq. (194) in [vanderstraeten2019]
- pfloat
The momentum of the state; for unit cells larger than 1, we already include the factor of the smaller Brillouin zone: p*L.
- Xlist of
Array
Current excitation tensors for each site of the unit cell.
- Returns:
L_sum – Array representing the left environment including the B tensors at each site.
- Return type:
Options
- config PlaneWaveExcitationEngine
option summary E_boost in MultiSitePlaneWaveExcitationEngine.run
uniform strength of the energy boosts (instead of specifying a list), defau [...]
init_env_data in PlaneWaveExcitationEngine
Dictionary as returned by ``self.env.get_initialization_data()`` from [...]
lambda_C1 in PlaneWaveExcitationEngine
Energy shift from contracting the infinite environments. If `None`, compute [...]
lanczos_params in PlaneWaveExcitationEngine
Lanczos parameters as described in :cfg:config:`KrylovBased`.
max_N_sites_per_ring (from Algorithm) in Algorithm
Threshold for raising errors on too many sites per ring. Default ``18``. [...]
sum_iterations in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Maximum number of iterations for the explicit summation (default sum_iterat [...]
sum_method in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Whether to explicitly sum the environment by applying the unit cell tensors [...]
sum_tol in MultiSitePlaneWaveExcitationEngine.infinite_sum_right
Convergence criterion for the explicit summation.
trunc_params (from Algorithm) in Algorithm
Truncation parameters as described in :cfg:config:`truncation`.
-
option sum_method:
explicit
|GMRES
Whether to explicitly sum the environment by applying the unit cell tensors until convergence (default) or solving the geometric series with the GMRES method.
- option sum_tol: float
Convergence criterion for the explicit summation.
- option sum_iterations: int
Maximum number of iterations for the explicit summation (default sum_iterations=100).
-
option sum_method:
- class Aligned_Effective_H(outer)[source]
Bases:
NpcLinearOperator
Class defining the effective Hamiltonian for the excitation tensors X.
Where the B tensors are in the same unit cell as the tensors we want to update.
For a single-site unit cell the effective Hamiltonian looks like this:
| .--- B ---. | | | | | LW----W0----RW | | | | | .---VL- --.
- Parameters:
outer (
PlaneWaveExcitationEngine
) – Parent engine for the plane wave excitation ansatz.
- matvec(vec)[source]
Calculate the action of the operator on a vector vec.
Note that we don’t require vec to be one-dimensional. However, for square operators we require that the result of matvec has the same legs (in the same order) as vec such that they can be added. Note that this excludes a non-trivial qtotal for square operators.
- class Unaligned_Effective_H(outer, p)[source]
Bases:
NpcLinearOperator
Class defining the effective Hamiltonian for the excitation tensors X, where the B tensors are in left (LB) or right (RB) environment.
For a single-site unit cell the effective Hamiltonian looks like this:
| .--- AR ---. .--- AL ---. | ip | | | -ip | | | | e LB----W0----RW + e LW----W0----RB | | | | | | | | .---VL- --. .---VL- --.
- Parameters:
outer (
PlaneWaveExcitationEngine
) – Parent engine for the plane wave excitation ansatz.p (float) – The momentum of the state; for unit cells larger than 1, we already include the factor of the smaller Brillouin zone: p*L.
- matvec(vec)[source]
Calculate the action of the operator on a vector vec.
Note that we don’t require vec to be one-dimensional. However, for square operators we require that the result of matvec has the same legs (in the same order) as vec such that they can be added. Note that this excludes a non-trivial qtotal for square operators.
- estimate_RAM(mem_saving_factor=None)[source]
Gives an approximate prediction for the required memory usage.
This calculation is based on the requested bond dimension, the local Hilbert space dimension, the number of sites, and the boundary conditions.
- Parameters:
mem_saving_factor (float) – Represents the amount of RAM saved due to conservation laws. By default, it is ‘None’ and is extracted from the model automatically. However, this is only possible in a few cases and needs to be estimated in most cases. This is due to the fact that it is dependent on the model parameters. If one has a better estimate, one can pass the value directly. This value can be extracted by building the initial state psi (usually by performing DMRG) and then calling
print(psi.get_B(0).sparse_stats())
TeNPy will automatically print the fraction of nonzero entries in the first line, for example,6 of 16 entries (=0.375) nonzero
. This fraction corresponds to the mem_saving_factor; in our example, it is 0.375.- Returns:
usage – Required RAM in MB.
- Return type:
See also
tenpy.simulations.simulation.estimate_simulation_RAM
global function calling this.
- get_resume_data(sequential_simulations=False)[source]
Return necessary data to resume a
run()
interrupted at a checkpoint.At a
checkpoint
, you can savepsi
,model
andoptions
along with the data returned by this function. When the simulation aborts, you can resume it using this saved data with:eng = AlgorithmClass(psi, model, options, resume_data=resume_data) eng.resume_run()
An algorithm which doesn’t support this should override resume_run to raise an Error.
- Parameters:
sequential_simulations (bool) – If True, return only the data for re-initializing a sequential simulation run, where we “adiabatically” follow the evolution of a ground state (for variational algorithms), or do series of quenches (for time evolution algorithms); see
run_seq_simulations()
.- Returns:
resume_data – Dictionary with necessary data (apart from copies of psi, model, options) that allows to continue the algorithm run from where we are now. It might contain an explicit copy of psi.
- Return type:
- classmethod switch_engine(other_engine, *, options=None, **kwargs)[source]
Initialize algorithm from another algorithm instance of a different class.
You can initialize one engine from another, not too different subclasses. Internally, this function calls
get_resume_data()
to extract data from the other_engine and then initializes the new class.Note that it transfers the data without making copies in most case; even the options! Thus, when you call run() on one of the two algorithm instances, it will modify the state, environment, etc. in the other. We recommend to make the switch as
engine = OtherSubClass.switch_engine(engine)
directly replacing the reference.- Parameters:
cls (class) – Subclass of
Algorithm
to be initialized.other_engine (
Algorithm
) – The engine from which data should be transferred. Another, but not too different algorithm subclass-class; e.g. you can switch from theTwoSiteDMRGEngine
to theOneSiteDMRGEngine
.options (None | dict-like) – If not None, these options are used for the new initialization. If None, take the options from the other_engine.
**kwargs – Further keyword arguments for class initialization. If not defined, resume_data is collected with
get_resume_data()
.