ALPS MPS Codes
Reference documentation.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ts_optimize< Matrix, SymmGroup, Storage > Class Template Reference

#include <ts_optimize.hpp>

Inheritance diagram for ts_optimize< Matrix, SymmGroup, Storage >:
optimizer_base< Matrix, SymmGroup, Storage >

Public Types

typedef optimizer_base< Matrix,
SymmGroup, Storage > 
base
 

Public Member Functions

 ts_optimize (MPS< Matrix, SymmGroup > &mps_, MPO< Matrix, SymmGroup > const &mpo_, BaseParameters &parms_, boost::function< bool()> stop_callback_, int initial_site_=0)
 
int to_site (const int L, const int i) const
 
void sweep (int sweep, OptimizeDirection d=Both)
 
results_collector const & iteration_results () const
 

Protected Member Functions

void boundary_left_step (MPO< Matrix, SymmGroup > const &mpo, int site)
 
void boundary_right_step (MPO< Matrix, SymmGroup > const &mpo, int site)
 
void init_left_right (MPO< Matrix, SymmGroup > const &mpo, int site)
 
double get_cutoff (int sweep) const
 
std::size_t get_Mmax (int sweep) const
 

Protected Attributes

results_collector iteration_results_
 
MPS< Matrix, SymmGroup > & mps
 
MPO< Matrix, SymmGroup > const & mpo
 
BaseParametersparms
 
boost::function< bool()> stop_callback
 
std::vector< Boundary
< typename
storage::constrained< Matrix >
::type, SymmGroup > > 
left_
 
std::vector< Boundary
< typename
storage::constrained< Matrix >
::type, SymmGroup > > 
right_
 
unsigned int northo
 
std::vector< std::vector
< block_matrix< typename
storage::constrained< Matrix >
::type, SymmGroup > > > 
ortho_left_
 
std::vector< std::vector
< block_matrix< typename
storage::constrained< Matrix >
::type, SymmGroup > > > 
ortho_right_
 
std::vector< MPS< Matrix,
SymmGroup > > 
ortho_mps
 

Detailed Description

template<class Matrix, class SymmGroup, class Storage>
class ts_optimize< Matrix, SymmGroup, Storage >

Definition at line 40 of file ts_optimize.hpp.

Member Typedef Documentation

template<class Matrix , class SymmGroup , class Storage >
typedef optimizer_base<Matrix, SymmGroup, Storage> ts_optimize< Matrix, SymmGroup, Storage >::base

Definition at line 44 of file ts_optimize.hpp.

Constructor & Destructor Documentation

template<class Matrix , class SymmGroup , class Storage >
ts_optimize< Matrix, SymmGroup, Storage >::ts_optimize ( MPS< Matrix, SymmGroup > &  mps_,
MPO< Matrix, SymmGroup > const &  mpo_,
BaseParameters parms_,
boost::function< bool()>  stop_callback_,
int  initial_site_ = 0 
)
inline

Definition at line 53 of file ts_optimize.hpp.

58  : base(mps_, mpo_, parms_, stop_callback_, to_site(mps_.length(), initial_site_))
59  , initial_site((initial_site_ < 0) ? 0 : initial_site_)
60  {
61  locale_shared l; // cache twosite mpo
62  make_ts_cache_mpo(mpo, ts_cache_mpo, mps);
63  }
int to_site(const int L, const int i) const
Definition: ts_optimize.hpp:65
std::size_t locale_shared
void make_ts_cache_mpo(MPO< MPOMatrix, SymmGroup > const &mpo_orig, MPO< MPSMatrix, SymmGroup > &mpo_out, MPS< MPSMatrix, SymmGroup > const &mps)
Definition: ts_ops.h:198
MPO< Matrix, SymmGroup > const & mpo
Definition: optimize.h:230
MPS< Matrix, SymmGroup > & mps
Definition: optimize.h:229
size_t length() const
Definition: mps.h:58
optimizer_base< Matrix, SymmGroup, Storage > base
Definition: ts_optimize.hpp:44

