Skip to content
Snippets Groups Projects
Commit 9d96216e authored by Joao Dos Santos Ferreira's avatar Joao Dos Santos Ferreira
Browse files

Delete auxiliar.h

parent cd83526e
Branches master
No related tags found
No related merge requests found
//
// My version of TEVOL, the args "Swipe"="left" or "right" defines the orientation.
//
#ifndef __ITENSOR_3TEBD_H
#define __ITENSOR_3TEBD_H
#include "itensor/all.h"
#include "itensor/mps/mpo.h"
#include "itensor/mps/TEvolObserver.h"
using namespace std;
using namespace itensor;
class BondGate3{
public:
enum Type { tReal, tImag};
BondGate3(SiteSet const& sites,int i1,int i2,int i3,Type type,Cplx tau,ITensor bondH,string which_swipe);
BondGate3(SiteSet const& sites,int i1,int i2,int i3,Type type,Real tau,ITensor bondH,string which_swipe);
int i1() const { return i1_; }
int i2() const { return i2_; }
int i3() const { return i3_; }
operator const ITensor&() const { return gate_; }
ITensor const&
gate() const { return gate_; }
Type
type() const { return type_; }
void
set_swipe(string tt) { swipe_=tt; }
string get_swipe() const {return swipe_;}
private:
Type type_;
int i1_,i2_,i3_; // The left, right indices of bond
ITensor gate_;
string swipe_; //
};
template <class Iterable>
Real
gateTEvolMixed(Iterable const& gatelist,Real ttotal,Real tstep,MPS & psi,Args args);
//
//
// Implementations
//
ITensor inline
operator*(BondGate3 const& G, ITensor T) { T *= G.gate(); return T; }
ITensor inline
operator*(ITensor T, BondGate3 const& G) { T *= G.gate(); return T; }
//Receives a three site tensor
inline BondGate3::
BondGate3(SiteSet const& sites,
int i1,
int i2,
int i3,
Type type,
Cplx tau,
ITensor bondH,
string which_swipe)
: type_(type)
{
swipe_=which_swipe;
vector<int> ax{i1,i2,i3};
std::sort(ax.begin(), ax.end());
i1_=ax[0];
i2_=ax[1];
i3_=ax[2];
if(!(type_ == tReal || type_ ==tImag))
{
Error("When providing bondH, type must be tReal or tImag");
}
bondH *= -tau;
ITensor unit = sites.op("Id",i1_)*sites.op("Id",i2_)*sites.op("Id",i3_);
if(type_ == tReal)
{
bondH *= Complex_i;
}
auto term = bondH;
bondH.replaceTags("1","2");
bondH.replaceTags("0","1");
// exp(x) = 1 + x + x^2/2! + x^3/3! ..
// = 1 + x * (1 + x/2 *(1 + x/3 * (...
// ~ ((x/3 + 1) * x/2 + 1) * x + 1
for(int ord = 100; ord >= 1; --ord)
{
term /= ord;
gate_ = unit + term;
term = gate_ * bondH;
term.replaceTags("2","1");
}
}
//Receives a three site tensor
inline BondGate3::
BondGate3(SiteSet const& sites,
int i1,
int i2,
int i3,
Type type,
Real tau,
ITensor bondH,
string which_swipe)
: type_(type)
{
swipe_=which_swipe;
vector<int> ax{i1,i2,i3};
std::sort(ax.begin(), ax.end());
i1_=ax[0];
i2_=ax[1];
i3_=ax[2];
if(!(type_ == tReal || type_ ==tImag))
{
Error("When providing bondH, type must be tReal or tImag");
}
bondH *= -tau;
ITensor unit = sites.op("Id",i1_)*sites.op("Id",i2_)*sites.op("Id",i3_);
if(type_ == tReal)
{
bondH *= Complex_i;
}
auto term = bondH;
bondH.replaceTags("1","2");
bondH.replaceTags("0","1");
// exp(x) = 1 + x + x^2/2! + x^3/3! ..
// = 1 + x * (1 + x/2 *(1 + x/3 * (...
// ~ ((x/3 + 1) * x/2 + 1) * x + 1
for(int ord = 100; ord >= 1; --ord)
{
term /= ord;
gate_ = unit + term;
term = gate_ * bondH;
term.replaceTags("2","1");
}
}
template <class Iterable>
Real
gateTEvolMixed(Iterable const& gatelist,Real ttotal,Real tstep,MPS & psi,Args args){
TEvolObserver obs(args);
const bool verbose = args.getBool("Verbose",false);
const bool do_normalize = args.getBool("Normalize",true);
const string which_swipe = args.getString("Swipe","error");
const int nt = int(ttotal/tstep+(1e-9*(ttotal/tstep)));
if(fabs(nt*tstep-ttotal) > 1E-9){
Error("Timestep not commensurate with total time");
}
if(verbose){
printfln("Taking %d steps of timestep %.5f, total time %.5f",nt,tstep,ttotal);
}
psi.position(gatelist.front().i1());
Real tot_norm = norm(psi);
Real tsofar = 0;
for(auto tt : range1(nt)){
auto g = gatelist.begin();
while(g != gatelist.end())
{
auto i1 = g->i1();
auto i2 = g->i2();
auto i3 = g->i3();
if(i3){
//print("*"+str(i1)+"-"+str(i2)+"-"+str(i3)+"*");
auto AA=psi(i1)*psi(i2)*psi(i3)*g->gate();
AA.replaceTags("Site,1","Site,0");
auto swp=g->get_swipe();
auto link_tags1 = tags(linkIndex(psi,i1));
auto link_tags2 = tags(linkIndex(psi,i2));
++g;
if(g != gatelist.end()){
//Look ahead to next gate position
auto nr = g->i1();
auto nl = g->i3();
if(!nl)
nl=g->i2();
if(swp=="right")
{
auto [U,S,V] = svd(AA,inds(psi(i1)),args);
if(args.getBool("DoNormalize",false))
S *= 1./itensor::norm(S);
auto ulink=commonIndex(U,S);
auto iset=IndexSet(ulink,siteIndex(psi,i2));
auto [U2,S2,V2] = svd(S*V,iset,args);
if(args.getBool("DoNormalize",false))
S2 *= 1./itensor::norm(S2);
psi.set(i1,U);
psi.set(i2,U2);
psi.set(i3,S2*V2);
//psi.svdBond(i1,AA,Fromleft,args);
auto ulink1=commonIndex(psi(i1),psi(i2));
auto ulink2=commonIndex(psi(i2),psi(i3));
psi.ref(i1).setTags(link_tags1,ulink1);
psi.ref(i2).setTags(link_tags1,ulink1);
psi.ref(i2).setTags(link_tags2,ulink2);
psi.ref(i3).setTags(link_tags2,ulink2);
psi.position(nr);
//does no work if position already ni1
}
else if(swp=="left")
{
auto [U,S,V] = svd(AA,inds(psi(i3)),args);
if(args.getBool("DoNormalize",false))
S *= 1./itensor::norm(S);
auto ulink=commonIndex(U,S);
auto iset=IndexSet(ulink,siteIndex(psi,i2));
auto [U2,S2,V2] = svd(S*V,iset,args);
if(args.getBool("DoNormalize",false))
S2 *= 1./itensor::norm(S2);
psi.set(i3,U);
psi.set(i2,U2);
psi.set(i1,S2*V2);
//psi.svdBond(i1,AA,Fromleft,args);
auto ulink1=commonIndex(psi(i1),psi(i2));
auto ulink2=commonIndex(psi(i2),psi(i3));
psi.ref(i1).setTags(link_tags1,ulink1);
psi.ref(i2).setTags(link_tags1,ulink1);
psi.ref(i2).setTags(link_tags2,ulink2);
psi.ref(i3).setTags(link_tags2,ulink2);
psi.position(nl);
}
else{
Error("You must define swipe direction - left or right");
}
}
else{
auto [U,S,V] = svd(AA,inds(psi(i3)),args);
if(args.getBool("DoNormalize",false))
S *= 1./itensor::norm(S);
auto ulink=commonIndex(U,S);
auto iset=IndexSet(ulink,siteIndex(psi,i2));
auto [U2,S2,V2] = svd(S*V,iset,args);
if(args.getBool("DoNormalize",false))
S2 *= 1./itensor::norm(S2);
psi.set(i3,U);
psi.set(i2,U2);
psi.set(i1,S2*V2);
auto ulink1=commonIndex(psi(i1),psi(i2));
auto ulink2=commonIndex(psi(i2),psi(i3));
psi.ref(i1).setTags(link_tags1,ulink1);
psi.ref(i2).setTags(link_tags1,ulink1);
psi.ref(i2).setTags(link_tags2,ulink2);
psi.ref(i3).setTags(link_tags2,ulink2);
psi.position(i1);
}
}
else{
auto AA = psi(i1)*psi(i2)*g->gate();
AA.replaceTags("Site,1","Site,0");
auto swp=g->get_swipe();
++g;
if(g != gatelist.end())
{
//Look ahead to next gate position
auto nr = g->i1();
auto nl = g->i3();
if(!nl)
nl=g->i2();
//SVD AA to restore MPS form
//before applying current gate
if(swp=="right")
{
psi.svdBond(i1,AA,Fromleft,args);
psi.position(nr); //does no work if position already ni1
}
else if(swp=="left")
{
psi.svdBond(i1,AA,Fromright,args);
psi.position(nl); //does no work if position already ni2
}
else{
Error("You must define swipe direction - left or right");
}
}
else
{
//No next gate to analyze, just restore MPS form
psi.svdBond(i1,AA,Fromright,args);
}
}
}
if(do_normalize)
{
tot_norm *= psi.normalize();
}
tsofar += tstep;
args.add("TimeStepNum",tt);
args.add("Time",tsofar);
args.add("TotalTime",ttotal);
obs.measure(args);
}
if(verbose)
{
printfln("\nTotal time evolved = %.5f\n",tsofar);
}
return tot_norm;
} // gateTEvol
#endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment