I have a system with two different types of degrees of freedom - one is high energy, and the other low energy. Therefore a good description for the low energy eigenstates of the system can probably be obtained by only keeping a few basis states for the high energy degrees of freedom, but keeping more states for the low energy degrees of freedom. I'm wondering how to implement such a thing? Here is my existing custom class. I'm envisioning having adim and edim, as opposed to just dim, and then constructing new operators based on those dimensions. Thanks!
namespace itensor {
class ExcitonSite;
using Exciton = BasicSiteSet<ExcitonSite>;
class ExcitonSite
{
IQIndex s;
int dim;
public:
ExcitonSite() { }
ExcitonSite(IQIndex I) : s(I) { }
ExcitonSite(int n, Args const& args = Args::global())
{
dim = args.getInt("on_site_dim", 5);
auto v = stdx::reserve_vector<IndexQN>(2 * dim + 1);
for (int j = -dim; j <= dim; ++j)
{
auto i = Index(std::to_string(j), 1, Site);
auto q = QN("Sz=", j);
v.emplace_back(i, q);
}
s = IQIndex(nameint("site=",n), std::move(v), Out, 0);
}
IQIndex
index() const { return s; }
IQIndexVal
state(std::string const& state)
{
int j = -dim;
while (j <= dim)
{
j++;
if (state == std::to_string(j))
{
return s(j + 1 + dim);
break;
}
}
if (j == dim+1)
{
Error("State " + state + " not recognized");
}
return IQIndexVal{};
}
IQTensor
op(std::string const& opname,
Args const& args) const
{
auto sP = prime(s);
auto Op = IQTensor(dag(s),sP);
if (opname == "n")
{
for (int j = -dim; j<= dim; ++j)
{
Op.set(s(j + 1 + dim), sP(j + 1 + dim), (float)j);
}
}
else
if (opname == "nsq")
{
for (int j = -dim; j<= dim; ++j)
{
Op.set(s(j + 1 + dim), sP(j + 1 + dim), (float)(j * j));
}
}
else
if (opname == "gp")
{
for (int j = -dim; j<= dim - 1; ++j)
{
Op.set(s(j + 1 + dim), sP(j + 2 + dim), +1.0);
}
}
else
if (opname == "igp")
{
for (int j = -dim; j<= dim-1; ++j)
{
Op.set(s(j + 1 + dim), sP(j + 2 + dim), +1.0_i);
}
}
else
if (opname == "gm")
{
for (int j = -dim + 1; j<= dim; ++j)
{
Op.set(s(j + 1 + dim), sP(j + dim), +1.0);
}
}
else
if (opname == "igm")
{
for (int j = -dim + 1; j<= dim; ++j)
{
Op.set(s(j + 1 + dim), sP(j + dim), +1.0_i);
}
}
else
{
Error("Operator \"" + opname + "\" name not recognized");
}
return Op;
}
};
} //namespace itensor
#endif