Member Function Documentation

template<class Matrix , class SymmGroup , class Storage >
void optimizer_base< Matrix, SymmGroup, Storage >::boundary_left_step ( MPO< Matrix, SymmGroup > const &  mpo,
int  site 
)
inlineprotectedinherited

Definition at line 135 of file optimize.h.

136  {
137  left_[site+1] = contraction::overlap_mpo_left_step(mps[site], mps[site], left_[site], mpo[site]);
138  Storage::pin(left_[site+1]);
139 
140  for (int n = 0; n < northo; ++n)
141  ortho_left_[n][site+1] = contraction::overlap_left_step(mps[site], ortho_mps[n][site], ortho_left_[n][site]);
142  }
std::vector< Boundary< typename storage::constrained< Matrix >::type, SymmGroup > > left_
Definition: optimize.h:235
std::vector< std::vector< block_matrix< typename storage::constrained< Matrix >::type, SymmGroup > > > ortho_left_
Definition: optimize.h:239
static block_matrix< OtherMatrix, SymmGroup > overlap_left_step(MPSTensor< Matrix, SymmGroup > const &bra_tensor, MPSTensor< Matrix, SymmGroup > const &ket_tensor, block_matrix< OtherMatrix, SymmGroup > const &left, block_matrix< OtherMatrix, SymmGroup > *localop=NULL)
Definition: contractions.h:41
MPS< Matrix, SymmGroup > & mps
Definition: optimize.h:229
static Boundary< OtherMatrix, SymmGroup > overlap_mpo_left_step(MPSTensor< Matrix, SymmGroup > const &bra_tensor, MPSTensor< Matrix, SymmGroup > const &ket_tensor, Boundary< OtherMatrix, SymmGroup > const &left, MPOTensor< Matrix, SymmGroup > const &mpo)
Definition: contractions.h:340
unsigned int northo
Definition: optimize.h:238
std::vector< MPS< Matrix, SymmGroup > > ortho_mps
Definition: optimize.h:240
template<class Matrix , class SymmGroup , class Storage >
void optimizer_base< Matrix, SymmGroup, Storage >::boundary_right_step ( MPO< Matrix, SymmGroup > const &  mpo,
int  site 
)
inlineprotectedinherited

Definition at line 144 of file optimize.h.

145  {
146  right_[site] = contraction::overlap_mpo_right_step(mps[site], mps[site], right_[site+1], mpo[site]);
147  Storage::pin(right_[site]);
148 
149  for (int n = 0; n < northo; ++n)
150  ortho_right_[n][site] = contraction::overlap_right_step(mps[site], ortho_mps[n][site], ortho_right_[n][site+1]);
151  }
static block_matrix< OtherMatrix, SymmGroup > overlap_right_step(MPSTensor< Matrix, SymmGroup > const &bra_tensor, MPSTensor< Matrix, SymmGroup > const &ket_tensor, block_matrix< OtherMatrix, SymmGroup > const &right, block_matrix< OtherMatrix, SymmGroup > *localop=NULL)
Definition: contractions.h:71
static Boundary< OtherMatrix, SymmGroup > overlap_mpo_right_step(MPSTensor< Matrix, SymmGroup > const &bra_tensor, MPSTensor< Matrix, SymmGroup > const &ket_tensor, Boundary< OtherMatrix, SymmGroup > const &right, MPOTensor< Matrix, SymmGroup > const &mpo)
Definition: contractions.h:386
std::vector< std::vector< block_matrix< typename storage::constrained< Matrix >::type, SymmGroup > > > ortho_right_
Definition: optimize.h:239
MPS< Matrix, SymmGroup > & mps
Definition: optimize.h:229
unsigned int northo
Definition: optimize.h:238
std::vector< Boundary< typename storage::constrained< Matrix >::type, SymmGroup > > right_
Definition: optimize.h:235
std::vector< MPS< Matrix, SymmGroup > > ortho_mps
Definition: optimize.h:240
template<class Matrix , class SymmGroup , class Storage >
double optimizer_base< Matrix, SymmGroup, Storage >::get_cutoff ( int  sweep) const
inlineprotectedinherited

