Hello everyone!
I want to measure the time-evolving wavefunction after every time-step.
So, I modified these two head files TEvolObserver.h and tevol.h for this purpose.
In short, I just added one more SiteSet parameter to TEvolObserver type and a MPS type parameter to measure, so that TEvolObserver can do measurements when it obtain both wavefunciton and sites information.
Here are the details.
For TEvolObserver.h, I added a SiteSet type parameter, so that site information can be passed in and I can construct operator op(site,"",i) and do measurements in the file.
c++
#ifndef __ITENSOR_TEVOLOBSERVER_H
#define __ITENSOR_TEVOLOBSERVER_H
#include "itensor/util/readwrite.h"
#include "itensor/mps/observer.h"
namespace itensor {
class TEvolObserver : public Observer
{
public:
using String = std::string;
TEvolObserver(SiteSet const& sites,
Args const& args = Args::global());
virtual ~TEvolObserver() { }
void virtual
measure(MPS const& psi,Args const& args = Args::global());
bool virtual
checkDone(Args const& args = Args::global());
private:
/////////////
//
// Data Members
SiteSet const& site;
bool done_,
show_percent_,
m_l_d;
//
/////////////
}; // class TEvolObserver
inline TEvolObserver::
TEvolObserver(SiteSet const& sites,Args const& args)
:
site(sites),
done_(false),
show_percent_(args.getBool("ShowPercent",true)),
m_l_d(args.getBool("MLD",true))
{
}
void inline TEvolObserver::
measure(MPS const& psi, const Args& args)
{
auto psi_=psi;
const Real t = args.getReal("Time");
if(show_percent_)
{
const Real ttotal = args.getReal("TotalTime");
Real percentdone = (100.*t)/ttotal;
if(percentdone < 99.5 || (std::fabs(t-ttotal) < 1E-10))
{
printf("\b\b\b%2.f%%",percentdone);
std::cout.flush();
}
}
auto N = length(site);
if(m_l_d)
{
String opname="N";
for(int b=1;b <= N; b++)
{
auto ket = psi_(b);
auto bra = dag(prime(ket,"Site"));
auto Numop = op(site,opname,b);
auto Numi = std::real(eltC(bra*Numop*ket));
printfln("Observable %s of site %i at time %d is %f\n",opname,b,t,Numi);
}
}
}
bool inline TEvolObserver::
checkDone(const Args& args)
{
const Real t = args.getReal("Time");
if(fileExists("STOP_TEVOL"))
{
println("File STOP_TEVOL found: stopping this time evolution run at time ",t);
std::remove("STOP_TEVOL");
return true;
}
//Set done_ flag to true so any outer callers using this Observer will also terminate.
if(fileExists("STOP_TEVOL_ALL"))
{
println("File STOP_TEVOL_ALL found: stopping this time evolution at time ",t);
std::remove("STOP_TEVOL_ALL");
done_ = true;
return done_;
}
return done_;
}
} //namespace itensor
#endif // __ITENSOR_TEVOLOBSERVER_H
Furthermore, I made a little bit modification in tevol.h so that wavefunction psi can be passed to measure method.
c++
#ifndef __ITENSOR_TEVOL_H
#define __ITENSOR_TEVOL_H
#include "itensor/mps/mpo.h"
#include "itensor/mps/bondgate.h"
#include "itensor/mps/TEvolObserver.h"
namespace itensor {
template <class Iterable>
Real
gateTEvol(Iterable const& gatelist,
Real ttotal,
Real tstep,
MPS & psi,
Args const& args = Args::global());
template <class Iterable>
Real
gateTEvol(Iterable const& gatelist,
Real ttotal,
Real tstep,
MPS & psi,
TEvolObserver& obs,
Args args = Args::global());
//
//
// Implementations
//
template <class Iterable>
Real
gateTEvol(Iterable const& gatelist,
Real ttotal,
Real tstep,
MPS & psi,
TEvolObserver& obs,
Args args)
{
const bool verbose = args.getBool("Verbose",false);
const bool do_normalize = args.getBool("Normalize",true);
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 AA = psi(i1)*psi(i2)*g->gate();
AA.replaceTags("Site,1","Site,0");
++g;
if(g != gatelist.end())
{
//Look ahead to next gate position
auto ni1 = g->i1();
auto ni2 = g->i2();
//SVD AA to restore MPS form
//before applying current gate
if(ni1 >= i2)
{
psi.svdBond(i1,AA,Fromleft,args);
psi.position(ni1); //does no work if position already ni1
}
else
{
psi.svdBond(i1,AA,Fromright,args);
psi.position(ni2); //does no work if position already ni2
}
}
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(psi,args);
}
if(verbose)
{
printfln("\nTotal time evolved = %.5f\n",tsofar);
}
return tot_norm;
} // gateTEvol
template <class Iterable>
Real
gateTEvol(Iterable const& gatelist,
Real ttotal,
Real tstep,
MPS & psi,
Args const& args)
{
TEvolObserver obs(psi,args);
return gateTEvol(gatelist,ttotal,tstep,psi,obs,args);
}
} //namespace itensor
#endif
These modification works well at least in the compliling state.
However, here comes many incorrect results in the output file, say
Observable N of site 3 at time 1.5 is 16.828193
Observable N of site 4 at time 1.5 is 39.418237
Observable N of site 5 at time 1.5 is 66.264575
Observable N of site 6 at time 1.5 is 60.269621
Observable N of site 7 at time 1.5 is 56.621804
Here is my site definition:
c++
auto sites=Boson(L,{"ConserveNb",true,"ConserveNb",true,"MaxOcc=",5});
The local density should less than 4.
So, does the problem attribute to the modification of head files?