Hi Miles,
I am calculating the Entanglement Entropy(EE) of a Bose-Fermi mixture chain. There are two degenerated ground states which belong to two different QN blocks ( The Bose particle number and Fermi particle number are both conservation.), one is located in the (Nb,Nf)=(22,0) Hilbert sub-space and the other is in the (Nb,Nf)=(21,1) space. The total number of sites is 36. Besides, I am using the IQMPS in the DMRG calculation.
I can do the EE calculation of the single ground state and it worked well. But I think the truly ground state |phi> is a linear combination of these two states and the coefficients can be determined by minimizing the EE (This way can be found from D.N.Sheng's paper). I try to evaluate these values and the key function used in my program is .plusEq(psi1,{"Maxm",3000,"Cutoff",1E-9}). The following is my code:
#include "itensor/all.h"
#include <fstream>
#include "itensor/decomp.h"
#include<math.h>
using namespace std;
using namespace itensor;
int main(int argc, char* argv[])
{
//Parse the input file //DMRG parameter
if(argc != 2) { printfln("Usage: %s inputfile_bsfm",argv[0]); return 0; }
auto input = InputGroup(argv[1],"input");
auto N = input.getInt("N"); // auto N=80;
auto sites2=BF(N);
readFromFile("sites_NbNf220", sites2);
IQMPS psi2(sites2);
readFromFile("psi_NbNf220", psi2);
auto sites3=BF(N);
readFromFile("sites_NbNf211", sites3);
IQMPS psi3(sites3);
readFromFile("psi_NbNf211", psi3);
ofstream SaveFile("EE_mixed.dat");
auto pi=4*atan(1.0);
std::complex<double> z {0, 1};
auto imag=z;
auto b=36;
auto c1=0.5;
auto c2=sqrt(1-pow(c1,2));
auto fai=pi;
psi2 *= c1;
psi3 *= c2*exp(imag*fai);
auto psi=toMPS(psi2);
auto psi1=toMPS(psi3);
println("Add two states ");
// psi=psi*c1+psi1*c2*exp(i*fai)
psi.plusEq(psi1,{"Maxm",3000,"Cutoff",1E-9});
//Entanglement Entropy
printfln("\n");
println("Entanglement Entropy");
//Given an MPS or IQMPS called "psi",
//and some particular bond "b" (1 <= b < psi.N())
//across which we want to compute the von Neumann entanglement
//"Gauge" the MPS to site b
if(b%2==0){
psi.position(b);
//Here assuming an MPS of ITensors, but same code works
//for IQMPS by replacing ITensor -> IQTensor
//Compute two-site wavefunction for sites (b,b+1)
// IQTensor wf = psi.A(b)*psi.A(b+1);
ITensor wf = psi.A(b)*psi.A(b+1);
//SVD this wavefunction to get the spectrum
//of density-matrix eigenvalues
auto U = psi.A(b);
// IQTensor S,V;
ITensor S,V;
auto spectrum = svd(wf,U,S,V);
//Apply von Neumann formula
//spectrum.eigs() is a Vector containing
//the density matrix eigenvalues
//(squares of the singular values)
Real SvN = 0.;
for(auto p : spectrum.eigs())
{
if(p > 1E-12) SvN += -p*log(p);
}
SaveFile<<b/2<<" "<<SvN<<endl;
}
println("\n");
SaveFile.close();
return 0;
}
It can be built but something wrong when I run the code. The error message is ITensor has different index structure. I believe that the error occurred at the psi.plusEq(psi1,{"Maxm",3000,"Cutoff",1E-9}); line, because the message "Add two states" has been output on the screen but there is no "Entanglement Entropy" string. So, my question is whether it's possible to do this or are there some tricks?
Thank you.
Best,
Chenrong