Definition at line 202 of file optimize.h.

203  {
204  double cutoff;
205  if (sweep >= parms.template get<int>("ngrowsweeps"))
206  cutoff = parms.template get<double>("truncation_final");
207  else
208  cutoff = log_interpolate(parms.template get<double>("truncation_initial"), parms.template get<double>("truncation_final"), parms.template get<int>("ngrowsweeps"), sweep);
209  return cutoff;
210  }
BaseParameters & parms
Definition: optimize.h:232
virtual void sweep(int sweep, OptimizeDirection d=Both)=0
double log_interpolate(double y0, double y1, int N, int i)
Definition: optimize.h:74
template<class Matrix , class SymmGroup , class Storage >
std::size_t optimizer_base< Matrix, SymmGroup, Storage >::get_Mmax ( int  sweep) const
inlineprotectedinherited

Definition at line 212 of file optimize.h.

213  {
214  std::size_t Mmax;
215  if (parms.is_set("sweep_bond_dimensions")) {
216  std::vector<std::size_t> ssizes = parms.template get<std::vector<std::size_t> >("sweep_bond_dimensions");
217  if (sweep >= ssizes.size())
218  Mmax = *ssizes.rbegin();
219  else
220  Mmax = ssizes[sweep];
221  } else
222  Mmax = parms.template get<std::size_t>("max_bond_dimension");
223  return Mmax;
224  }
BaseParameters & parms
Definition: optimize.h:232
bool is_set(std::string const &key) const
virtual void sweep(int sweep, OptimizeDirection d=Both)=0
template<class Matrix , class SymmGroup , class Storage >
void optimizer_base< Matrix, SymmGroup, Storage >::init_left_right ( MPO< Matrix, SymmGroup > const &  mpo,
int  site 
)
inlineprotectedinherited

Definition at line 153 of file optimize.h.

154  {
156  std::size_t L = mps.length();
157 
158  left_.resize(mpo.length()+1);
159  right_.resize(mpo.length()+1);
160 
161  ortho_left_.resize(northo);
162  ortho_right_.resize(northo);
163  for (int n = 0; n < northo; ++n) {
164  ortho_left_[n].resize(L+1);
165  ortho_right_[n].resize(L+1);
166 
167  ortho_left_[n][0] = mps.left_boundary()[0];
168  ortho_right_[n][L] = mps.right_boundary()[0];
169  }
170 
171  //Timer tlb("Init left boundaries"); tlb.begin();
172  Storage::drop(left_[0]);
173  left_[0] = mps.left_boundary();
174  Storage::pin(left_[0]);
175 
176  for (int i = 0; i < site; ++i) {
177  Storage::drop(left_[i+1]);
178  boundary_left_step(mpo, i);
179  Storage::evict(left_[i]);
180  }
181  Storage::evict(left_[site]);
182  //tlb.end();
183 
184  maquis::cout << "Boundaries are partially initialized...\n";
185 
186  //Timer trb("Init right boundaries"); trb.begin();
187  Storage::drop(right_[L]);
188  right_[L] = mps.right_boundary();
189  Storage::pin(right_[L]);
190 
191  for (int i = L-1; i >= site; --i) {
192  Storage::drop(right_[i]);
193  boundary_right_step(mpo, i);
194  Storage::evict(right_[i+1]);
195  }
196  Storage::evict(right_[site]);
197  //trb.end();
198 
199  maquis::cout << "Boundaries are fully initialized...\n";
200  }
std::vector< std::vector< int > > construct_placements(const MPO< Matrix, SymmGroup > &mpo)
Definition: placement.h:101
std::vector< Boundary< typename storage::constrained< Matrix >::type, SymmGroup > > left_
Definition: optimize.h:235
std::vector< std::vector< block_matrix< typename storage::constrained< Matrix >::type, SymmGroup > > > ortho_left_
Definition: optimize.h:239
std::vector< std::vector< block_matrix< typename storage::constrained< Matrix >::type, SymmGroup > > > ortho_right_
Definition: optimize.h:239
MPS< Matrix, SymmGroup > & mps
Definition: optimize.h:229
unsigned int northo
Definition: optimize.h:238
std::vector< Boundary< typename storage::constrained< Matrix >::type, SymmGroup > > right_
Definition: optimize.h:235
void boundary_right_step(MPO< Matrix, SymmGroup > const &mpo, int site)
Definition: optimize.h:144
void boundary_left_step(MPO< Matrix, SymmGroup > const &mpo, int site)
Definition: optimize.h:135
template<class Matrix , class SymmGroup , class Storage >
results_collector const& optimizer_base< Matrix, SymmGroup, Storage >::iteration_results ( ) const
inlineinherited

