VMPS++
|
Solver that calculates the ground state of a UMPS. Analogue of the DmrgSolver class.
Symmetry | : A class providing all relevant functions and infos that are determined by the Symmetry. Click here for more information. |
Scalar | : double or complex<double> |
Definition at line 32 of file VumpsSolver.h.
#include <VumpsSolver.h>
Public Member Functions | |
VumpsSolver (DMRG::VERBOSITY::OPTION VERBOSITY=DMRG::VERBOSITY::SILENT) | |
void | set_verbosity (DMRG::VERBOSITY::OPTION VERBOSITY) |
DMRG::VERBOSITY::OPTION | get_verbosity () |
void | userSetGlobParam () |
void | userSetDynParam () |
void | userSetLocParam () |
void | edgeState (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout, qarray< Symmetry::Nq > Qtot, LANCZOS::EDGE::OPTION EDGE=LANCZOS::EDGE::GROUND, bool USE_STATE=false) |
const double & | errVar () |
const double & | errState () |
const double & | errEigval () |
Mps< Symmetry, Scalar > | create_Mps (size_t Ncells, const Eigenstate< Umps< Symmetry, Scalar > > &V, const MpHamiltonian &H, size_t x0) |
vector< Mps< Symmetry, Scalar > > | create_Mps (size_t Ncells, const Eigenstate< Umps< Symmetry, Scalar > > &V, const MpHamiltonian &H, const Mpo< Symmetry, Scalar > &O, const vector< Mpo< Symmetry, Scalar > > &Omult, double tol_OxV=2.) |
Mps< Symmetry, Scalar > | create_Mps (size_t Ncells, const Eigenstate< Umps< Symmetry, Scalar > > &V, const MpHamiltonian &H, const Mpo< Symmetry, Scalar > &O, const Mpo< Symmetry, Scalar > &Omult, double tol_OxV=2.) |
void | set_boundary (const Umps< Symmetry, Scalar > &Vin, Mps< Symmetry, Scalar > &Vout, bool LEFT=false, bool RIGHT=true) |
void | prepare (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout, qarray< Symmetry::Nq > Qtot, bool USE_STATE=false) |
void | build_cellEnv (const MpHamiltonian &H, const Eigenstate< Umps< Symmetry, Scalar > > &Vout, size_t power=1) |
double | get_err_eigval () const |
double | get_err_state () const |
double | get_err_var () const |
void | cleanup (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout) |
void | solve_linear (VMPS::DIRECTION::OPTION gauge, size_t ab, const vector< vector< Biped< Symmetry, MatrixType > > > &A, const Tripod< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > &Y_LR, const Biped< Symmetry, MatrixType > &LReigen, const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > &W, const vector< vector< qarray< Symmetry::Nq > > > &qloc, const vector< vector< qarray< Symmetry::Nq > > > &qOp, Scalar LRdotY, const Tripod< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > &LRguess, Tripod< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > &LRres) |
void | solve_linear (VMPS::DIRECTION::OPTION DIR, const vector< vector< Biped< Symmetry, MatrixType > > > &A, const Biped< Symmetry, MatrixType > &hLR, const Biped< Symmetry, MatrixType > &LReigen, const vector< vector< qarray< Symmetry::Nq > > > &qloc, Scalar hLRdotLR, Biped< Symmetry, MatrixType > &LRres) |
void | set_LanczosTolerances (double &tolLanczosEigval, double &tolLanczosState) |
void | calc_errors (const MpHamiltonian &H, const Eigenstate< Umps< Symmetry, Scalar > > &Vout, bool CALC_ERR_STATE=true, VUMPS::TWOSITE_A::OPTION option=VUMPS::TWOSITE_A::ALxAC) |
string | test_LReigen (const Eigenstate< Umps< Symmetry, Scalar > > &Vout) const |
Mps< Symmetry, Scalar > | assemble_Mps (size_t Ncells, const Umps< Symmetry, Scalar > &V, const vector< vector< Biped< Symmetry, MatrixType > > > &AL, const vector< vector< Biped< Symmetry, MatrixType > > > &AR, const vector< vector< qarray< Symmetry::Nq > > > &qloc_input, const Tripod< Symmetry, MatrixType > &L, const Tripod< Symmetry, MatrixType > &R, int x0) |
string | info () const |
string | eigeninfo () const |
double | memory (MEMUNIT memunit=GB) const |
size_t | iterations () |
void | set_log (int N_log_input, string file_e_input, string file_err_eigval_input, string file_err_var_input, string file_err_state_input) |
void | prepare_h2site (const TwoSiteHamiltonian &h2site, const vector< qarray< Symmetry::Nq > > &qloc_input, Eigenstate< Umps< Symmetry, Scalar > > &Vout, qarray< Symmetry::Nq > Qtot_input, size_t M, size_t Nqmax, bool USE_STATE=false) |
void | iteration_h2site (Eigenstate< Umps< Symmetry, Scalar > > &Vout) |
void | iteration_parallel (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout) |
void | iteration_sequential (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout) |
void | build_LR (const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > &AL, const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > &AR, const vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > &C, const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > &W, const vector< vector< qarray< Symmetry::Nq > > > &qloc, const vector< vector< qarray< Symmetry::Nq > > > &qOp, Tripod< Symmetry, MatrixType > &L, Tripod< Symmetry, MatrixType > &R) |
void | build_R (const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > &AR, const Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > &Cintercell, const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > &W, const vector< vector< qarray< Symmetry::Nq > > > &qloc, const vector< vector< qarray< Symmetry::Nq > > > &qOp, Tripod< Symmetry, MatrixType > &R) |
void | build_L (const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > &AL, const Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > &Cintercell, const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > &W, const vector< vector< qarray< Symmetry::Nq > > > &qloc, const vector< vector< qarray< Symmetry::Nq > > > &qOp, Tripod< Symmetry, MatrixType > &L) |
void | expand_basis (size_t loc, size_t DeltaD, const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout, VUMPS::TWOSITE_A::OPTION option=VUMPS::TWOSITE_A::ALxAC) |
void | expand_basis (size_t DeltaD, const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout, VUMPS::TWOSITE_A::OPTION option=VUMPS::TWOSITE_A::ALxAC) |
void | expand_basis2 (size_t DeltaD, const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout, VUMPS::TWOSITE_A::OPTION option=VUMPS::TWOSITE_A::ALxAC) |
void | calc_B2 (size_t loc, const MpHamiltonian &H, const Umps< Symmetry, Scalar > &Psi, VUMPS::TWOSITE_A::OPTION option, Biped< Symmetry, MatrixType > &B2, vector< Biped< Symmetry, MatrixType > > &NL, vector< Biped< Symmetry, MatrixType > > &NR) const |
void | calc_B2 (size_t loc, const MpHamiltonian &H, const Umps< Symmetry, Scalar > &Psi, VUMPS::TWOSITE_A::OPTION option, Biped< Symmetry, MatrixType > &B2) const |
Public Attributes | |
VUMPS::CONTROL::GLOB | GlobParam |
VUMPS::CONTROL::DYN | DynParam |
VUMPS::CONTROL::LOC | LocParam |
bool | FORCE_DO_SOMETHING = false |
size_t | N_sites |
bool | USER_SET_GLOBPARAM = false |
bool | USER_SET_DYNPARAM = false |
bool | USER_SET_LOCPARAM = false |
size_t | N_iterations =0ul |
size_t | N_iterations_without_expansion =0ul |
double | err_eigval |
double | err_eigval_old |
double | err_var |
double | err_var_old |
double | err_state =std::nan("1") |
double | err_state_old =std::nan("1") |
vector< PivumpsMatrix1< Symmetry, Scalar, Scalar > > | Heff |
vector< PivotMatrix1< Symmetry, Scalar, Scalar > > | HeffA |
vector< PivotMatrix1< Symmetry, Scalar, Scalar > > | HeffC |
vector< qarray< Symmetry::Nq > > | qloc |
TwoSiteHamiltonian | h2site |
size_t | D |
size_t | M |
size_t | dW |
size_t | dW_singlet |
vector< pair< qarray< Symmetry::Nq >, size_t > > | basis_order |
std::unordered_map< pair< qarray< Symmetry::Nq >, size_t >, size_t > | basis_order_map |
double | eL |
double | eR |
double | eoldR |
double | eoldL |
DMRG::VERBOSITY::OPTION | CHOSEN_VERBOSITY |
Tripod< Symmetry, MatrixType > | YLlast |
Tripod< Symmetry, MatrixType > | YRfrst |
Private Types | |
typedef Matrix< Scalar, Dynamic, Dynamic > | MatrixType |
typedef Matrix< complex< Scalar >, Dynamic, Dynamic > | ComplexMatrixType |
typedef Matrix< Scalar, Dynamic, 1 > | VectorType |
typedef boost::multi_array< Scalar, 4 > | TwoSiteHamiltonian |
double | Eold = std::nan("1") |
void | iteration_idmrg (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout) |
void | prepare_idmrg (const MpHamiltonian &H, Eigenstate< Umps< Symmetry, Scalar > > &Vout, qarray< Symmetry::Nq > Qtot, bool USE_STATE=false) |
size_t | N_log = 0 |
string | file_e |
string | file_err_eigval |
string | file_err_var |
string | file_err_state |
vector< double > | eL_mem |
vector< double > | eR_mem |
vector< double > | err_eigval_mem |
vector< double > | err_var_mem |
vector< double > | err_state_mem |
void | write_log (bool FORCE=false) |
|
private |
Definition at line 35 of file VumpsSolver.h.
|
private |
Definition at line 34 of file VumpsSolver.h.
|
private |
Definition at line 37 of file VumpsSolver.h.
|
private |
Definition at line 36 of file VumpsSolver.h.
|
inline |
Definition at line 41 of file VumpsSolver.h.
Mps< Symmetry, Scalar > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::assemble_Mps | ( | size_t | Ncells, |
const Umps< Symmetry, Scalar > & | V, | ||
const vector< vector< Biped< Symmetry, MatrixType > > > & | AL, | ||
const vector< vector< Biped< Symmetry, MatrixType > > > & | AR, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qloc_input, | ||
const Tripod< Symmetry, MatrixType > & | L, | ||
const Tripod< Symmetry, MatrixType > & | R, | ||
int | x0 | ||
) |
Assembles an Mps from the VUMPS solution with a heterogeneous section and infinite boundary conditions with known environments. This is a low-lever work function that is called from the higher-level create_Mps
.
Ncells | : amount of cells to generate the heterogeneous section, the total length becomes Lcell*Ncells |
V | : converged ground state to generate from |
AL | : This left-orthogonal tensor is put into the environment (can contain Jordan-Wigner string) |
AR | : This right-orthogonal tensor is put into the environment (can contain quantum number shift) |
qloc_input | : This local basis of the unit cell is put into the environment |
L | : left environment with the Hamiltonian |
R | : right environment with the Hamiltonian |
x0 | : put the pivot site here |
Definition at line 2422 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::build_cellEnv | ( | const MpHamiltonian & | H, |
const Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
size_t | power = 1 |
||
) |
Builds environments for each site of the unit cell.
Definition at line 1029 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::build_L | ( | const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > & | AL, |
const Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > & | Cintercell, | ||
const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > & | W, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qloc, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qOp, | ||
Tripod< Symmetry, MatrixType > & | L | ||
) |
Builds the environment of a unit cell.
Definition at line 790 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::build_LR | ( | const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > & | AL, |
const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > & | AR, | ||
const vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > & | C, | ||
const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > & | W, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qloc, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qOp, | ||
Tripod< Symmetry, MatrixType > & | L, | ||
Tripod< Symmetry, MatrixType > & | R | ||
) |
Builds the environment of a unit cell.
Definition at line 880 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::build_R | ( | const vector< vector< Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > > > & | AR, |
const Biped< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > & | Cintercell, | ||
const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > & | W, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qloc, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qOp, | ||
Tripod< Symmetry, MatrixType > & | R | ||
) |
Builds the environment of a unit cell.
Definition at line 835 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::calc_B2 | ( | size_t | loc, |
const MpHamiltonian & | H, | ||
const Umps< Symmetry, Scalar > & | Psi, | ||
VUMPS::TWOSITE_A::OPTION | option, | ||
Biped< Symmetry, MatrixType > & | B2 | ||
) | const |
A wrapper, if you want to discard the nullspaces when calculating B2.
Definition at line 2162 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::calc_B2 | ( | size_t | loc, |
const MpHamiltonian & | H, | ||
const Umps< Symmetry, Scalar > & | Psi, | ||
VUMPS::TWOSITE_A::OPTION | option, | ||
Biped< Symmetry, MatrixType > & | B2, | ||
vector< Biped< Symmetry, MatrixType > > & | NL, | ||
vector< Biped< Symmetry, MatrixType > > & | NR | ||
) | const |
Calculates the two-site B-tensor (from double tangent space). It is relevant for orthogonal information in the Umps as well as for the global state error.
loc | : Calculate the two-site B-tensor at sites loc and loc+1 |
H | : Mpo |
Vout | : Umps |
option | : how to calculate the two-site A-tensor (ALxAC, ARxAC or ALxCxAR) |
B2 | : The two-site B-tensor as the return value. |
NL | : The left nullspace, which is calculated during the routine and returned for later use. (Needed when performing the enrichment) |
NR | : You can guess what this is, probably. |
Definition at line 2086 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::calc_errors | ( | const MpHamiltonian & | H, |
const Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
bool | CALC_ERR_STATE = true , |
||
VUMPS::TWOSITE_A::OPTION | option = VUMPS::TWOSITE_A::ALxAC |
||
) |
Calculates the errors.
Definition at line 518 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::cleanup | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout | ||
) |
Cleans up after the iteration process.
Mps< Symmetry, Scalar > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::create_Mps | ( | size_t | Ncells, |
const Eigenstate< Umps< Symmetry, Scalar > > & | V, | ||
const MpHamiltonian & | H, | ||
const Mpo< Symmetry, Scalar > & | O, | ||
const Mpo< Symmetry, Scalar > & | Omult, | ||
double | tol_OxV = 2. |
||
) |
Variant of create_Mps without a vector, with a single MPO/MPS.
Definition at line 2339 of file VumpsSolver.h.
vector< Mps< Symmetry, Scalar > > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::create_Mps | ( | size_t | Ncells, |
const Eigenstate< Umps< Symmetry, Scalar > > & | V, | ||
const MpHamiltonian & | H, | ||
const Mpo< Symmetry, Scalar > & | O, | ||
const vector< Mpo< Symmetry, Scalar > > & | Omult, | ||
double | tol_OxV = 2. |
||
) |
Creates an Mps from the VUMPS solution with a heterogeneous section and infinite boundary conditions for a local operator. Already performs mutliplication with the operator. A possible Jordan-Wigner string on the left and a shifted quantum number on the right are absorbed into the environments.
Ncells | : amount of cells to generate the heterogeneous section, the total length becomes Lcell*Ncells |
V | : converged ground state to generate from |
H | : Hamiltonian of the VUMPS ground state, needed to recalculate environment |
O | : Operator to get the boundaries from. Should be local, with the excitation centre away from the boundaries. |
Omult | : Operator to multiply the state with |
Definition at line 2198 of file VumpsSolver.h.
Mps< Symmetry, Scalar > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::create_Mps | ( | size_t | Ncells, |
const Eigenstate< Umps< Symmetry, Scalar > > & | V, | ||
const MpHamiltonian & | H, | ||
size_t | x0 | ||
) |
Creates an Mps from the VUMPS solution with a heterogeneous section and infinite boundary conditions.
Ncells | : amount of cells to generate the heterogeneous section, the total length becomes Lcell*Ncells |
V | : converged ground state to generate from |
H | : Hamiltonian of the VUMPS ground state, needed to recalculate environment |
x0 | : Puts on this site |
Definition at line 2170 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::edgeState | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
qarray< Symmetry::Nq > | Qtot, | ||
LANCZOS::EDGE::OPTION | EDGE = LANCZOS::EDGE::GROUND , |
||
bool | USE_STATE = false |
||
) |
Calculates the highest or lowest eigenstate with an MPO (algorithm 6).
Definition at line 1784 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eigeninfo |
Prints some info about this class.
Definition at line 391 of file VumpsSolver.h.
|
inline |
Definition at line 91 of file VumpsSolver.h.
|
inline |
Definition at line 90 of file VumpsSolver.h.
|
inline |
Definition at line 89 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::expand_basis | ( | size_t | DeltaD, |
const MpHamiltonian & | H, | ||
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
VUMPS::TWOSITE_A::OPTION | option = VUMPS::TWOSITE_A::ALxAC |
||
) |
This function adds orthogonal information to the UMPS and enlarge therewith the bond dimension and the number of symmetry blocks. For information see appendix B in Zauner-Stauber et al. 2018. Just calls expand_basis() for all positions in the unit cell.
Definition at line 2075 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::expand_basis | ( | size_t | loc, |
size_t | DeltaD, | ||
const MpHamiltonian & | H, | ||
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
VUMPS::TWOSITE_A::OPTION | option = VUMPS::TWOSITE_A::ALxAC |
||
) |
This function adds orthogonal information to the UMPS in the unit cell at site loc and enlarge therewith the bond dimension and the number of symmetry blocks. For information see appendix B in Zauner-Stauber et al. 2018.
Definition at line 1981 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::expand_basis2 | ( | size_t | DeltaD, |
const MpHamiltonian & | H, | ||
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
VUMPS::TWOSITE_A::OPTION | option = VUMPS::TWOSITE_A::ALxAC |
||
) |
This function adds orthogonal information to the UMPS in the unit cell at site loc and enlarge therewith the bond dimension and the number of symmetry blocks. For information see appendix B in Zauner-Stauber et al. 2018.
Definition at line 2518 of file VumpsSolver.h.
|
inline |
Definition at line 131 of file VumpsSolver.h.
|
inline |
Definition at line 132 of file VumpsSolver.h.
|
inline |
Definition at line 133 of file VumpsSolver.h.
|
inline |
Returnx the verbosity level.
Definition at line 49 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::info |
Prints some info about this class.
Definition at line 380 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::iteration_h2site | ( | Eigenstate< Umps< Symmetry, Scalar > > & | Vout | ) |
Performs an iteration with 1-site unit cell. Used with an explicit 2-site Hamiltonian.
Definition at line 1539 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::iteration_idmrg | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout | ||
) |
Performs an IDMRG iteration with a 2-site unit cell.
Definition at line 1649 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::iteration_parallel | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout | ||
) |
Performs an iteration with an n-site unit cell (in parallel, algorithm 3). Used with an MPO.
Definition at line 1100 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::iteration_sequential | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout | ||
) |
Performs an iteration with an n-site unit cell (sequentially, algorithm 4). Used with an MPO.
Definition at line 1351 of file VumpsSolver.h.
|
inline |
Prints some info about this class.
Definition at line 70 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::memory | ( | MEMUNIT | memunit = GB | ) | const |
Calculates the (theoretically) allocated memory (note: by default in GB).
Definition at line 407 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::prepare | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
qarray< Symmetry::Nq > | Qtot, | ||
bool | USE_STATE = false |
||
) |
Prepares the class setting up the environments. Used with an Mpo.
Definition at line 608 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::prepare_h2site | ( | const TwoSiteHamiltonian & | h2site, |
const vector< qarray< Symmetry::Nq > > & | qloc_input, | ||
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
qarray< Symmetry::Nq > | Qtot_input, | ||
size_t | M, | ||
size_t | Nqmax, | ||
bool | USE_STATE = false |
||
) |
Prepares the class, setting up the environments. Used with an explicit 2-site Hamiltonian.
Definition at line 565 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::prepare_idmrg | ( | const MpHamiltonian & | H, |
Eigenstate< Umps< Symmetry, Scalar > > & | Vout, | ||
qarray< Symmetry::Nq > | Qtot, | ||
bool | USE_STATE = false |
||
) |
Prepares the class, setting up the environments for IDMRG.
Definition at line 744 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::set_boundary | ( | const Umps< Symmetry, Scalar > & | Vin, |
Mps< Symmetry, Scalar > & | Vout, | ||
bool | LEFT = false , |
||
bool | RIGHT = true |
||
) |
Definition at line 2471 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::set_LanczosTolerances | ( | double & | tolLanczosEigval, |
double & | tolLanczosState | ||
) |
Sets the Lanczos tolerances adaptively, depending on the current errors.
Definition at line 492 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::set_log | ( | int | N_log_input, |
string | file_e_input, | ||
string | file_err_eigval_input, | ||
string | file_err_var_input, | ||
string | file_err_state_input | ||
) |
Setup a logfile of the iterations.
N_log_input | : save the log every N_log_input half-sweeps |
file_e_input | : file for the ground-state energy in the format [min(eL,eR), eL, eR] |
file_err_eigval_input | : file for the energy error |
file_err_var_input | : file for the variational error |
file_err_state_input | : file for the global state error |
Definition at line 426 of file VumpsSolver.h.
|
inline |
Resets the verbosity level.
Definition at line 46 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::solve_linear | ( | VMPS::DIRECTION::OPTION | DIR, |
const vector< vector< Biped< Symmetry, MatrixType > > > & | A, | ||
const Biped< Symmetry, MatrixType > & | hLR, | ||
const Biped< Symmetry, MatrixType > & | LReigen, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qloc, | ||
Scalar | hLRdotLR, | ||
Biped< Symmetry, MatrixType > & | LRres | ||
) |
Solves the linear system (eq. 15) using GMRES.
DIR | : L or R |
A | : contracted A-tensor of the cell |
hLR | : (h_L|, |h_R) for eq. 15 |
LReigen | : (L| or |R) |
qloc | : local basis |
hLRdotLR | : (h_L|R), (L|h_R) for eq. 15 |
LRres | : resulting (H_L| or |H_R) |
Definition at line 1943 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::solve_linear | ( | VMPS::DIRECTION::OPTION | gauge, |
size_t | ab, | ||
const vector< vector< Biped< Symmetry, MatrixType > > > & | A, | ||
const Tripod< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > & | Y_LR, | ||
const Biped< Symmetry, MatrixType > & | LReigen, | ||
const vector< vector< vector< vector< Biped< Symmetry, SparseMatrix< Scalar > > > > > > & | W, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qloc, | ||
const vector< vector< qarray< Symmetry::Nq > > > & | qOp, | ||
Scalar | LRdotY, | ||
const Tripod< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > & | LRguess, | ||
Tripod< Symmetry, Matrix< Scalar, Dynamic, Dynamic > > & | LRres | ||
) |
Solves the linear system (eq. C25ab) using GMRES.
gauge | : L or R |
ab | : fixed index a (rows) or b (cols) |
A | : contracted A-tensor of the cell |
Y_LR | : |Y_Ra), (Y_La| for eq. C25ab |
LReigen | : (L| or |R) |
W | : Mpo tensor for the transfer matrix |
qloc | : local basis |
qOp | : operator basis |
LRdotY | : (Y_La|R), (L|Y_Ra) for eq. C25ab |
LRguess | : the starting guess for the linear solver |
LRres | : resulting (H_L| or |H_R) |
Definition at line 1908 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::test_LReigen | ( | const Eigenstate< Umps< Symmetry, Scalar > > & | Vout | ) | const |
Tests if and are the eigenvectors of the transfer matrices.
Definition at line 1763 of file VumpsSolver.h.
|
inline |
Definition at line 53 of file VumpsSolver.h.
|
inline |
Definition at line 52 of file VumpsSolver.h.
|
inline |
Definition at line 54 of file VumpsSolver.h.
void VumpsSolver< Symmetry, MpHamiltonian, Scalar >::write_log | ( | bool | FORCE = false | ) |
Function to write out the logfiles.
FORCE | : if true , forced write without checking any conditions |
Definition at line 442 of file VumpsSolver.h.
vector<pair<qarray<Symmetry::Nq>,size_t> > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::basis_order |
std::unordered_map<pair<qarray<Symmetry::Nq>,size_t>,size_t> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::basis_order_map |
Definition at line 274 of file VumpsSolver.h.
DMRG::VERBOSITY::OPTION VumpsSolver< Symmetry, MpHamiltonian, Scalar >::CHOSEN_VERBOSITY |
control of verbosity and algorithms
Definition at line 322 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::D |
bond dimension per subspace, bond dimension per site, Mpo bond dimension, Mpo bond dimension in the singlet sector
Definition at line 270 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::dW |
Definition at line 270 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::dW_singlet |
Definition at line 270 of file VumpsSolver.h.
VUMPS::CONTROL::DYN VumpsSolver< Symmetry, MpHamiltonian, Scalar >::DynParam |
Definition at line 57 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eL |
left and right error (eq. 18) and old errors from previous half-sweep
Definition at line 277 of file VumpsSolver.h.
vector<double> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eL_mem |
log data
Definition at line 348 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::Eold = std::nan("1") |
old energy for comparison in IDMRG.
Definition at line 170 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eoldL |
Definition at line 277 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eoldR |
Definition at line 277 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eR |
Definition at line 277 of file VumpsSolver.h.
vector<double> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::eR_mem |
Save log every N_log
optimization steps.
Definition at line 348 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_eigval |
errors
Definition at line 252 of file VumpsSolver.h.
vector<double> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_eigval_mem |
Save log every N_log
optimization steps.
Definition at line 348 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_eigval_old |
Definition at line 252 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_state =std::nan("1") |
Definition at line 252 of file VumpsSolver.h.
vector<double> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_state_mem |
Save log every N_log
optimization steps.
Definition at line 348 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_state_old =std::nan("1") |
Definition at line 252 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_var |
Definition at line 252 of file VumpsSolver.h.
vector<double> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_var_mem |
Save log every N_log
optimization steps.
Definition at line 348 of file VumpsSolver.h.
double VumpsSolver< Symmetry, MpHamiltonian, Scalar >::err_var_old |
Definition at line 252 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::file_e |
log filenames
Definition at line 345 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::file_err_eigval |
Save log every N_log
optimization steps.
Definition at line 345 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::file_err_state |
Save log every N_log
optimization steps.
Definition at line 345 of file VumpsSolver.h.
string VumpsSolver< Symmetry, MpHamiltonian, Scalar >::file_err_var |
Save log every N_log
optimization steps.
Definition at line 345 of file VumpsSolver.h.
bool VumpsSolver< Symmetry, MpHamiltonian, Scalar >::FORCE_DO_SOMETHING = false |
Definition at line 93 of file VumpsSolver.h.
VUMPS::CONTROL::GLOB VumpsSolver< Symmetry, MpHamiltonian, Scalar >::GlobParam |
Definition at line 56 of file VumpsSolver.h.
TwoSiteHamiltonian VumpsSolver< Symmetry, MpHamiltonian, Scalar >::h2site |
stored 2-site Hamiltonian
Definition at line 267 of file VumpsSolver.h.
vector<PivumpsMatrix1<Symmetry,Scalar,Scalar> > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::Heff |
environment for the 2-site Hamiltonian version
Definition at line 255 of file VumpsSolver.h.
vector<PivotMatrix1<Symmetry,Scalar,Scalar> > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::HeffA |
environment of AL
and AR
for the Mpo version
Definition at line 258 of file VumpsSolver.h.
vector<PivotMatrix1<Symmetry,Scalar,Scalar> > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::HeffC |
environment of C
for the Mpo version
Definition at line 261 of file VumpsSolver.h.
VUMPS::CONTROL::LOC VumpsSolver< Symmetry, MpHamiltonian, Scalar >::LocParam |
Definition at line 58 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::M |
Definition at line 270 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::N_iterations =0ul |
keeping track of iterations
Definition at line 249 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::N_iterations_without_expansion =0ul |
Definition at line 249 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::N_log = 0 |
Save log every N_log
optimization steps.
Definition at line 342 of file VumpsSolver.h.
size_t VumpsSolver< Symmetry, MpHamiltonian, Scalar >::N_sites |
chain length
Definition at line 242 of file VumpsSolver.h.
vector<qarray<Symmetry::Nq> > VumpsSolver< Symmetry, MpHamiltonian, Scalar >::qloc |
local base
Definition at line 264 of file VumpsSolver.h.
bool VumpsSolver< Symmetry, MpHamiltonian, Scalar >::USER_SET_DYNPARAM = false |
Definition at line 245 of file VumpsSolver.h.
bool VumpsSolver< Symmetry, MpHamiltonian, Scalar >::USER_SET_GLOBPARAM = false |
Definition at line 244 of file VumpsSolver.h.
bool VumpsSolver< Symmetry, MpHamiltonian, Scalar >::USER_SET_LOCPARAM = false |
Definition at line 246 of file VumpsSolver.h.
Tripod<Symmetry,MatrixType> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::YLlast |
saved , see eq. (C26), (C27)
Definition at line 332 of file VumpsSolver.h.
Tripod<Symmetry,MatrixType> VumpsSolver< Symmetry, MpHamiltonian, Scalar >::YRfrst |
saved , see eq. (C26), (C27)
Definition at line 335 of file VumpsSolver.h.