Commit b5ad647c authored by Alvarez, Gonzalo's avatar Alvarez, Gonzalo
Browse files

KroneckerDumper: unhooked

parent 776137fd
......@@ -317,13 +317,12 @@ private:
if (tstStruct_.omega().first != PsimagLite::FREQ_REAL)
throw PsimagLite::RuntimeError("Matsubara only with KRYLOV\n");
RealType fakeTime = 0;
const RealType fakeTime = 0;
typename ModelHelperType::Aux aux(p, lrs_);
typename ModelType::HamiltonianConnectionType hc(lrs_,
model_.superGeometry(),
ModelType::modelLinks(),
fakeTime,
0);
model_.superOpHelper());
LanczosMatrixType h(model_, hc, aux);
RealType E0 = energy_;
CorrectionVectorFunctionType cvft(h,tstStruct_,E0);
......
......@@ -111,7 +111,6 @@ public:
typedef typename SparseMatrixType::value_type ComplexOrRealType;
typedef typename ModelType::ModelHelperType ModelHelperType;
typedef typename ModelHelperType::LeftRightSuperType LeftRightSuperType;
typedef typename LeftRightSuperType::ParamsForKroneckerDumperType ParamsForKroneckerDumperType;
typedef typename TargetingType::MatrixVectorType MatrixVectorType;
typedef typename ModelType::InputValidatorType InputValidatorType;
typedef typename PsimagLite::Vector<RealType>::Type VectorRealType;
......@@ -459,20 +458,11 @@ private:
SizeType loopIndex)
{
const OptionsType& options = parameters_.options;
const bool dumperEnabled = options.isSet("KroneckerDumper");
ParamsForKroneckerDumperType paramsKrDumper(dumperEnabled,
parameters_.dumperBegin,
parameters_.dumperEnd,
parameters_.precision);
ParamsForKroneckerDumperType* paramsKrDumperPtr = 0;
if (lrs.super().block().size() == model_.superGeometry().numberOfSites())
paramsKrDumperPtr = &paramsKrDumper;
HamiltonianConnectionType hc(lrs,
model_.superGeometry(),
ModelType::modelLinks(),
targetTime,
paramsKrDumperPtr);
model_.superOpHelper());
const SizeType saveOption = parameters_.finiteLoop[loopIndex].saveOption;
typename ModelHelperType::Aux aux(partitionIndex, lrs);
......
......@@ -88,17 +88,19 @@ DISCLOSED WOULD NOT INFRINGE PRIVATELY OWNED RIGHTS.
#include "OperatorStorage.h"
#include "OperatorsCached.h"
#include "ManyToTwoConnection.h"
#include "SuperOpHelperBase.h"
namespace Dmrg {
// Keep this class independent of x and y in x = H*y
// For things that depend on x and y use ParallelHamiltonianConnection.h
template<typename ModelLinksType, typename ModelHelperType_>
template<typename ModelLinksType, typename ModelHelperType_, typename ParamsForSolverType_>
class HamiltonianConnection {
public:
typedef ModelHelperType_ ModelHelperType;
typedef ParamsForSolverType_ ParamsForSolverType;
typedef typename ModelLinksType::SuperGeometryType SuperGeometryType;
typedef HamiltonianAbstract<SuperGeometryType> HamiltonianAbstractType;
typedef typename ModelHelperType::RealType RealType;
......@@ -113,24 +115,23 @@ public:
typedef typename PsimagLite::Concurrency ConcurrencyType;
typedef PsimagLite::Vector<SizeType>::Type VectorSizeType;
typedef typename ModelHelperType::LeftRightSuperType LeftRightSuperType;
typedef typename LeftRightSuperType::ParamsForKroneckerDumperType ParamsForKroneckerDumperType;
typedef typename LeftRightSuperType::KroneckerDumperType KroneckerDumperType;
typedef typename PsimagLite::Vector<LinkType>::Type VectorLinkType;
typedef typename ModelHelperType::Aux AuxType;
typedef OperatorsCached<LeftRightSuperType> OperatorsCachedType;
typedef typename ModelLinksType::TermType::OneLinkType OneLinkType;
typedef ManyToTwoConnection<ModelLinksType, LeftRightSuperType> ManyToTwoConnectionType;
typedef SuperOpHelperBase<SuperGeometryType, ParamsForSolverType> SuperOpHelperBaseType;
typedef ManyToTwoConnection<ModelLinksType, LeftRightSuperType, SuperOpHelperBaseType>
ManyToTwoConnectionType;
HamiltonianConnection(const LeftRightSuperType& lrs,
const SuperGeometryType& superGeometry,
const ModelLinksType& lpb,
RealType targetTime,
const ParamsForKroneckerDumperType* pKroneckerDumper)
const SuperOpHelperBaseType& superOpHelper)
: modelHelper_(lrs),
superGeometry_(superGeometry),
superGeometry_(superOpHelper.superGeometry()),
modelLinks_(lpb),
targetTime_(targetTime),
kroneckerDumper_(pKroneckerDumper, lrs),
superOpHelper_(superOpHelper),
operatorsCached_(lrs),
progress_("HamiltonianConnection"),
systemBlock_(modelHelper_.leftRightSuper().left().block()),
......@@ -151,7 +152,7 @@ public:
SizeType last = lrs.super().block().size();
assert(last > 0);
--last;
SizeType numberOfSites = superGeometry.numberOfSites();
SizeType numberOfSites = superGeometry_.numberOfSites();
assert(numberOfSites > 0);
bool superIsReallySuper = (lrs.super().block()[0] == 0 &&
lrs.super().block()[last] == numberOfSites - 1);
......@@ -251,11 +252,6 @@ public:
return link2;
}
KroneckerDumperType& kroneckerDumper() const
{
return kroneckerDumper_;
}
const ModelHelperType& modelHelper() const { return modelHelper_; }
SizeType tasks() const {return lps_.size(); }
......@@ -309,7 +305,8 @@ private:
ManyToTwoConnectionType manyToTwo(hItems,
type,
oneLink,
modelHelper_.leftRightSuper());
modelHelper_.leftRightSuper(),
superOpHelper_);
LinkType link2(manyToTwo.finalIndices(),
type,
......@@ -355,7 +352,7 @@ private:
const SuperGeometryType& superGeometry_;
const ModelLinksType& modelLinks_;
RealType targetTime_;
mutable KroneckerDumperType kroneckerDumper_;
const SuperOpHelperBaseType& superOpHelper_;
OperatorsCachedType operatorsCached_;
PsimagLite::ProgressIndicator progress_;
VectorLinkType lps_;
......
......@@ -9,7 +9,7 @@
namespace Dmrg {
template<typename LeftRightSuperType>
template<typename LeftRightSuperType, typename SolverParamsType>
class KroneckerDumper {
typedef PsimagLite::Concurrency ConcurrencyType;
......@@ -42,20 +42,29 @@ public:
SizeType nOfQns;
}; // struct ParamsForKroneckerDumper
KroneckerDumper(const ParamsForKroneckerDumper* p,
const LeftRightSuperType& lrs)
: enabled_(p && p->enabled),pairCount_(0),disable_(false)
KroneckerDumper(const SolverParamsType& params,
const LeftRightSuperType& lrs,
ProgramGlobals::DirectionEnum dir)
: enabled_(false),pairCount_(0),disable_(false)
{
if (!enabled_) return;
if (dir == ProgramGlobals::DirectionEnum::INFINITE)
return;
bool b = (p->end > 0 && counter_ >= p->end);
if (counter_ < p->begin || b) {
counter_++;
enabled_ = params.options.isSet("KroneckerDumper");
if (!enabled_) return;
ParamsForKroneckerDumper p(enabled_,
params.dumperBegin,
params.dumperEnd,
params.precision);
bool b = (p.end > 0 && counter_ >= p.end);
if (counter_ < p.begin || b) {
++counter_;
enabled_ = false;
return;
}
if (p->nOfQns == 0) {
if (p.nOfQns == 0) {
PsimagLite::String msg("KroneckerDumper::ctor(): internal error ");
throw PsimagLite::RuntimeError(msg + "nOfQns\n");
}
......@@ -64,13 +73,13 @@ public:
PsimagLite::String filename = "kroneckerDumper" + ttos(counter_) + ".txt";
fout_.open(filename.c_str());
fout_.precision(p->precision);
fout_.precision(p.precision);
fout_<<"KroneckerDumper for DMRG++ version "<<DMRGPP_VERSION<<"\n";
fout_<<"Instance="<<counter_<<"\n";
fout_<<"EncodingOfQuantumNumbers="<<(2*ProgramGlobals::maxElectronsOneSpin)<<"\n";
printOneBasis("Left",lrs.left(),p->nOfQns);
printOneBasis("Right",lrs.right(),p->nOfQns);
printOneBasis("Left",lrs.left(),p.nOfQns);
printOneBasis("Right",lrs.right(),p.nOfQns);
fout_<<"SuperBasisPermutation\n";
fout_<<lrs.super().permutationVector();
......@@ -78,7 +87,7 @@ public:
//fout_<<qtarget<<"\n";
signs_ = lrs.left().signs();
counter_++;
++counter_;
}
~KroneckerDumper()
......@@ -217,8 +226,8 @@ private:
ConcurrencyType::MutexType mutex_;
}; // class KroneckerDumpter
template<typename SparseMatrixType>
SizeType KroneckerDumper<SparseMatrixType>::counter_ = 0;
template<typename T1, typename T2>
SizeType KroneckerDumper<T1, T2>::counter_ = 0;
} // namespace Dmrg
#endif // KRONECKERDUMPER_H
......@@ -80,7 +80,6 @@ DISCLOSED WOULD NOT INFRINGE PRIVATELY OWNED RIGHTS.
#define LEFT_RIGHT_SUPER_H
#include "ProgressIndicator.h"
#include "KroneckerDumper.h"
#include "Io/IoNg.h"
namespace Dmrg {
......@@ -100,8 +99,6 @@ public:
typedef typename BasisType::BlockType BlockType;
typedef PsimagLite::ProgressIndicator ProgressIndicatorType;
typedef LeftRightSuper<BasisWithOperatorsType_,SuperBlockType> ThisType;
typedef KroneckerDumper<ThisType> KroneckerDumperType;
typedef typename KroneckerDumperType::ParamsForKroneckerDumper ParamsForKroneckerDumperType;
typedef typename BasisType::QnType QnType;
template<typename IoInputter>
......@@ -365,7 +362,11 @@ private:
typedef LeftRightSuper<BasisWithOperatorsType, BasisType> LeftRightSuper2Type;
Xbasis.setOneSite(X, model, time);
leftOrRight.setToProduct(pS, Xbasis, model.superOpHelper(pS.block(), X, dir));
assert(X.size() == 1);
SizeType lastS = pS.block().size();
assert(lastS > 0);
model.superOpHelper().setToProduct(pS.block()[--lastS], X[0], dir);
leftOrRight.setToProduct(pS, Xbasis, model.superOpHelper());
SparseMatrixType matrix = leftOrRight.hamiltonian().getCRS();
......
......@@ -4,7 +4,9 @@
namespace Dmrg {
template<typename ModelLinksType, typename LeftRightSuperType>
template<typename ModelLinksType,
typename LeftRightSuperType,
typename SuperOpHelperType>
class ManyToTwoConnection {
public:
......@@ -18,10 +20,11 @@ public:
ManyToTwoConnection(const VectorSizeType& hItems,
ProgramGlobals::ConnectionEnum type,
const ModelTermLinkType& oneLink,
const LeftRightSuperType& lrs)
const LeftRightSuperType& lrs,
const SuperOpHelperType& superOpHelper)
: oneLink_(oneLink), lrs_(lrs)
{
finalIndices_ = finalIndices(hItems, type);
finalIndices_ = finalIndices(hItems, type, superOpHelper);
assert(oneLink.mods.size() == 2);
mods_ = PairCharType(oneLink.mods[0], oneLink.mods[1]);
}
......@@ -70,11 +73,12 @@ private:
}
PairSizeType finalIndices(const VectorSizeType& hItems,
ProgramGlobals::ConnectionEnum type) const
ProgramGlobals::ConnectionEnum type,
const SuperOpHelperType& superOpHelper) const
{
if (hItems.size() == 4) {
assert(type == ProgramGlobals::ConnectionEnum::SYSTEM_ENVIRON);
return finalIndices4sites(hItems);
return finalIndices4sites(hItems, superOpHelper);
}
const ProgramGlobals::SysOrEnvEnum sysOrEnv =
......@@ -115,7 +119,8 @@ private:
lrs_.right().localOperatorIndex(i, sigma);
}
PairSizeType finalIndices4sites(const VectorSizeType& hItems) const
PairSizeType finalIndices4sites(const VectorSizeType& hItems,
const SuperOpHelperType& superOpHelper) const
{
SizeType last = lrs_.left().block().size();
assert(last > 0);
......@@ -133,9 +138,9 @@ private:
assert(sysSites.size() > 0);
assert(envSites.size() > 0);
PairSizeType finalIndex(lrs_.left(). superOperatorIndices(sysSites, oneLink_.indices[0]),
lrs_.right().superOperatorIndices(envSites, oneLink_.indices[1]));
return finalIndex;
const SizeType finalLeft = superOpHelper.leftIndex(sysSites, oneLink_.indices[0]);
const SizeType finalRight = superOpHelper.leftIndex(sysSites, oneLink_.indices[1]);
return PairSizeType(finalLeft, finalRight);
}
const ModelTermLinkType& oneLink_;
......
......@@ -157,7 +157,7 @@ public:
typedef typename ModelLinksType::TermType ModelTermType;
typedef OpaqueOp OpForLinkType;
typedef typename ModelLinksType::AtomKindBase AtomKindBaseType;
typedef SuperOpHelperBase SuperOpHelperType;
typedef SuperOpHelperBase<SuperGeometryType, ParametersType> SuperOpHelperBaseType;
ModelBase(const ParametersType& params,
const SuperGeometryType& superGeometry,
......@@ -181,6 +181,9 @@ public:
customOperators();
modelLinks_.postCtor2();
superOpHelper_ = setSuperOpHelper();
assert(superOpHelper_);
ProgramGlobals::init(maxElectronsOneSpin());
}
......@@ -338,19 +341,6 @@ for (SizeType dof = 0; dof < numberOfDofs; ++dof) {
virtual PsimagLite::String oracle() const { return ""; }
virtual const SuperOpHelperType& superOpHelper(const VectorSizeType& bigBlock,
const VectorSizeType& smallBlock,
ProgramGlobals::DirectionEnum dir) const
{
if (superOpHelper_) {
delete superOpHelper_;
superOpHelper_ = nullptr;
}
superOpHelper_ = new SuperOpHelperType(bigBlock, smallBlock, dir);
return *superOpHelper_;
}
/**
The function \cppFunction{addHamiltonianConnection} implements
the Hamiltonian connection (e.g. tight-binding links in the case of the Hubbard Model
......@@ -375,10 +365,9 @@ for (SizeType dof = 0; dof < numberOfDofs; ++dof) {
VectorSizeType nzs(total, 0);
HamiltonianConnectionType hc(lrs,
modelCommon_.superGeometry(),
modelLinks_,
currentTime,
0);
superOpHelper());
for (SizeType m = 0; m < total; ++m) {
SizeType offset = lrs.super().partition(m);
......@@ -601,6 +590,11 @@ for (SizeType dof = 0; dof < numberOfDofs; ++dof) {
InputValidatorType_& ioIn() const { return ioIn_; }
SuperOpHelperBaseType& superOpHelper() const
{
return *superOpHelper_;
}
protected:
PsimagLite::String oracle(const RealType& energy,
......@@ -629,6 +623,11 @@ protected:
return modelLinks_.createTerm(name);
}
virtual SuperOpHelperBaseType* setSuperOpHelper()
{
return new SuperOpHelperBaseType(modelCommon_.superGeometry());
}
private:
static void offsetsFromSizes(std::unordered_map<QnType, SizeType>& offsets,
......@@ -684,7 +683,7 @@ private:
TargetQuantumElectronsType targetQuantum_;
InputValidatorType_& ioIn_;
AtomKindBaseType* atomKind_;
mutable SuperOpHelperType* superOpHelper_;
mutable SuperOpHelperBaseType* superOpHelper_;
static LabeledOperatorsType labeledOperators_;
static ModelLinksType modelLinks_;
static VectorQnType qns_;
......
......@@ -114,7 +114,8 @@ public:
typedef typename ModelHelperType::BasisWithOperatorsType BasisWithOperatorsType;
typedef LabeledOperators<OperatorType> LabeledOperatorsType;
typedef ModelLinks<LabeledOperatorsType, SuperGeometryType> ModelLinksType;
typedef HamiltonianConnection<ModelLinksType, ModelHelperType> HamiltonianConnectionType;
typedef HamiltonianConnection<ModelLinksType, ModelHelperType, ParametersType>
HamiltonianConnectionType;
typedef typename HamiltonianConnectionType::VectorLinkType VectorLinkType;
typedef typename ModelHelperType::LeftRightSuperType LeftRightSuperType;
typedef typename PsimagLite::Vector<OperatorType>::Type VectorOperatorType;
......
......@@ -39,17 +39,17 @@ public:
if (taskNumber == 0) {
hc_.modelHelper().hamiltonianLeftProduct(xtemp_[threadNum], y_, aux_);
const SparseMatrixType& hamiltonian = hc_.modelHelper().leftRightSuper().
left().hamiltonian().getCRS();
hc_.kroneckerDumper().push(true, hamiltonian, y_);
// const SparseMatrixType& hamiltonian = hc_.modelHelper().leftRightSuper().
// left().hamiltonian().getCRS();
// hc_.kroneckerDumper().push(true, hamiltonian, y_);
return;
}
if (taskNumber == 1) {
hc_.modelHelper().hamiltonianRightProduct(xtemp_[threadNum], y_, aux_);
const SparseMatrixType& hamiltonian = hc_.modelHelper().leftRightSuper().
right().hamiltonian().getCRS();
hc_.kroneckerDumper().push(false, hamiltonian, y_);
// const SparseMatrixType& hamiltonian = hc_.modelHelper().leftRightSuper().
// right().hamiltonian().getCRS();
// hc_.kroneckerDumper().push(false, hamiltonian, y_);
return;
}
......@@ -66,11 +66,11 @@ public:
link2,
aux_);
hc_.kroneckerDumper().push(A->getCRS(),
B->getCRS(),
link2.value,
link2.fermionOrBoson,
y_);
// hc_.kroneckerDumper().push(A->getCRS(),
// B->getCRS(),
// link2.value,
// link2.fermionOrBoson,
// y_);
}
SizeType tasks() const { return hc_.tasks() + 2; }
......
......@@ -134,13 +134,12 @@ private:
MatrixComplexOrRealType& V,
SizeType i0)
{
SizeType p = lrs_.super().findPartitionNumber(phi.offset(i0));
const SizeType p = lrs_.super().findPartitionNumber(phi.offset(i0));
typename ModelHelperType::Aux aux(p, lrs_);
typename ModelType::HamiltonianConnectionType hc(lrs_,
model_.superGeometry(),
ModelType::modelLinks(),
currentTime_,
0);
model_.superOpHelper());
typename LanczosSolverType::MatrixType lanczosHelper(model_, hc, aux);
typename LanczosSolverType::ParametersSolverType params(io_,"Tridiag");
......
......@@ -5,6 +5,7 @@
namespace Dmrg {
template<typename SuperGeometryType, typename ParametersType>
class SuperOpHelperBase {
public:
......@@ -12,13 +13,17 @@ public:
typedef PsimagLite::Vector<SizeType>::Type VectorSizeType;
typedef std::pair<bool, SizeType> PairBoolSizeType;
SuperOpHelperBase(const VectorSizeType&,
const VectorSizeType&,
ProgramGlobals::DirectionEnum dir)
: dir_(dir) {}
SuperOpHelperBase(const SuperGeometryType& superGeometry)
: superGeometry_(superGeometry)
{}
virtual ~SuperOpHelperBase() {}
virtual void setToProduct(SizeType, SizeType, ProgramGlobals::DirectionEnum dir)
{
dir_ = dir;
}
// This below is for a plaquette, and will have to be
// written somewhere else
// testing devel FIXME TODO
......@@ -34,10 +39,23 @@ public:
return PairBoolSizeType(false, 0);
}
virtual SizeType leftIndex(VectorSizeType&, SizeType) const
{
throw PsimagLite::RuntimeError("SuperOpHelperBase::leftIndex\n");
}
// non virtual below
const SuperGeometryType& superGeometry() const
{
return superGeometry_;
}
const ProgramGlobals::DirectionEnum& dir() const { return dir_; }
private:
const SuperGeometryType& superGeometry_;
ProgramGlobals::DirectionEnum dir_;
};
}
......
......@@ -282,13 +282,12 @@ private:
const VectorType& sv,
SizeType p)
{
RealType fakeTime = 0;
const RealType fakeTime = 0;
typename ModelHelperType::Aux aux(p, BaseType::lrs());
typename ModelType::HamiltonianConnectionType hc(BaseType::lrs(),
BaseType::model().superGeometry(),
ModelType::modelLinks(),
fakeTime,
0);
BaseType::model().superOpHelper());
typename LanczosSolverType::MatrixType h(BaseType::model(), hc, aux);
paramsForSolver_.lotaMemory = true;
LanczosSolverType lanczosSolver(h,paramsForSolver_);
......
......@@ -815,10 +815,9 @@ private:
SizeType p = this->lrs().super().findPartitionNumber(phi.offset(i0));
typename ModelHelperType::Aux aux(p, BaseType::lrs());
typename ModelType::HamiltonianConnectionType hc(BaseType::lrs(),
BaseType::model().superGeometry(),
BaseType::ModelType::modelLinks(),
this->common().aoe().time(),
0);
model_.superOpHelper());
typename LanczosSolverType::MatrixType lanczosHelper(BaseType::model(),
hc,
aux);
......
......@@ -313,18 +313,17 @@ private:
SizeType whatTarget,
SizeType i0) const
{
SizeType p = this->lrs().super().findPartitionNumber(phi.offset(i0));
const SizeType p = this->lrs().super().findPartitionNumber(phi.offset(i0));
typename ModelHelperType::Aux aux(p, BaseType::lrs());
typename ModelType::HamiltonianConnectionType hc(BaseType::lrs(),
BaseType::model().superGeometry(),
ModelType::modelLinks(),
this->common().aoe().time(),
0);
BaseType::model().superOpHelper());
typename LanczosSolverType::MatrixType lanczosHelper(BaseType::model(),
hc,
aux);
SizeType total = phi.effectiveSize(i0);
const SizeType total = phi.effectiveSize(i0);
TargetVectorType phi2(total);
phi.extract(phi2,i0);
TargetVectorType x(total);
......
......@@ -284,10 +284,9 @@ private:
SizeType p = lrs_.super().findPartitionNumber(phi.offset(i0));
typename ModelHelperType::Aux aux(p, lrs_);
typename ModelType::HamiltonianConnectionType hc(lrs_,
model_.superGeometry(),
ModelType::modelLinks(),
time(),
0);
model_.superOpHelper());
MatrixLanczosType lanczosHelper(model_, hc, aux);
ProgramGlobals::VerboseEnum verbose = (model_.params().options.isSet("VerboseCheby"))
......@@ -380,10 +379,9 @@ private:
const SizeType p = lrs_.super().findPartitionNumber(phi.offset(i0));
typename ModelHelperType::Aux aux(p, lrs_);
typename ModelType::HamiltonianConnectionType hc(lrs_,
model_.superGeometry(),
ModelType::modelLinks(),
time(),
0);
model_.superOpHelper());
MatrixLanczosType lanczosHelper(model_, hc, aux);
ProgramGlobals::VerboseEnum verbose = (model_.params().options.isSet("VerboseCheby"))
? ProgramGlobals::VerboseEnum::YES
......
......@@ -188,7 +188,7 @@ private:
timeDirection_(timeDirection),
p_(lrs.super().findPartitionNumber(phi.offset(i0))),
aux_(p_, lrs),
hc_(lrs, model.superGeometry(), ModelType::modelLinks(), currentTime, 0),
hc_(lrs, ModelType::modelLinks(), currentTime, model.superOpHelper()),
lanczosHelper_(model, hc_, aux_)
{}
......