Definition at line 131 of file optimize.h.

131 { return iteration_results_; }
results_collector iteration_results_
Definition: optimize.h:227
template<class Matrix , class SymmGroup , class Storage >
void ts_optimize< Matrix, SymmGroup, Storage >::sweep ( int  sweep,
OptimizeDirection  d = Both 
)
inlinevirtual

Implements optimizer_base< Matrix, SymmGroup, Storage >.

Definition at line 71 of file ts_optimize.hpp.

72  {
73  boost::chrono::high_resolution_clock::time_point sweep_now = boost::chrono::high_resolution_clock::now();
74 
76 
77  std::size_t L = mps.length();
78 
79  int _site = 0, site = 0;
80  if (initial_site != -1) {
81  _site = initial_site;
82  site = to_site(L, _site);
83  }
84 
85  if (_site < L-1) {
86  Storage::prefetch(left_[site]);
87  Storage::prefetch(right_[site+2]);
88  } else {
89  Storage::prefetch(left_[site-1]);
90  Storage::prefetch(right_[site+1]);
91  }
92 
93 #ifndef NDEBUG
94  maquis::cout << mps.description() << std::endl;
95 #endif
96  for (; _site < 2*L-2; ++_site) {
97  /* (0,1), (1,2), ... , (L-1,L), (L-1,L), (L-2, L-1), ... , (0,1)
98  | | |
99  site 1 |
100  | left to right | right to left, lr = -1
101  site 2 | */
102 
103  int lr, site1, site2;
104  if (_site < L-1) {
105  site = to_site(L, _site);
106  lr = 1;
107  site1 = site;
108  site2 = site+1;
109  ts_cache_mpo[site1].placement_l = mpo[site1].placement_l;
110  ts_cache_mpo[site1].placement_r = get_right_placement(ts_cache_mpo[site1], mpo[site1].placement_l, mpo[site2].placement_r);
111  } else {
112  site = to_site(L, _site);
113  lr = -1;
114  site1 = site-1;
115  site2 = site;
116  ts_cache_mpo[site1].placement_l = get_left_placement(ts_cache_mpo[site1], mpo[site1].placement_l, mpo[site2].placement_r);
117  ts_cache_mpo[site1].placement_r = mpo[site2].placement_r;
118  }
119 
120  //if (lr == +1) mps.canonize(site1);
121  //else mps.canonize(site2);
122 
123  maquis::cout << std::endl;
124  maquis::cout << "Sweep " << sweep << ", optimizing sites " << site1 << " and " << site2 << std::endl;
125 
126  // MD: some changes needed to re-enable it.
127 // if (parms.template get<bool>("beta_mode")) {
128 // if (sweep == 0 && lr == 1) {
129 // mpo = zero_after(mpo_orig, 0);
130 // if (site == 0)
131 // this->init_left_right(mpo, 0);
132 // } else if (sweep == 0 && lr == -1 && site == L-1) {
133 // mpo = mpo_orig;
134 // }
135 // }
136 
137  if (_site != L-1)
138  {
139  Storage::fetch(left_[site1]);
140  Storage::fetch(right_[site2+1]);
141  }
142 
143  if (lr == +1) {
144  if (site2+2 < right_.size()){
145  Storage::prefetch(right_[site2+2]);
146  }
147  } else {
148  if (site1 > 0){
149  Storage::prefetch(left_[site1-1]);
150  }
151  }
152 
153 
154  boost::chrono::high_resolution_clock::time_point now, then;
155 
156  // Create TwoSite objects
157  TwoSiteTensor<Matrix, SymmGroup> tst(mps[site1], mps[site2]);
158  MPSTensor<Matrix, SymmGroup> twin_mps = tst.make_mps();
159  SiteProblem<Matrix, SymmGroup> sp(left_[site1], right_[site2+1], ts_cache_mpo[site1]);
160 
161  /// Compute orthogonal vectors
162  std::vector<MPSTensor<Matrix, SymmGroup> > ortho_vecs(base::northo);
163  for (int n = 0; n < base::northo; ++n) {
165  ortho_vecs[n] = contraction::site_ortho_boundaries(twin_mps, ts_ortho.make_mps(),
166  base::ortho_left_[n][site1], base::ortho_right_[n][site2+1]);
167  }
168 
169  std::pair<double, MPSTensor<Matrix, SymmGroup> > res;
170 
171  if (d == Both ||
172  (d == LeftOnly && lr == -1) ||
173  (d == RightOnly && lr == +1))
174  {
175  if (parms["eigensolver"] == std::string("IETL")) {
176  BEGIN_TIMING("IETL")
177  res = solve_ietl_lanczos(sp, twin_mps, parms);
178  END_TIMING("IETL")
179  } else if (parms["eigensolver"] == std::string("IETL_JCD")) {
180  BEGIN_TIMING("JCD")
181  res = solve_ietl_jcd(sp, twin_mps, parms, ortho_vecs);
182  END_TIMING("JCD")
183  } else {
184  throw std::runtime_error("I don't know this eigensolver.");
185  }
186 
187  tst << res.second;
188  }
189 
190 #ifndef NDEBUG
191  // Caution: this is an O(L) operation, so it really should be done only in debug mode
192  for (int n = 0; n < base::northo; ++n)
193  maquis::cout << "MPS overlap: " << overlap(mps, base::ortho_mps[n]) << std::endl;
194 #endif
195 
196  maquis::cout << "Energy " << lr << " " << res.first << std::endl;
197  iteration_results_["Energy"] << res.first;
198 
199  double cutoff = this->get_cutoff(sweep);
200  std::size_t Mmax = this->get_Mmax(sweep);
201  truncation_results trunc;
202 
203  if (lr == +1)
204  {
205  // Write back result from optimization
206  boost::tie(mps[site1], mps[site2], trunc) = tst.split_mps_l2r(Mmax, cutoff);
207 
209 
210  //t = mps[site1].normalize_left(DefaultSolver());
211  //mps[site2].multiply_from_left(t);
212  //mps[site2].divide_by_scalar(mps[site2].scalar_norm());
213 
214  t = mps[site2].normalize_left(DefaultSolver());
215  // MD: DEBUGGING OUTPUT
216  maquis::cout << "Propagating t with norm " << t.norm() << std::endl;
217  if (site2 < L-1) mps[site2+1].multiply_from_left(t);
218 
219  if (site1 != L-2)
220  Storage::drop(right_[site2+1]);
221 
222  this->boundary_left_step(mpo, site1); // creating left_[site2]
223 
224  if (site1 != L-2){
225  Storage::evict(mps[site1]);
226  Storage::evict(left_[site1]);
227  }
228  }
229  if (lr == -1){
230  // Write back result from optimization
231  boost::tie(mps[site1], mps[site2], trunc) = tst.split_mps_r2l(Mmax, cutoff);
232 
234 
235  //t = mps[site2].normalize_right(DefaultSolver());
236  //mps[site1].multiply_from_right(t);
237  //mps[site1].divide_by_scalar(mps[site1].scalar_norm());
238 
239  t = mps[site1].normalize_right(DefaultSolver());
240  // MD: DEBUGGING OUTPUT
241  maquis::cout << "Propagating t with norm " << t.norm() << std::endl;
242  if (site1 > 0) mps[site1-1].multiply_from_right(t);
243 
244  if(site1 != 0)
245  Storage::drop(left_[site1]);
246 
247  this->boundary_right_step(mpo, site2); // creating right_[site2]
248 
249  if(site1 != 0){
250  Storage::evict(mps[site2]);
251  Storage::evict(right_[site2+1]);
252  }
253  }
254 
255  iteration_results_["BondDimension"] << trunc.bond_dimension;
256  iteration_results_["TruncatedWeight"] << trunc.truncated_weight;
257  iteration_results_["TruncatedFraction"] << trunc.truncated_fraction;
258  iteration_results_["SmallestEV"] << trunc.smallest_ev;
259 
260 
261  boost::chrono::high_resolution_clock::time_point sweep_then = boost::chrono::high_resolution_clock::now();
262  double elapsed = boost::chrono::duration<double>(sweep_then - sweep_now).count();
263  maquis::cout << "Sweep has been running for " << elapsed << " seconds." << std::endl;
264 
265  if (stop_callback())
266  throw dmrg::time_limit(sweep, _site+1);
267 
268  } // for sites
269  initial_site = -1;
270  } // sweep
std::size_t get_Mmax(int sweep) const
Definition: optimize.h:212
BaseParameters & parms
Definition: optimize.h:232
void sweep(int sweep, OptimizeDirection d=Both)
Definition: ts_optimize.hpp:71
std::pair< double, MPSTensor< Matrix, SymmGroup > > solve_ietl_lanczos(SiteProblem< Matrix, SymmGroup > &sp, MPSTensor< Matrix, SymmGroup > const &initial, BaseParameters &params)
std::vector< Boundary< typename storage::constrained< Matrix >::type, SymmGroup > > left_
Definition: optimize.h:235
int to_site(const int L, const int i) const
Definition: ts_optimize.hpp:65
Definition: optimize.h:84
static MPSTensor< Matrix, SymmGroup > site_ortho_boundaries(MPSTensor< Matrix, SymmGroup > const &mps, MPSTensor< Matrix, SymmGroup > const &ortho_mps, block_matrix< OtherMatrix, SymmGroup > const &ortho_left, block_matrix< OtherMatrix, SymmGroup > const &ortho_right)
Definition: contractions.h:479
std::vector< std::vector< block_matrix< typename storage::constrained< Matrix >::type, SymmGroup > > > ortho_left_
Definition: optimize.h:239
std::vector< int > get_right_placement(const MPOTensor< Matrix, SymmGroup > &mpo, const std::vector< int > &placement_l, const std::vector< int > &placement_r_old)
Definition: placement.h:85
results_collector iteration_results_
Definition: optimize.h:227
MPS< Matrix, SymmGroup >::scalar_type overlap(MPS< Matrix, SymmGroup > const &mps1, MPS< Matrix, SymmGroup > const &mps2)
Definition: mps_mpo_ops.h:136
std::vector< std::vector< block_matrix< typename storage::constrained< Matrix >::type, SymmGroup > > > ortho_right_
Definition: optimize.h:239
std::pair< double, MPSTensor< Matrix, SymmGroup > > solve_ietl_jcd(SiteProblem< Matrix, SymmGroup > &sp, MPSTensor< Matrix, SymmGroup > const &initial, BaseParameters &params, std::vector< MPSTensor< Matrix, SymmGroup > > ortho_vecs=std::vector< MPSTensor< Matrix, SymmGroup > >())
boost::function< bool()> stop_callback
Definition: optimize.h:233
MPO< Matrix, SymmGroup > const & mpo
Definition: optimize.h:230
MPS< Matrix, SymmGroup > & mps
Definition: optimize.h:229
double get_cutoff(int sweep) const
Definition: optimize.h:202
unsigned int northo
Definition: optimize.h:238
real_type norm() const
#define BEGIN_TIMING(name)
#define END_TIMING(name)
std::vector< Boundary< typename storage::constrained< Matrix >::type, SymmGroup > > right_
Definition: optimize.h:235
void boundary_right_step(MPO< Matrix, SymmGroup > const &mpo, int site)
Definition: optimize.h:144
std::vector< MPS< Matrix, SymmGroup > > ortho_mps
Definition: optimize.h:240
std::vector< int > get_left_placement(const MPOTensor< Matrix, SymmGroup > &mpo, const std::vector< int > &placement_l_old, const std::vector< int > &placement_r)
Definition: placement.h:69
void boundary_left_step(MPO< Matrix, SymmGroup > const &mpo, int site)
Definition: optimize.h:135
template<class Matrix , class SymmGroup , class Storage >
int ts_optimize< Matrix, SymmGroup, Storage >::to_site ( const int  L,
const int  i 
) const
inline

