Commit 99854983 by Francois Gygi

Fixed loading of gfdata into slater_det with non-zero deltaspin


git-svn-id: http://qboxcode.org/svn/qb/trunk@1316 cba15fb0-1239-40c8-b417-11db7ca47a34
parent 3db3b8a0
......@@ -28,7 +28,7 @@ using namespace std;
////////////////////////////////////////////////////////////////////////////////
SampleHandler::SampleHandler(Sample& s, DoubleMatrix& gfdata) :
s_(s), gfdata_(gfdata) {}
s_(s), gfdata_(gfdata), current_gfdata_pos(0) {}
////////////////////////////////////////////////////////////////////////////////
SampleHandler::~SampleHandler() {}
......@@ -67,14 +67,14 @@ StructureHandler* SampleHandler::startSubHandler(const XMLCh* const uri,
else if ( qnm == "wavefunction" )
{
read_wf = true;
return new WavefunctionHandler(s_.wf,gfdata_);
return new WavefunctionHandler(s_.wf,gfdata_,current_gfdata_pos);
}
else if ( qnm == "wavefunction_velocity" )
{
read_wfv = true;
assert(read_wf);
s_.wfv = new Wavefunction(s_.wf);
return new WavefunctionHandler(*s_.wfv,gfdata_);
return new WavefunctionHandler(*s_.wfv,gfdata_,current_gfdata_pos);
}
else
{
......
......@@ -33,6 +33,7 @@ class SampleHandler : public StructureHandler
Sample& s_;
DoubleMatrix& gfdata_;
int current_gfdata_pos;
public:
......
......@@ -35,7 +35,7 @@ using namespace std;
////////////////////////////////////////////////////////////////////////////////
WavefunctionHandler::WavefunctionHandler(Wavefunction& wf,
DoubleMatrix& gfdata) : wf_(wf), gfdata_(gfdata) {}
DoubleMatrix& gfdata, int& current_gfdata_pos ) : wf_(wf), gfdata_(gfdata) {}
////////////////////////////////////////////////////////////////////////////////
WavefunctionHandler::~WavefunctionHandler() {}
......@@ -335,7 +335,8 @@ void WavefunctionHandler::startElement(const XMLCh* const uri,
void WavefunctionHandler::endElement(const XMLCh* const uri,
const XMLCh* const localname, const XMLCh* const qname, string& content)
{
bool onpe0 = wf_.context().onpe0();
const Context& ctxt = wf_.context();
bool onpe0 = ctxt.onpe0();
string locname(XMLString::transcode(localname));
//cout << " WavefunctionHandler::endElement " << locname << endl;
if ( locname == "density_matrix")
......@@ -361,19 +362,19 @@ void WavefunctionHandler::endElement(const XMLCh* const uri,
SlaterDet* sd = wf_.sd(current_ispin,current_ikp);
#if DEBUG
cout << ctxt_.mype() << ": mapping gfdata matrix..."
cout << ctxt.mype() << ": mapping gfdata matrix..."
<< endl;
cout << ctxt_.mype() << ": gfdata: (" << gfdata.m() << "x" << gfdata.n()
<< ") (" << gfdata.mb() << "x" << gfdata.nb() << ") blocks" << endl;
cout << ctxt_.mype() << ": gfdata.mloc()=" << gfdata.mloc()
<< " gfdata.nloc()=" << gfdata.nloc() << endl;
cout << ctxt.mype() << ": gfdata: (" << gfdata_.m() << "x" << gfdata_.n()
<< ") (" << gfdata_.mb() << "x" << gfdata_.nb() << ") blocks" << endl;
cout << ctxt.mype() << ": gfdata.mloc()=" << gfdata_.mloc()
<< " gfdata.nloc()=" << gfdata_.nloc() << endl;
#endif
sd->set_occ(dmat_);
const Basis& basis = sd->basis();
#if DEBUG
cout << ctxt_.mype() << ": ft->np012loc()=" << ft->np012loc() << endl;
cout << ctxt_.mype() << ": ft->context(): " << ft->context();
cout << ctxt.mype() << ": ft->np012loc()=" << ft->np012loc() << endl;
cout << ctxt.mype() << ": ft->context(): " << ft->context();
#endif
ComplexMatrix& c = sd->c();
......@@ -391,15 +392,26 @@ void WavefunctionHandler::endElement(const XMLCh* const uri,
wftmpr_block_size = 2*ft->np012loc(0);
}
#if DEBUG
cout << ctxt_.mype() << ": wftmpr_size: " << wftmpr_size << endl;
cout << ctxt_.mype() << ": wftmpr_block_size: "
cout << ctxt.mype() << ": wftmpr_size: " << wftmpr_size << endl;
cout << ctxt.mype() << ": wftmpr_block_size: "
<< wftmpr_block_size << endl;
#endif
DoubleMatrix wftmpr(sd->context(),wftmpr_size,sd->nst(),
wftmpr_block_size,c.nb());
wftmpr.getsub(gfdata_,wftmpr_size,sd->nst(),0,
current_ikp*sd->nst()+current_ispin*wf_.nkp()*sd->nst());
#if DEBUG
// parameters of the getsub operation
cout << "WavefunctionHandler: current_ikp= " << current_ikp << endl;
cout << "WavefunctionHandler: current_ispin=" << current_ispin << endl;
cout << "WavefunctionHandler: sd->nst()= " << sd->nst() << endl;
cout << "WavefunctionHandler: wf_.nkp()= " << wf_.nkp() << endl;
cout << "WavefunctionHandler: current_gfdata_pos= "
<< current_gfdata_pos << endl;
#endif
assert(current_gfdata_pos < gfdata_.n());
wftmpr.getsub(gfdata_,wftmpr_size,sd->nst(),0,current_gfdata_pos);
current_gfdata_pos += sd->nst();
#if DEBUG
// Check orthogonality by computing overlap matrix
......@@ -451,6 +463,7 @@ void WavefunctionHandler::endElement(const XMLCh* const uri,
for ( int ispin = 0; ispin < wf_.nspin(); ispin++ )
wf_.nst_[ispin] = wf_.sd_[ispin][0]->nst();
#if 1
// if nspin==2, adjust deltaspin_ to reflect the number of states
// of each spin that were read
if ( wf_.nspin() == 2 )
......@@ -461,6 +474,7 @@ void WavefunctionHandler::endElement(const XMLCh* const uri,
// even and odd number of electrons
wf_.deltaspin_ = (wf_.nst(0)-wf_.nst(1))/2;
}
#endif
}
}
......
......@@ -38,7 +38,7 @@ class WavefunctionHandler : public StructureHandler
int nx_, ny_, nz_;
int current_gf_nx,current_gf_ny,current_gf_nz;
std::string current_gf_encoding;
int current_ispin,current_ikp;
int current_ispin,current_ikp, current_gfdata_pos;
std::vector<double> dmat_;
double current_kx, current_ky, current_kz, current_weight;
int current_size;
......@@ -68,7 +68,8 @@ class WavefunctionHandler : public StructureHandler
const XMLCh* const localname, const XMLCh* const qname,
const StructureHandler* const subHandler);
WavefunctionHandler(Wavefunction& wf, DoubleMatrix& gfdata);
WavefunctionHandler(Wavefunction& wf, DoubleMatrix& gfdata,
int& current_gfdata_pos);
~WavefunctionHandler();
};
#endif
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment