VMPS++
Loading...
Searching...
No Matches
TDVPPropagator.h
Go to the documentation of this file.
1#ifndef TDVP_PROPAGATOR
2#define TDVP_PROPAGATOR
3
4#include "Stopwatch.h" // from TOOLS
5#include "LanczosPropagator.h" // from ALGS
6
7#include "pivot/DmrgPivotMatrix0.h"
8#include "pivot/DmrgPivotMatrix2.h"
9//include "tensors/DmrgContractions.h"
10#include "DmrgJanitor.h" // for turnaround
11#include <gsl/gsl_math.h> // for GSL_SIGN
12
13template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
14class TDVPPropagator
15{
16public:
17
18 TDVPPropagator(){};
19 TDVPPropagator (const Hamiltonian &H, VectorType &Vinout);
20
21 string info() const;
22 double memory (MEMUNIT memunit=GB) const;
23 double overhead (MEMUNIT memunit=GB) const;
24
25 void t_step (const Hamiltonian &H, const VectorType &Vin, VectorType &Vout, TimeScalar dt, int N_stages=1, double tol_Lanczos=1e-8);
26
27 void t_step (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, int N_stages=1, double tol_Lanczos=1e-8);
28 void t_step0 (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, int N_stages=1, double tol_Lanczos=1e-8);
29
30 void t_step_adaptive (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, const vector<bool> &TWO_STEP_AT, int N_stages=1, double tol_Lanczos=1e-8);
31
32 inline VectorXd get_deltaE() const {return deltaE;};
33 inline double get_t_tot() const {return t_tot;};
34 inline vector<size_t> get_dimK2_log() const {return dimK2_log;};
35 inline vector<size_t> get_dimK1_log() const {return dimK1_log;};
36 inline vector<size_t> get_dimK0_log() const {return dimK0_log;};
37
38private:
39
40 void t_step_pivot (double x, const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, double tol_Lanczos=1e-8);
41 void t0_step_pivot (bool BACK, double x, const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, double tol_Lanczos=1e-8, bool TURN_FIRST=true);
42
43 void test_edge_eigenvector (const PivotVector<Symmetry,TimeScalar> &Asingle);
44 VectorXd deltaE;
45 VectorXd dimKlog;
46
47 vector<PivotMatrix1<Symmetry,TimeScalar,MpoScalar> > Heff;
48 PivotMatrix1<Symmetry,TimeScalar,MpoScalar> HeffLast;
49 PivotMatrix1<Symmetry,TimeScalar,MpoScalar> HeffFrst;
50
51 double x (int alg, size_t l, int N_stages);
52
53 void set_blocks (const Hamiltonian &H, VectorType &Vinout);
54
55 size_t N_sites;
56 int pivot;
57 DMRG::DIRECTION::OPTION CURRENT_DIRECTION;
58
59 void build_L (const Hamiltonian &H, const VectorType &Vinout, int loc);
60 void build_R (const Hamiltonian &H, const VectorType &Vinout, int loc);
61
62 double dist_max = 0.;
63 double dimK_max = 0.;
64 vector<size_t> dimK0_log;
65 vector<size_t> dimK1_log;
66 vector<size_t> dimK2_log;
67 int N_stages_last = 0;
68
69 double t_0site = 0;
70 double t_1site = 0;
71 double t_2site = 0;
72 double t_ohead = 0; // precalc
73 double t_contr = 0; // contract & sweep
74 double t_tot = 0; // full time step
75
76 TimeScalar last_dt;
77};
78
79template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
80string TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
81info() const
82{
83 stringstream ss;
84 ss << "TDVPPropagator: ";
85 ss << " dt=" << last_dt;
86 ss << ", max(dist)=" << dist_max << ", ";
87// ss << "max(dimK)=" << dimK_max << ", ";
88 if (dimK2_log.size() > 0)
89 {
90 ss << "dimK2=" << *min_element(dimK2_log.begin(), dimK2_log.end()) << "~" << *max_element(dimK2_log.begin(), dimK2_log.end()) << ", ";
91 }
92 if (dimK1_log.size() > 0)
93 {
94 ss << "dimK1=" << *min_element(dimK1_log.begin(), dimK1_log.end()) << "~" << *max_element(dimK1_log.begin(), dimK1_log.end()) << ", ";
95 }
96 if (dimK0_log.size() > 0)
97 {
98 ss << "dimK0=" << *min_element(dimK0_log.begin(), dimK0_log.end()) << "~" << *max_element(dimK0_log.begin(), dimK0_log.end()) << ", ";
99 }
100// ss << "N_stages=" << N_stages_last << ", ";
101// ss << "mem=" << round(memory(GB),3) << "GB, ";
102 ss << "δE@edge: L=" << deltaE(0) << ", R=" << deltaE(deltaE.rows()-1) << ", ";
103 ss << "overhead=" << round(overhead(MB),3) << "MB, ";
104 ss << "t[s]=" << t_tot << ", "
105 << "t0=" << round(t_0site/t_tot*100.,0) << "%, "
106 << "t1=" << round(t_1site/t_tot*100.) << "%, "
107 << "t2=" << round(t_2site/t_tot*100.) << "%, "
108 << "t_ohead=" << round(t_ohead/t_tot*100.) << "%, "
109 << "t_contr=" << round(t_contr/t_tot*100.) << "%";
110
111 return ss.str();
112}
113
114template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
115double TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
116memory (MEMUNIT memunit) const
117{
118 double res = 0.;
119// for (size_t l=0; l<N_sites; ++l)
120// {
121// res += Heff[l].L.memory(memunit);
122// res += Heff[l].R.memory(memunit);
123// for (size_t s1=0; s1<Heff[l].W.size(); ++s1)
124// for (size_t s2=0; s2<Heff[l].W[s1].size(); ++s2)
125// {
126// res += calc_memory(Heff[l].W[s1][s2],memunit);
127// }
128// }
129 return res;
130}
131
132template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
133double TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
134overhead (MEMUNIT memunit) const
135{
136 double res = 0.;
137// for (size_t l=0; l<N_sites; ++l)
138// {
139// res += Heff[l].L.overhead(memunit);
140// res += Heff[l].R.overhead(memunit);
141// res += 2. * calc_memory<size_t>(Heff[l].qlhs.size(),memunit);
142// res += 4. * calc_memory<size_t>(Heff[l].qrhs.size(),memunit);
143// }
144 return res;
145}
146
147template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
148TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
149TDVPPropagator (const Hamiltonian &H, VectorType &Vinout)
150{
151 set_blocks(H,Vinout);
152 assert(H.HAS_TWO_SITE_DATA() and "You need to call H.precalc_TwoSiteData() before dynamics!");
153}
154
155template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
156void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
157set_blocks (const Hamiltonian &H, VectorType &Vinout)
158{
159 N_sites = H.length();
160
161 // set edges
162 Heff.clear();
163 Heff.resize(N_sites);
164// Heff[0].L.setVacuum();
165// Heff[N_sites-1].R.setTarget(qarray3<Symmetry::Nq>{Vinout.Qtarget(), Vinout.Qtarget(), Symmetry::qvacuum()});
166 for (int l=0; l<N_sites; ++l) Heff[l].Terms.resize(1);
167 Heff[0].Terms[0].L = Vinout.get_boundaryTensor(DMRG::DIRECTION::LEFT);
168 Heff[N_sites-1].Terms[0].R = Vinout.get_boundaryTensor(DMRG::DIRECTION::RIGHT);
169
170 for (size_t l=0; l<N_sites; ++l)
171 {
172 Heff[l].Terms[0].W = H.W[l];
173 }
174
175 //----- workaround sign bug >>>>>
176 auto Vtmp = Vinout;
177 //<<<<< workaround sign bug -----
178
179 // initial sweep, right-to-left:
180 for (size_t l=N_sites-1; l>0; --l)
181 {
182 Vinout.sweepStep(DMRG::DIRECTION::LEFT, l, DMRG::BROOM::QR);
183 build_R(H,Vinout,l-1);
184// cout << "l=" << l << ", dot=" << Vtmp.dot(Vinout) << endl;
185 }
186 CURRENT_DIRECTION = DMRG::DIRECTION::RIGHT;
187 pivot = 0;
188
189 //----- workaround sign bug >>>>>
190 double overlap = isReal(Vtmp.dot(Vinout));
191 if (GSL_SIGN(overlap) == -1)
192 {
193 Vinout *= -1.;
194 lout << termcolor::yellow << "dot=" << overlap << ", minus overlap, compensating..." << termcolor::reset << endl;
195 }
196 //<<<<< workaround sign bug -----
197
198 deltaE.resize(N_sites);
199
200 // initial sweep, left-to-right:
201// for (size_t l=0; l<N_sites-1; ++l)
202// {
203// Vinout.sweepStep(DMRG::DIRECTION::RIGHT, l, DMRG::BROOM::QR);
204// build_L(H,Vinout,l+1);
205// }
206// CURRENT_DIRECTION = DMRG::DIRECTION::LEFT;
207// pivot = N_sites-1;
208}
209
210template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
211double TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
212x (int alg, size_t l, int N_stages)
213{
214 double out;
215 int N_updates = (alg==2)? N_sites-1 : N_sites;
216
217 if (N_stages == 1)
218 {
219 out = 0.5;
220 }
221 else if (N_stages == 2)
222 {
223 if (l<N_updates) {out = 0.25;}
224 else if (l>=N_updates and l<2*N_updates) {out = 0.25;}
225 else if (l>=2*N_updates and l<3*N_updates) {out = 0.25;}
226 else if (l>=3*N_updates) {out = 0.25;}
227 }
228 else if (N_stages == 3)
229 {
230 double tripleJump1 = 1./(2.-pow(2.,1./3.));
231 double tripleJump2 = 1.-2.*tripleJump1;
232
233 if (l< 2*N_updates) {out = 0.5*tripleJump1;}
234 else if (l>=2*N_updates and l<4*N_updates) {out = 0.5*tripleJump2;}
235 else if (l>=4*N_updates) {out = 0.5*tripleJump1;}
236 }
237 return out;
238}
239
240template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
241void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
242t_step (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, int N_stages, double tol_Lanczos)
243{
244 last_dt = dt;
245 assert(N_stages==1 or N_stages==2 or N_stages==3 and "Only N_stages=1,2,3 implemented for TDVPPropagator::t_step!");
246 dist_max = 0.;
247 dimK_max = 0;
248 N_stages_last = N_stages;
249
250 t_0site = 0;
251 t_1site = 0;
252 t_2site = 0;
253 t_ohead = 0;
254 t_contr = 0;
255 t_tot = 0;
256
257 Stopwatch<> Wtot;
258
259 dimK2_log.clear();
260 dimK1_log.clear();
261 dimK0_log.clear();
262
263 for (size_t l=0; l<2*N_stages*(N_sites-1); ++l)
264 {
265 Stopwatch<> Chronos;
266 turnaround(pivot, N_sites, CURRENT_DIRECTION);
267
268 // 2-site propagation
269 size_t loc1 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot : pivot-1;
270 size_t loc2 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot+1 : pivot;
271// cout << "2site between: " << loc1 << "," << loc2 << ", pivot=" << pivot << endl;
272
273 Stopwatch<> Wc;
274 PivotVector<Symmetry,TimeScalar> Apair(Vinout.A[loc1], Vinout.locBasis(loc1),
275 Vinout.A[loc2], Vinout.locBasis(loc2),
276 Vinout.QoutTop[loc1], Vinout.QoutBot[loc1]);
277 t_contr += Wc.time(SECONDS);
278 PivotMatrix2<Symmetry,TimeScalar,MpoScalar> Heff2(Heff[loc1].Terms[0].L, Heff[loc2].Terms[0].R,
279 H.W_at(loc1), H.W_at(loc2),
280 H.locBasis(loc1), H.locBasis(loc2),
281 H.opBasis(loc1), H.opBasis(loc2));
282
283 Stopwatch<> Woh2;
284 // WRONG FOR MULTISITE TERMS:
285// precalc_blockStructure (Heff[loc1].L, Apair.data, Heff2.W12, Heff2.W34, Apair.data, Heff[loc2].R,
286// H.locBasis(loc1), H.locBasis(loc2), H.opBasis(loc1), H.opBasis(loc2),
287// H.TSD[loc1],
288// Heff2.qlhs, Heff2.qrhs, Heff2.factor_cgcs);
289 precalc_blockStructure (Heff2.Terms[0].L, Apair.data, Heff2.Terms[0].W12, Heff2.Terms[0].W34, Apair.data, Heff2.Terms[0].R,
290 H.locBasis(loc1), H.locBasis(loc2), H.opBasis(loc1), H.opBasis(loc2),
291 Heff2.qlhs, Heff2.qrhs, Heff2.factor_cgcs);
292 t_ohead += Woh2.time(SECONDS);
293
294// LanczosPropagator<PivotMatrix2<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar> > Lutz2(tol_Lanczos);
295 LanczosPropagator<PivotMatrix2<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz2(tol_Lanczos);
296 Stopwatch<> W2;
297// Lutz2.t_step(Heff2, Apair, -x(2,l,N_stages)*dt.imag()); // 2-site algorithm
298 Lutz2.t_step(Heff2, Apair, x(2,l,N_stages)*dt); // 2-site algorithm
299// cout << Lutz2.info() << endl;
300 t_2site += W2.time(SECONDS);
301
302 dimK2_log.push_back(Lutz2.get_dimK());
303 if (Lutz2.get_dist() > dist_max) {dist_max = Lutz2.get_dist();}
304 if (Lutz2.get_dimK() > dimK_max) {dimK_max = Lutz2.get_dimK();}
305
306 Stopwatch<> Ws;
307 Vinout.sweepStep2(CURRENT_DIRECTION, loc1, Apair.data);
308 t_contr += Ws.time(SECONDS);
309 (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)? build_L(H,Vinout,loc2) : build_R(H,Vinout,loc1);
310 pivot = Vinout.get_pivot();
311
312 // 1-site propagation
313 if ((CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT and pivot != N_sites-1) or
314 (CURRENT_DIRECTION==DMRG::DIRECTION::LEFT and pivot != 0))
315 {
316 Stopwatch<> Woh1;
317 precalc_blockStructure (Heff[pivot].Terms[0].L, Vinout.A[pivot], Heff[pivot].Terms[0].W, Vinout.A[pivot], Heff[pivot].Terms[0].R,
318 H.locBasis(pivot), H.opBasis(pivot),
319 Heff[pivot].qlhs, Heff[pivot].qrhs, Heff[pivot].factor_cgcs);
320 t_ohead += Woh1.time(SECONDS);
321
322 PivotVector<Symmetry,TimeScalar> Asingle(Vinout.A[pivot]);
323
324// lout << "1site at: " << pivot << endl;
325// LanczosPropagator<PivotMatrix1<Symmetry,TimeScalar,MpoScalar>, PivotVector<Symmetry,TimeScalar> > Lutz(tol_Lanczos);
326 LanczosPropagator<PivotMatrix1<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz(tol_Lanczos);
327 Stopwatch<> W1;
328// Lutz.t_step(Heff[pivot], Asingle, +x(2,l,N_stages)*dt.imag()); // 1-site algorithm
329 Lutz.t_step(Heff[pivot], Asingle, -x(2,l,N_stages)*dt); // 1-site algorithm
330// cout << Lutz.info() << endl;
331 t_1site += W1.time(SECONDS);
332
333 dimK1_log.push_back(Lutz2.get_dimK());
334 if (Lutz.get_dist() > dist_max) {dist_max = Lutz2.get_dist();}
335 if (Lutz.get_dimK() > dimK_max) {dimK_max = Lutz2.get_dimK();}
336
337 Vinout.A[pivot] = Asingle.data;
338 }
339 }
340
341// lout << "final pivot=" << pivot << endl;
342
343 t_tot = Wtot.time(SECONDS);
344
345// double norm_Psi_t = Vref.squaredNorm();
346// double norm_Psi_dt = Vinout.squaredNorm();
347// double overlap = dot(Vref,Vinout).real();
348// double eta = (norm_Psi_t+norm_Psi_dt-2.*overlap)/pow(dt.imag(),2);
349//
350// double avgH = avg(Vinout,H,Vinout).real();
351// double avgHsq = avg(Vinout,H,Vinout,true).real();
352//
353// cout << "error: " << pow(avgHsq-avgH,2)-pow(eta,2) << "\t" << avgHsq-avgH << "\t" << eta << endl;
354}
355
356//template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
357//void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
358//t_step3 (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, int N_stages, double tol_Lanczos)
359//{
360// assert(N_stages==1 or N_stages==2 or N_stages==3 and "Only N_stages=1,2,3 implemented for TDVPPropagator::t_step!");
361// dist_max = 0.;
362// dimK_max = 0;
363// N_stages_last = N_stages;
364//
365// t_0site = 0;
366// t_1site = 0;
367// t_2site = 0;
368// t_ohead = 0;
369// t_contr = 0;
370// t_tot = 0;
371//
372// Stopwatch<> Wtot;
373//
374// dimK2_log.clear();
375// dimK1_log.clear();
376// dimK0_log.clear();
377//
378// for (size_t l=0; l<2*N_stages*(N_sites-2); ++l)
379// {
380// Stopwatch<> Chronos;
381// turnaround(pivot, N_sites, CURRENT_DIRECTION);
382//
383// // 3-site propagation
384// size_t loc1 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot : pivot-2;
385// size_t loc2 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot+1 : pivot-1;
386// size_t loc3 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot+2 : pivot;
387//
388// Stopwatch<> Wc;
389// PivotVector<Symmetry,TimeScalar> Apair(Vinout.A[loc1], Vinout.locBasis(loc1),
390// Vinout.A[loc2], Vinout.locBasis(loc2),
391// Vinout.QoutTop[loc1], Vinout.QoutBot[loc1]);
392// t_contr += Wc.time(SECONDS);
393// PivotMatrix2<Symmetry,TimeScalar,MpoScalar> Heff2(Heff[loc1].L, Heff[loc2].R,
394// H.W_at(loc1), H.W_at(loc2),
395// H.locBasis(loc1), H.locBasis(loc2),
396// H.opBasis(loc1), H.opBasis(loc2));
397//
398// Stopwatch<> Woh2;
399// WRONG FOR MULTISITE TERMS:
400// precalc_blockStructure (Heff[loc1].L, Apair.data, Heff2.W12, Heff2.W34, Apair.data, Heff[loc2].R,
401// H.locBasis(loc1), H.locBasis(loc2), H.opBasis(loc1), H.opBasis(loc2),
402// H.TSD[loc1],
403// Heff2.qlhs, Heff2.qrhs, Heff2.factor_cgcs);
404// t_ohead += Woh2.time(SECONDS);
405//
406// LanczosPropagator<PivotMatrix2<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz2(tol_Lanczos);
407// Stopwatch<> W2;
408// Lutz2.t_step(Heff2, Apair, x(2,l,N_stages)*dt); // 2-site algorithm
409// t_2site += W2.time(SECONDS);
410//
411// dimK2_log.push_back(Lutz2.get_dimK());
412// if (Lutz2.get_dist() > dist_max) {dist_max = Lutz2.get_dist();}
413// if (Lutz2.get_dimK() > dimK_max) {dimK_max = Lutz2.get_dimK();}
414//
415// Stopwatch<> Ws;
416// Vinout.sweepStep2(CURRENT_DIRECTION, loc1, Apair.data);
417// t_contr += Ws.time(SECONDS);
418// (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)? build_L(H,Vinout,loc2) : build_R(H,Vinout,loc1);
419// pivot = Vinout.get_pivot();
420//
421// // 1-site propagation
422// if ((CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT and pivot != N_sites-1) or
423// (CURRENT_DIRECTION==DMRG::DIRECTION::LEFT and pivot != 0))
424// {
425// Stopwatch<> Woh1;
426// precalc_blockStructure (Heff[pivot].L, Vinout.A[pivot], Heff[pivot].W, Vinout.A[pivot], Heff[pivot].R,
427// H.locBasis(pivot), H.opBasis(pivot),
428// Heff[pivot].qlhs, Heff[pivot].qrhs, Heff[pivot].factor_cgcs);
429// t_ohead += Woh1.time(SECONDS);
430//
431// PivotVector<Symmetry,TimeScalar> Asingle(Vinout.A[pivot]);
432//
433// LanczosPropagator<PivotMatrix1<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz(tol_Lanczos);
434// Stopwatch<> W1;
435// Lutz.t_step(Heff[pivot], Asingle, -x(2,l,N_stages)*dt); // 1-site algorithm
436// t_1site += W1.time(SECONDS);
437//
438// dimK1_log.push_back(Lutz2.get_dimK());
439// if (Lutz.get_dist() > dist_max) {dist_max = Lutz2.get_dist();}
440// if (Lutz.get_dimK() > dimK_max) {dimK_max = Lutz2.get_dimK();}
441//
442// Vinout.A[pivot] = Asingle.data;
443// }
444// }
445//
446// t_tot = Wtot.time(SECONDS);
447//}
448
449template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
450void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
451t_step0 (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, int N_stages, double tol_Lanczos)
452{
453 last_dt = dt;
454 dist_max = 0.;
455 dimK_max = 0.;
456 N_stages_last = N_stages;
457
458 t_0site = 0;
459 t_1site = 0;
460 t_2site = 0;
461 t_ohead = 0;
462 t_contr = 0;
463 t_tot = 0;
464
465 Stopwatch<> Wtot;
466
467 dimK2_log.clear();
468 dimK1_log.clear();
469 dimK0_log.clear();
470
471// VectorType Vref = Vinout;
472
473 for (size_t l=0; l<2*N_stages*N_sites; ++l)
474 {
475 turnaround(pivot, N_sites, CURRENT_DIRECTION);
476
477 // 1-site propagation
478 PivotVector<Symmetry,TimeScalar> Asingle(Vinout.A[pivot]);
479 Stopwatch<> Woh1;
480 precalc_blockStructure (Heff[pivot].Terms[0].L, Vinout.A[pivot], Heff[pivot].Terms[0].W, Vinout.A[pivot], Heff[pivot].Terms[0].R,
481 H.locBasis(pivot), H.opBasis(pivot), Heff[pivot].qlhs, Heff[pivot].qrhs, Heff[pivot].factor_cgcs);
482 t_ohead += Woh1.time(SECONDS);
483
484// lout << "1site at: " << pivot << endl;
485 LanczosPropagator<PivotMatrix1<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz(tol_Lanczos);
486
487 Stopwatch<> W1;
488// Lutz.t_step(Heff[pivot], Asingle, -x(1,l,N_stages)*dt.imag()); // 1-site algorithm
489 Lutz.t_step(Heff[pivot], Asingle, x(1,l,N_stages)*dt); // 1-site algorithm
490 t_1site += W1.time(SECONDS);
491
492 dimK1_log.push_back(Lutz.get_dimK());
493 if (Lutz.get_dist() > dist_max) {dist_max = Lutz.get_dist();}
494 if (Lutz.get_dimK() > dimK_max) {dimK_max = Lutz.get_dimK();}
495 Vinout.A[pivot] = Asingle.data;
496
497 // 0-site propagation
498 if ((l+1)%N_sites != 0)
499 {
500 PivotVector<Symmetry,TimeScalar> Azero;
501 int old_pivot = pivot;
502 (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)? Vinout.rightSplitStep(pivot,Azero.data[0]) : Vinout.leftSplitStep(pivot,Azero.data[0]);
503 pivot = Vinout.get_pivot();
504 (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)? build_L(H,Vinout,pivot) : build_R(H,Vinout,pivot);
505
506// if (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)
507// {
508// lout << "0site between " << old_pivot << "," << old_pivot+1 << endl;
509// }
510// else
511// {
512// lout << "0site between " << old_pivot-1 << "," << old_pivot << endl;
513// }
514 LanczosPropagator<PivotMatrix0<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz0(tol_Lanczos);
515
516 PivotMatrix0<Symmetry,TimeScalar,MpoScalar> Heff0;
517 (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)?
518 Heff0 = PivotMatrix0<Symmetry,TimeScalar,MpoScalar>(Heff[old_pivot+1].Terms[0].L, Heff[old_pivot].Terms[0].R):
519 Heff0 = PivotMatrix0<Symmetry,TimeScalar,MpoScalar>(Heff[old_pivot].Terms[0].L, Heff[old_pivot-1].Terms[0].R);
520
521 Stopwatch<> W0;
522// Lutz0.t_step(Heff0, Azero, +x(1,l,N_stages)*dt.imag()); // 0-site algorithm
523 Lutz0.t_step(Heff0, Azero, -x(1,l,N_stages)*dt); // 0-site algorithm
524 t_0site += W0.time(SECONDS);
525
526 dimK0_log.push_back(Lutz0.get_dimK());
527 if (Lutz0.get_dist() > dist_max) {dist_max = Lutz0.get_dist();}
528 if (Lutz0.get_dimK() > dimK_max) {dimK_max = Lutz0.get_dimK();}
529
530 Vinout.absorb(pivot, CURRENT_DIRECTION, Azero.data[0]);
531 }
532 }
533
534// lout << "final pivot=" << pivot << endl;
535
536 t_tot = Wtot.time(SECONDS);
537
538// double norm_Psi_t = Vref.squaredNorm();
539// double norm_Psi_dt = Vinout.squaredNorm();
540// double overlap = dot(Vref,Vinout).real();
541// double eta = (norm_Psi_t+norm_Psi_dt-2.*overlap)/dt.imag();
542//
543// double avgH = avg(Vinout,H,Vinout).real();
544// double avgHsq = avg(Vinout,H,Vinout,true).real();
545//
546// cout << "error: " << pow(avgH-avgHsq,2)-pow(eta,2) << endl;
547}
548
549template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
550void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
551t_step (const Hamiltonian &H, const VectorType &Vin, VectorType &Vout, TimeScalar dt, int N_stages, double tol_Lanczos)
552{
553 Vout = Vin;
554 set_blocks(H,Vout);
555 t_step(H,Vout,dt,N_stages,tol_Lanczos);
556}
557
558template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
559void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
560t_step_pivot (double x, const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, double tol_Lanczos)
561{
562 last_dt = dt;
563 turnaround(pivot, N_sites, CURRENT_DIRECTION);
564
565 // 2-site propagation
566 size_t loc1 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot : pivot-1;
567 size_t loc2 = (CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT)? pivot+1 : pivot;
568// cout << "2site between: " << loc1 << "," << loc2 << ", pivot=" << pivot << endl;
569
570 Stopwatch<> Wc;
571 PivotVector<Symmetry,TimeScalar> Apair(Vinout.A[loc1], Vinout.locBasis(loc1),
572 Vinout.A[loc2], Vinout.locBasis(loc2),
573 Vinout.QoutTop[loc1], Vinout.QoutBot[loc1]);
574 t_contr += Wc.time(SECONDS);
575 PivotMatrix2<Symmetry,TimeScalar,MpoScalar> Heff2(Heff[loc1].Terms[0].L, Heff[loc2].Terms[0].R,
576 H.W_at(loc1), H.W_at(loc2),
577 H.locBasis(loc1), H.locBasis(loc2),
578 H.opBasis(loc1), H.opBasis(loc2));
579
580 Stopwatch<> Woh2;
581 // WRONG FOR MULTISITE TERMS:
582// precalc_blockStructure (Heff[loc1].L, Apair.data, Heff2.W12, Heff2.W34, Apair.data, Heff[loc2].R,
583// H.locBasis(loc1), H.locBasis(loc2), H.opBasis(loc1), H.opBasis(loc2),
584// H.TSD[loc1],
585// Heff2.qlhs, Heff2.qrhs, Heff2.factor_cgcs);
586 precalc_blockStructure (Heff2.Terms[0].L, Apair.data, Heff2.Terms[0].W12, Heff2.Terms[0].W34, Apair.data, Heff2.Terms[0].R,
587 H.locBasis(loc1), H.locBasis(loc2), H.opBasis(loc1), H.opBasis(loc2),
588 Heff2.qlhs, Heff2.qrhs, Heff2.factor_cgcs);
589 t_ohead += Woh2.time(SECONDS);
590
591 LanczosPropagator<PivotMatrix2<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz2(tol_Lanczos);
592 Stopwatch<> W2;
593// Lutz2.t_step(Heff2, Apair, -x*dt.imag()); // 2-site algorithm
594 Lutz2.t_step(Heff2, Apair, x*dt); // 2-site algorithm
595// cout << Lutz2.info() << endl;
596 t_2site += W2.time(SECONDS);
597
598 dimK2_log.push_back(Lutz2.get_dimK());
599 if (Lutz2.get_dist() > dist_max) {dist_max = Lutz2.get_dist();}
600 if (Lutz2.get_dimK() > dimK_max) {dimK_max = Lutz2.get_dimK();}
601
602 Stopwatch<> Ws;
603 Vinout.sweepStep2(CURRENT_DIRECTION, loc1, Apair.data);
604 t_contr += Ws.time(SECONDS);
605
606 (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)? build_L(H,Vinout,loc2) : build_R(H,Vinout,loc1);
607 pivot = Vinout.get_pivot();
608
609 // 1-site propagation
610 if ((CURRENT_DIRECTION==DMRG::DIRECTION::RIGHT and pivot != N_sites-1) or
611 (CURRENT_DIRECTION==DMRG::DIRECTION::LEFT and pivot != 0))
612 {
613 Stopwatch<> Woh1;
614 precalc_blockStructure (Heff[pivot].Terms[0].L, Vinout.A[pivot], Heff[pivot].Terms[0].W, Vinout.A[pivot], Heff[pivot].Terms[0].R,
615 H.locBasis(pivot), H.opBasis(pivot),
616 Heff[pivot].qlhs, Heff[pivot].qrhs, Heff[pivot].factor_cgcs);
617 t_ohead += Woh1.time(SECONDS);
618
619 PivotVector<Symmetry,TimeScalar> Asingle(Vinout.A[pivot]);
620// test_edge_eigenvector(Asingle);
621
622// cout << "1site at: " << pivot << endl;
623 LanczosPropagator<PivotMatrix1<Symmetry,TimeScalar,MpoScalar>, PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz(tol_Lanczos);
624 Stopwatch<> W1;
625// Lutz.t_step(Heff[pivot], Asingle, +x*dt.imag()); // 1-site algorithm
626 Lutz.t_step(Heff[pivot], Asingle, -x*dt); // 1-site algorithm
627 deltaE(pivot) = Lutz.get_deltaE();
628// cout << Lutz.info() << endl;
629 t_1site += W1.time(SECONDS);
630
631 dimK1_log.push_back(Lutz.get_dimK());
632 if (Lutz.get_dist() > dist_max) {dist_max = Lutz2.get_dist();}
633 if (Lutz.get_dimK() > dimK_max) {dimK_max = Lutz2.get_dimK();}
634
635 Vinout.A[pivot] = Asingle.data;
636 }
637}
638
639template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
640void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
641t0_step_pivot (bool BACK, double x, const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, double tol_Lanczos, bool TURN_FIRST)
642{
643 last_dt = dt;
644 if (TURN_FIRST) turnaround(pivot, N_sites, CURRENT_DIRECTION);
645
646 // 1-site propagation
647 PivotVector<Symmetry,TimeScalar> Asingle(Vinout.A[pivot]);
648 Stopwatch<> Woh1;
649 precalc_blockStructure (Heff[pivot].Terms[0].L, Vinout.A[pivot], Heff[pivot].Terms[0].W, Vinout.A[pivot], Heff[pivot].Terms[0].R,
650 H.locBasis(pivot), H.opBasis(pivot),
651 Heff[pivot].qlhs, Heff[pivot].qrhs, Heff[pivot].factor_cgcs);
652 t_ohead += Woh1.time(SECONDS);
653
654// test_edge_eigenvector(Asingle);
655
656// cout << "1site at: " << pivot << endl;
657 LanczosPropagator<PivotMatrix1<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz(tol_Lanczos);
658
659 Stopwatch<> W1;
660// Lutz.t_step(Heff[pivot], Asingle, -x*dt.imag()); // 1-site algorithm
661 Lutz.t_step(Heff[pivot], Asingle, x*dt); // 1-site algorithm
662 deltaE(pivot) = Lutz.get_deltaE();
663 t_1site += W1.time(SECONDS);
664
665 dimK1_log.push_back(Lutz.get_dimK());
666 if (Lutz.get_dist() > dist_max) {dist_max = Lutz.get_dist();}
667 if (Lutz.get_dimK() > dimK_max) {dimK_max = Lutz.get_dimK();}
668 Vinout.A[pivot] = Asingle.data;
669
670 // 0-site propagation
671 if (BACK)
672 {
673 PivotVector<Symmetry,TimeScalar> Azero;
674 int old_pivot = pivot;
675 if (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)
676 {
677// cout << "split to RIGHT at pivot=" << pivot << endl;
678 Vinout.rightSplitStep(pivot,Azero.data[0]);
679 }
680 else
681 {
682// cout << "split to LEFT at pivot=" << pivot << endl;
683 Vinout.leftSplitStep(pivot,Azero.data[0]);
684 }
685 pivot = Vinout.get_pivot();
686 if (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)
687 {
688 if (old_pivot+1 == N_sites)
689 {
690 build_L(H,Vinout,N_sites);
691 }
692 else
693 {
694 build_L(H,Vinout,pivot);
695 }
696 }
697 else
698 {
699 if (old_pivot-1 == -1)
700 {
701 build_R(H,Vinout,-1);
702 }
703 else
704 {
705 build_R(H,Vinout,pivot);
706 }
707 }
708
709// if (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)
710// {
711// cout << "0site between " << old_pivot << "," << old_pivot+1 << endl;
712// }
713// else
714// {
715// cout << "0site between " << old_pivot-1 << "," << old_pivot << endl;
716// }
717 LanczosPropagator<PivotMatrix0<Symmetry,TimeScalar,MpoScalar>,PivotVector<Symmetry,TimeScalar>,TimeScalar> Lutz0(tol_Lanczos);
718
719 PivotMatrix0<Symmetry,TimeScalar,MpoScalar> Heff0;
720 if (CURRENT_DIRECTION == DMRG::DIRECTION::RIGHT)
721 {
722 if (old_pivot+1 == N_sites)
723 {
724 Heff0 = PivotMatrix0<Symmetry,TimeScalar,MpoScalar>(HeffLast.Terms[0].L, Heff[old_pivot].Terms[0].R);
725 }
726 else
727 {
728 Heff0 = PivotMatrix0<Symmetry,TimeScalar,MpoScalar>(Heff[old_pivot+1].Terms[0].L, Heff[old_pivot].Terms[0].R);
729 }
730 }
731 else
732 {
733 if (old_pivot-1 == -1)
734 {
735 Heff0 = PivotMatrix0<Symmetry,TimeScalar,MpoScalar>(Heff[old_pivot].Terms[0].L, HeffFrst.Terms[0].R);
736 }
737 else
738 {
739 Heff0 = PivotMatrix0<Symmetry,TimeScalar,MpoScalar>(Heff[old_pivot].Terms[0].L, Heff[old_pivot-1].Terms[0].R);
740 }
741 }
742
743 Stopwatch<> W0;
744// Lutz0.t_step(Heff0, Azero, +x*dt.imag()); // 0-site algorithm
745 Lutz0.t_step(Heff0, Azero, -x*dt); // 0-site algorithm
746// cout << Lutz0.info() << endl;
747 t_0site += W0.time(SECONDS);
748
749 dimK0_log.push_back(Lutz0.get_dimK());
750 if (Lutz0.get_dist() > dist_max) {dist_max = Lutz0.get_dist();}
751 if (Lutz0.get_dimK() > dimK_max) {dimK_max = Lutz0.get_dimK();}
752
753 if (!TURN_FIRST) turnaround(pivot, N_sites, CURRENT_DIRECTION);
754
755// cout << "absorb at pivot=" << pivot << " going " << CURRENT_DIRECTION << endl;
756 Vinout.absorb(pivot, CURRENT_DIRECTION, Azero.data[0]);
757 }
758}
759
760template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
761void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
762t_step_adaptive (const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, const vector<bool> &TWO_STEP_AT, int N_stages, double tol_Lanczos)
763{
764 last_dt = dt;
765 assert(N_stages==1 and "Only N_stages=1 implemented for TDVPPropagator::t_step_adaptive!");
766 dist_max = 0.;
767 dimK_max = 0;
768 N_stages_last = N_stages;
769
770 t_0site = 0;
771 t_1site = 0;
772 t_2site = 0;
773 t_ohead = 0;
774 t_contr = 0;
775 t_tot = 0;
776
777 Stopwatch<> Wtot;
778
779 dimK2_log.clear();
780 dimK1_log.clear();
781 dimK0_log.clear();
782
783 for (size_t l=0; l<N_sites-1; ++l)
784 {
785 if (TWO_STEP_AT[l] == true)
786 {
787 t_step_pivot(x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
788 }
789 else
790 {
791 t0_step_pivot(true,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
792 }
793 }
794
795 if (TWO_STEP_AT[N_sites-2])
796 {
797 t_step_pivot(x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
798 }
799 else
800 {
801 // at N_sites-1 RIGHT
802 // at N_sites-1 LEFT
803// if (Vinout.Boundaries.IS_TRIVIAL())
804// {
805 t0_step_pivot(false,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
806 t0_step_pivot(true,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
807// }
808// else
809// {
810// t0_step_pivot(true,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos,false);
812// t0_step_pivot(true,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
813// }
814 }
815
816 for (int l=N_sites-3; l>=0; --l)
817 {
818 if (TWO_STEP_AT[l] == true)
819 {
820 t_step_pivot(x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
821 }
822 else
823 {
824 t0_step_pivot(true,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
825 }
826 }
827
828 if (!TWO_STEP_AT[0])
829 {
830// if (Vinout.Boundaries.IS_TRIVIAL())
831// {
832 t0_step_pivot(false,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos);
833// }
834// else
835// {
836// t0_step_pivot(true,x(1,0,N_stages),H,Vinout,dt,tol_Lanczos,false);
837// }
838 }
839
840// cout << "final pivot=" << pivot << endl << endl;
841
842 t_tot = Wtot.time(SECONDS);
843}
844
845template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
846void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
847test_edge_eigenvector (const PivotVector<Symmetry,TimeScalar> &Asingle)
848{
849 if (pivot == 0 or pivot == N_sites-1)
850 {
851 auto V = Asingle;
852 normalize(V);
853 PivotVector<Symmetry,TimeScalar> HV;
854 HxV(Heff[pivot],V,HV);
855 double res = abs(dot(HV,HV).real()-pow(dot(V,HV).real(),2));
856// if (pivot==0)
857// {
858// EvarL = res;
859// }
860// else
861// {
862// EvarR = res;
863// }
864// lout << "pivot=" << pivot << ", eigenstate test=" << res << endl;
865 }
866}
867
868template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
869inline void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
870build_L (const Hamiltonian &H, const VectorType &Vinout, int loc)
871{
872 if (loc != 0)
873 {
874 if (loc == N_sites)
875 {
876 contract_L(Heff[loc-1].Terms[0].L, Vinout.A[loc-1], H.W[loc-1], Vinout.A[loc-1], H.locBasis(loc-1), H.opBasis(loc-1), HeffLast.Terms[0].L);
877 }
878 else
879 {
880 contract_L(Heff[loc-1].Terms[0].L, Vinout.A[loc-1], H.W[loc-1], Vinout.A[loc-1], H.locBasis(loc-1), H.opBasis(loc-1), Heff[loc].Terms[0].L);
881 }
882 }
883}
884
885template<typename Hamiltonian, typename Symmetry, typename MpoScalar, typename TimeScalar, typename VectorType>
886inline void TDVPPropagator<Hamiltonian,Symmetry,MpoScalar,TimeScalar,VectorType>::
887build_R (const Hamiltonian &H, const VectorType &Vinout, int loc)
888{
889 if (loc != N_sites-1)
890 {
891 if (loc == -1)
892 {
893 contract_R(Heff[loc+1].Terms[0].R, Vinout.A[loc+1], H.W[loc+1], Vinout.A[loc+1], H.locBasis(loc+1), H.opBasis(loc+1), HeffFrst.Terms[0].R);
894 }
895 else
896 {
897 contract_R(Heff[loc+1].Terms[0].R, Vinout.A[loc+1], H.W[loc+1], Vinout.A[loc+1], H.locBasis(loc+1), H.opBasis(loc+1), Heff[loc].Terms[0].R);
898 }
899 }
900}
901
902#endif
precalc_blockStructure
void precalc_blockStructure(const Tripod< Symmetry, Eigen::Matrix< Scalar, Dynamic, Dynamic > > &L, const vector< Biped< Symmetry, Eigen::Matrix< Scalar, Dynamic, Dynamic > > > &Abra, const vector< vector< vector< Biped< Symmetry, MpoMatrixType > > > > &W, const vector< Biped< Symmetry, Eigen::Matrix< Scalar, Dynamic, Dynamic > > > &Aket, const Tripod< Symmetry, Eigen::Matrix< Scalar, Dynamic, Dynamic > > &R, const vector< qarray< Symmetry::Nq > > &qloc, const vector< qarray< Symmetry::Nq > > &qOp, vector< std::array< size_t, 2 > > &qlhs, vector< vector< std::array< size_t, 6 > > > &qrhs, vector< vector< Scalar > > &factor_cgcs)
Definition: DmrgIndexGymnastics.h:488
turnaround
void turnaround(int pivot, size_t L, DMRG::DIRECTION::OPTION &DIR)
Flips the sweep direction when the edge is reached.
Definition: DmrgJanitor.h:7
dot
Scalar dot(const Mps< Symmetry, Scalar > &Vbra, const Mps< Symmetry, Scalar > &Vket)
Definition: DmrgLinearAlgebra.h:24
HxV
void HxV(const Mpo< Symmetry, MpoScalar > &H, const Mps< Symmetry, Scalar > &Vin, Mps< Symmetry, Scalar > &Vout, bool VERBOSE=true, double tol_compr=1e-4, int Mincr=100, int Mlimit=10000, int max_halfsweeps=100)
Definition: DmrgLinearAlgebra.h:500
normalize
void normalize(PivotVector< Symmetry, Scalar_ > &V)
Definition: DmrgPivotVector.h:355
isReal
double isReal(double x)
Definition: DmrgTypedefs.h:21
TDVPPropagator::get_deltaE
VectorXd get_deltaE() const
Definition: TDVPPropagator.h:32
TDVPPropagator::CURRENT_DIRECTION
DMRG::DIRECTION::OPTION CURRENT_DIRECTION
Definition: TDVPPropagator.h:57
TDVPPropagator::t_step_pivot
void t_step_pivot(double x, const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, double tol_Lanczos=1e-8)
Definition: TDVPPropagator.h:560
TDVPPropagator::last_dt
TimeScalar last_dt
Definition: TDVPPropagator.h:76
TDVPPropagator::t_step0
void t_step0(const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, int N_stages=1, double tol_Lanczos=1e-8)
Definition: TDVPPropagator.h:451
TDVPPropagator::dimKlog
VectorXd dimKlog
Definition: TDVPPropagator.h:45
TDVPPropagator::get_dimK0_log
vector< size_t > get_dimK0_log() const
Definition: TDVPPropagator.h:36
TDVPPropagator::t_step
void t_step(const Hamiltonian &H, const VectorType &Vin, VectorType &Vout, TimeScalar dt, int N_stages=1, double tol_Lanczos=1e-8)
Definition: TDVPPropagator.h:551
TDVPPropagator::build_R
void build_R(const Hamiltonian &H, const VectorType &Vinout, int loc)
Definition: TDVPPropagator.h:887
TDVPPropagator::memory
double memory(MEMUNIT memunit=GB) const
Definition: TDVPPropagator.h:116
TDVPPropagator::get_t_tot
double get_t_tot() const
Definition: TDVPPropagator.h:33
TDVPPropagator::x
double x(int alg, size_t l, int N_stages)
Definition: TDVPPropagator.h:212
TDVPPropagator::info
string info() const
Definition: TDVPPropagator.h:81
TDVPPropagator::build_L
void build_L(const Hamiltonian &H, const VectorType &Vinout, int loc)
Definition: TDVPPropagator.h:870
TDVPPropagator::overhead
double overhead(MEMUNIT memunit=GB) const
Definition: TDVPPropagator.h:134
TDVPPropagator::HeffFrst
PivotMatrix1< Symmetry, TimeScalar, MpoScalar > HeffFrst
Definition: TDVPPropagator.h:49
TDVPPropagator::test_edge_eigenvector
void test_edge_eigenvector(const PivotVector< Symmetry, TimeScalar > &Asingle)
Definition: TDVPPropagator.h:847
TDVPPropagator::get_dimK2_log
vector< size_t > get_dimK2_log() const
Definition: TDVPPropagator.h:34
TDVPPropagator::Heff
vector< PivotMatrix1< Symmetry, TimeScalar, MpoScalar > > Heff
Definition: TDVPPropagator.h:47
TDVPPropagator::t_step_adaptive
void t_step_adaptive(const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, const vector< bool > &TWO_STEP_AT, int N_stages=1, double tol_Lanczos=1e-8)
Definition: TDVPPropagator.h:762
TDVPPropagator::set_blocks
void set_blocks(const Hamiltonian &H, VectorType &Vinout)
Definition: TDVPPropagator.h:157
TDVPPropagator::dimK2_log
vector< size_t > dimK2_log
Definition: TDVPPropagator.h:66
TDVPPropagator::t0_step_pivot
void t0_step_pivot(bool BACK, double x, const Hamiltonian &H, VectorType &Vinout, TimeScalar dt, double tol_Lanczos=1e-8, bool TURN_FIRST=true)
Definition: TDVPPropagator.h:641
TDVPPropagator::HeffLast
PivotMatrix1< Symmetry, TimeScalar, MpoScalar > HeffLast
Definition: TDVPPropagator.h:48
TDVPPropagator::dimK0_log
vector< size_t > dimK0_log
Definition: TDVPPropagator.h:64
TDVPPropagator::get_dimK1_log
vector< size_t > get_dimK1_log() const
Definition: TDVPPropagator.h:35
TDVPPropagator::dimK1_log
vector< size_t > dimK1_log
Definition: TDVPPropagator.h:65
contract_R
void contract_R(const Tripod< Symmetry, MatrixType2 > &Rold, const vector< Biped< Symmetry, MatrixType > > &Abra, const vector< vector< vector< Biped< Symmetry, MpoMatrixType > > > > &W, const vector< Biped< Symmetry, MatrixType > > &Aket, const vector< qarray< Symmetry::Nq > > &qloc, const vector< qarray< Symmetry::Nq > > &qOp, Tripod< Symmetry, MatrixType2 > &Rnew, bool RANDOMIZE=false, tuple< CONTRACT_LR_MODE, size_t > MODE_input=make_pair(FULL, 0), const std::unordered_map< pair< qarray< Symmetry::Nq >, size_t >, size_t > &basis_order_map={})
Definition: DmrgContractions.h:180
contract_L
void contract_L(const Tripod< Symmetry, MatrixType2 > &Lold, const vector< Biped< Symmetry, MatrixType > > &Abra, const vector< vector< vector< Biped< Symmetry, MpoMatrixType > > > > &W, const vector< Biped< Symmetry, MatrixType > > &Aket, const vector< qarray< Symmetry::Nq > > &qloc, const vector< qarray< Symmetry::Nq > > &qOp, Tripod< Symmetry, MatrixType2 > &Lnew, bool RANDOMIZE=false, tuple< CONTRACT_LR_MODE, size_t > MODE_input=make_pair(FULL, 0), const std::unordered_map< pair< qarray< Symmetry::Nq >, size_t >, size_t > &basis_order_map={})
Definition: DmrgContractions.h:35
PivotMatrix2::qrhs
vector< vector< std::array< size_t, 12 > > > qrhs
Definition: DmrgPivotMatrix2.h:68
PivotMatrix2::qlhs
vector< std::array< size_t, 2 > > qlhs
Definition: DmrgPivotMatrix2.h:67
PivotMatrix2::factor_cgcs
vector< vector< Scalar > > factor_cgcs
Definition: DmrgPivotMatrix2.h:69
PivotMatrix2::Terms
vector< PivotMatrix2Terms< Symmetry, Scalar, MpoScalar > > Terms
Definition: DmrgPivotMatrix2.h:64
PivotVector::data
vector< Biped< Symmetry, Matrix< Scalar_, Dynamic, Dynamic > > > data
Definition: DmrgPivotVector.h:116