Definition at line 65 of file ts_optimize.hpp.

66  {
67  if (i < 0) return 0;
68  /// i, or (L-1) - (i - (L-1))
69  return (i < L-1) ? i : 2*L - 2 - i;
70  }

Member Data Documentation

template<class Matrix , class SymmGroup , class Storage >
results_collector optimizer_base< Matrix, SymmGroup, Storage >::iteration_results_
protectedinherited

Definition at line 227 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
std::vector<Boundary<typename storage::constrained<Matrix>::type, SymmGroup> > optimizer_base< Matrix, SymmGroup, Storage >::left_
protectedinherited

Definition at line 235 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
MPO<Matrix, SymmGroup> const& optimizer_base< Matrix, SymmGroup, Storage >::mpo
protectedinherited

Definition at line 230 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
MPS<Matrix, SymmGroup>& optimizer_base< Matrix, SymmGroup, Storage >::mps
protectedinherited

Definition at line 229 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
unsigned int optimizer_base< Matrix, SymmGroup, Storage >::northo
protectedinherited

Definition at line 238 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
std::vector< std::vector<block_matrix<typename storage::constrained<Matrix>::type, SymmGroup> > > optimizer_base< Matrix, SymmGroup, Storage >::ortho_left_
protectedinherited

Definition at line 239 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
std::vector<MPS<Matrix, SymmGroup> > optimizer_base< Matrix, SymmGroup, Storage >::ortho_mps
protectedinherited

Definition at line 240 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
std::vector< std::vector<block_matrix<typename storage::constrained<Matrix>::type, SymmGroup> > > optimizer_base< Matrix, SymmGroup, Storage >::ortho_right_
protectedinherited

Definition at line 239 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
BaseParameters& optimizer_base< Matrix, SymmGroup, Storage >::parms
protectedinherited

Definition at line 232 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
std::vector<Boundary<typename storage::constrained<Matrix>::type, SymmGroup> > optimizer_base< Matrix, SymmGroup, Storage >::right_
protectedinherited

Definition at line 235 of file optimize.h.

template<class Matrix , class SymmGroup , class Storage >
boost::function<bool ()> optimizer_base< Matrix, SymmGroup, Storage >::stop_callback
protectedinherited

Definition at line 233 of file optimize.h.


The documentation for this class was generated from the following file: