scnc
/
qmcpack
mirror of https://github.com/QMCPACK/qmcpack.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Fixing introduced bug on PPconvert 

git-svn-id: https://subversion.assembla.com/svn/qmcdev/trunk@6814 e5b18d87-469d-4833-9cc0-8cdfa06e9491
This commit is contained in:
Anouar Benali 2016-03-05 02:04:04 +00:00
parent d18be3b579
commit 533ea838b2
5 changed files with 786 additions and 67 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
src/QMCTools/CMakeLists.txt
View File

@ -21,6 +21,7 @@ SET(MOSRCS
BMakeFunc.cpp
GamesFMOParser.cpp
VSVBParser.cpp
QPParser.cpp
)
# create libmocommon

94
src/QMCTools/GamesFMOParser.cpp
View File

@ -20,13 +20,12 @@ GamesFMOParser::GamesFMOParser()
basisName = "Gaussian-G2";
Normalized = "no";
Mono=false;
BohrUnit=false;
BohrUnit=true;
MOtype="Canonical";
angular_type="cartesian";
readtype=0;
psi_tag="psi0";
ion_tag="ion0";
SpinMultiplicity = 1;
SpinRestricted=true;
TotNumAtom=0;
NumMonomer=0;
@ -35,7 +34,7 @@ GamesFMOParser::GamesFMOParser()
FMOMethod=0;
FMO=true;
MonomerID=0;
FixValence=true;
}
GamesFMOParser::GamesFMOParser(int argc, char** argv):
@ -44,13 +43,12 @@ GamesFMOParser::GamesFMOParser(int argc, char** argv):
basisName = "Gaussian-G2";
Normalized = "no";
Mono=false;
BohrUnit=false;
BohrUnit=true;
MOtype="Canonical";
angular_type="cartesian";
readtype=0;
psi_tag="psi0";
ion_tag="ion0";
SpinMultiplicity = 1;
SpinRestricted=true;
TotNumAtom=0;
NumMonomer=0;
@ -59,6 +57,7 @@ GamesFMOParser::GamesFMOParser(int argc, char** argv):
FMOMethod=0;
FMO=true;
MonomerID=0;
FixValence=true;
}
@ -79,6 +78,10 @@ void GamesFMOParser::parse(const std::string& fname)
fin.seekg(pivot_begin);
search(fin,"NUMBER OF CARTESIAN GAUSSIAN BASIS FUNCTIONS",aline);
parsewords(aline.c_str(),currentWords);
SizeOfBasisSet = atoi(currentWords[6].c_str());
search(fin,"TOTAL NUMBER OF ATOMS",aline);
parsewords(aline.c_str(),currentWords);
TotNumAtom = atoi(currentWords[4].c_str());
@ -110,6 +113,7 @@ void GamesFMOParser::parse(const std::string& fname)
readtype=0;
TotNumMonomer=NumMonomer;
IDMonomer = new string[NumMonomer];
@ -132,6 +136,7 @@ void GamesFMOParser::parse(const std::string& fname)
pivot= fin.tellg();
for(int i=0;i<NumMonomer;i++)
{
fin.seekg(pivot);
ostringstream convert;
convert <<i+1;
IDMonomer[i]="CURRENT N-MER COORDINATES, I= "+convert.str()+" J= 0 K= 0";
@ -151,12 +156,15 @@ void GamesFMOParser::parse(const std::string& fname)
IDDimer[count].IndexI=i-1;
IDDimer[count].IndexJ=j-1;
IDDimer[count].IndexK=0;
cout<<IDDimer[count].MyId<<endl;
count++;
}
}
if(FMOMethod==3){
fin.seekg(pivot);
count=0;
IDTrimer = new IDmer[NumTrimer];
for(int i=3; i<=NumMonomer; i++){
@ -168,10 +176,13 @@ void GamesFMOParser::parse(const std::string& fname)
convertK <<k;
IDTrimer[count].MyId="CURRENT N-MER COORDINATES, I= "+convertI.str()+" J= "+convertJ.str()+" K= "+convertK.str();
cout<<IDTrimer[count].MyId<<endl;
count++;
IDTrimer[count].IndexI=i-1;
IDTrimer[count].IndexJ=j-1;
IDTrimer[count].IndexK=k-1;
count++;
}
}
}
@ -210,8 +221,11 @@ void GamesFMOParser::parse(const std::string& fname)
for (int i=0; i<This; i++)
{
streampos pivot;
std::ifstream fin(fname.c_str());
pivot_begin= fin.tellg();
cout<<"Working On "<<Mer<<"-"<<i<<endl;
string temp_tag,Frag_Tag;
fin.seekg(pivot_begin);
if (Mer==0){
@ -247,21 +261,22 @@ void GamesFMOParser::parse(const std::string& fname)
FMOIndexK=IDTrimer[i].IndexK;
}
}
if(lookFor(fin,temp_tag,aline))
{
search(fin,"MUMBER OF ATOMS IN FRAGMENT",aline);
cout<<"temp_tag="<<temp_tag<<endl;
search(fin,"NUMBER OF ATOMS IN FRAGMENT",aline);
parsewords(aline.c_str(),currentWords);
NumberOfAtoms= atoi(currentWords[5].c_str());
cout<<"Total Number of Atoms="<<NumberOfAtoms<<endl;
search(fin,"NUMBER OF CARTESIAN ATOMIC ORBITALS",aline);
parsewords(aline.c_str(),currentWords);
SizeOfBasisSet= atoi(currentWords[5].c_str());
cout<<"Number of Cartesian atomic orbitals="<<SizeOfBasisSet<<endl;
search(fin,"TOTAL NUMBER OF MOS IN VARIATION SPACE",aline);
parsewords(aline.c_str(),currentWords);
numMO= atoi(currentWords[7].c_str());
cout<<"Number of Molecular Orbitals="<<numMO<<endl;
search(fin,"SPIN MULTIPLICITY",aline);
parsewords(aline.c_str(),currentWords);
SpinMultiplicity= atoi(currentWords[2].c_str());
@ -279,12 +294,18 @@ void GamesFMOParser::parse(const std::string& fname)
NumberOfBeta = atoi(currentWords[5].c_str());
cout<<"Number of beta electrons: " <<NumberOfBeta <<endl;
IonSystem.setName(ion_tag);
IonSystem.create(NumberOfAtoms);
GroupName.resize(NumberOfAtoms);
getGeometry(fin);
search(fin,"TOTAL NUMBER OF MOS IN VARIATION SPACE",aline);
parsewords(aline.c_str(),currentWords);
numMO= atoi(currentWords[7].c_str());
cout<<"Number of Molecular Orbitals="<<numMO<<endl;
fin.seekg(pivot_begin);
@ -296,9 +317,11 @@ void GamesFMOParser::parse(const std::string& fname)
getMO(fin,temp_tag);
Fmodump(psi_tag, ion_tag,Frag_Tag);
IonSystem.clear();
GroupName.clear();
fin.close();
}
else
{
@ -322,11 +345,13 @@ void GamesFMOParser::parse(const std::string& fname)
delete [] ESPIonChargeIndex;
delete [] ESPValenceChargeIndex;
delete [] ESPAtomicNumberIndex;
delete [] ESPSystem;
OHMMS::Controller->finalize();
}
void GamesFMOParser::getAllMonomer(std::istream& is, int Index)
{
const double ang_to_bohr=1.0/0.529177e0;
std::string aline;
//atomic numbers
vector<int> atomic_number,core;
@ -341,7 +366,7 @@ void GamesFMOParser::getAllMonomer(std::istream& is, int Index)
{
getwords(currentWords,is);
MyAtoms= atoi(currentWords[5].c_str());
ostringstream temp;
std::string ESP_ion_tag;
temp<<Index;
@ -358,13 +383,13 @@ void GamesFMOParser::getAllMonomer(std::istream& is, int Index)
currentWords[2] == "X" )
{
getwords(currentWords,is);
while(currentWords[0] != "------------------------------------------------------------")
while(currentWords[0]!="------------------------------------------------------------")
{
natms++;
double z=atof(currentWords[1].c_str());
int zint = (int)z; // is this dangerous???
int zint = (int)z;
atomic_number.push_back(zint);
q.push_back(z); // if using ECPs, change below
q.push_back(z);
tags.push_back(currentWords[0]);
pos.push_back(atof(currentWords[2].c_str()));
pos.push_back(atof(currentWords[3].c_str()));
@ -384,7 +409,7 @@ void GamesFMOParser::getAllMonomer(std::istream& is, int Index)
cerr<<"Could not find atomic coordinates for all atoms. \n";
abort();
}
SpeciesSet& ESPspecies(ESPSystem[Index].getSpeciesSet());
for(int i=0, ii=0; i<MyAtoms; i++)
{
@ -402,12 +427,22 @@ void GamesFMOParser::getAllMonomer(std::istream& is, int Index)
ESPspecies(ESPAtomicNumberIndex[Index],speciesID)=1;
ESPspecies(ESPIonChargeIndex[Index],speciesID)=ESP[Index][i];
}
ESPSystem[Index].R *= ang_to_bohr;
}
else
{
cerr<<" Error Finding Monomer coordinates!!! Impossible to create ESP charges!!"<<endl;
abort();
}
}
void GamesFMOParser::getGeometry(std::istream& is)
{
//atomic numbers
const double ang_to_bohr=1.0/0.529177e0;
vector<int> atomic_number,core;
vector<double> q,pos;
int natms=0;
@ -418,18 +453,20 @@ void GamesFMOParser::getGeometry(std::istream& is)
currentWords[2] == "X" )
{
getwords(currentWords,is);
while(currentWords[0] != "------------------------------------------------------------")
while(currentWords[0]!="------------------------------------------------------------")
{
natms++;
double z=atof(currentWords[1].c_str());
int zint = (int)z; // is this dangerous???
atomic_number.push_back(zint);
q.push_back(z); // if using ECPs, change below
q.push_back(z);
tags.push_back(currentWords[0]);
pos.push_back(atof(currentWords[2].c_str()));
pos.push_back(atof(currentWords[3].c_str()));
pos.push_back(atof(currentWords[4].c_str()));
//cout<<"Atom number "<<natms<<" Words:"<<currentWords[0]<<" "<< currentWords[2]<<" "<< currentWords[3]<<" "<< currentWords[4]<<" "<<endl;
getwords(currentWords,is);
}
}
else
@ -438,6 +475,7 @@ void GamesFMOParser::getGeometry(std::istream& is)
abort();
}
cout<<"natms="<<natms<<" NumberOfAtoms="<<NumberOfAtoms<<endl;
if(natms != NumberOfAtoms)
{
cerr<<"Could not find atomic coordinates for all atoms. \n";
@ -456,6 +494,7 @@ void GamesFMOParser::getGeometry(std::istream& is)
species(IonChargeIndex,speciesID)=q[i];
}
IonSystem.R *= ang_to_bohr;
}
@ -465,6 +504,7 @@ void GamesFMOParser::getGaussianCenters(std::istream& is)
{
gBound.resize(NumberOfAtoms+1);
int ng,nx;
int OldShell=1;
string aline;
std::map<std::string,int> basisDataMap;
int nUniqAt=0;
@ -511,16 +551,20 @@ void GamesFMOParser::getGaussianCenters(std::istream& is)
found=true;
}
getwords(currentWords,is); // empty line
int currPos=-1;
while(true)
{
getwords(currentWords,is);
if(currentWords.empty()) continue;
if(currentWords[0] == "------------------------------------------------------------")
if(currentWords[0] == "TOTAL" && currentWords[1] == "NUMBER" &&
currentWords[2] == "OF" && currentWords[3] == "BASIS")
{
ng=atoi(currentWords.back().c_str());
break;
}
if(currentWords.size() == 1) //found Species
@ -570,6 +614,7 @@ void GamesFMOParser::getGaussianCenters(std::istream& is)
coef[currPos].push_back(atof(currentWords[4].c_str()));
ncoeffpershell[currPos][nshll[currPos]]++;
shID[currPos][nshll[currPos]] = currentWords[1];
if(gsMap[currentWords[1]] == 2)
{
cerr<<"Can't handle SP basis states yet. Fix later.\n";
@ -616,7 +661,7 @@ void GamesFMOParser::getGaussianCenters(std::istream& is)
gC0.push_back(coef[indx][k]);
}
gBound[NumberOfAtoms] = gtot;
ng=gShell.size();
// ng=gShell.size();
}
@ -630,7 +675,8 @@ void GamesFMOParser::getMO(std::istream& is,std::string temp_tag)
for(int i=0;i<NumberOfAtoms;i++)
getwords(currentWords,is);
getwords(currentWords,is); // ----------------------
getwords(currentWords,is); // empty line
getwords(currentWords,is); // TOTAL NUMBER OF MOS IN VARIATION SPACE
getwords(currentWords,is); // Empty Line

619
src/QMCTools/QPParser.cpp Normal file
View File

@ -0,0 +1,619 @@
#include "QMCTools/QPParser.h"
#include <fstream>
#include <iterator>
#include <algorithm>
#include <set>
#include <map>
using namespace std;
QPParser::QPParser()
{
basisName = "Gaussian-G2";
Normalized = "no";
usingECP=false;
BohrUnit=true;
MOtype="Canonical";
angular_type="cartesian";
readtype=0;
NFZC=0;
numAO=0;
FixValence=true;
}
QPParser::QPParser(int argc, char** argv):
QMCGaussianParserBase(argc,argv)
{
basisName = "Gaussian-G2";
Normalized = "no";
usingECP=false;
BohrUnit=true;
MOtype="Canonical";
angular_type="cartesian";
SpinRestricted=true;
readtype=0;
NFZC=0;
numAO=0;
FixValence=true;
}
void QPParser::parse(const std::string& fname)
{
std::ifstream fin(fname.c_str());
pivot_begin= fin.tellg();
std::string aline;
search(fin,"do_pseudo",aline);
parsewords(aline.c_str(),currentWords);
if(currentWords[1]=="True")
usingECP=true;
else
usingECP=false;
cout<<"usingECP: " <<(usingECP?("yes"):("no")) <<endl;
cout.flush();
search(fin,"multi_det",aline);
parsewords(aline.c_str(),currentWords);
if(currentWords[1]=="True")
multideterminant=true;
else
multideterminant=false;
cout<<"Multideterminant: " <<(multideterminant?("yes"):("no")) <<endl;
search(fin,"ao_num",aline);
parsewords(aline.c_str(),currentWords);
numAO = atoi(currentWords[1].c_str());
cout<<"NUMBER OF AOs: " <<numAO <<endl;
SizeOfBasisSet = numAO;
cout<<"Size of Basis Set: " <<SizeOfBasisSet <<endl;
search(fin,"mo_num",aline);
parsewords(aline.c_str(),currentWords);
numMO = atoi(currentWords[1].c_str());
cout<<"NUMBER OF MOs: " <<numMO <<endl;
search(fin,"elec_alpha_num",aline);
parsewords(aline.c_str(),currentWords);
NumberOfAlpha = atoi(currentWords[1].c_str());
cout<<"Number of alpha electrons: " <<NumberOfAlpha <<endl;
search(fin,"elec_beta_num",aline);
parsewords(aline.c_str(),currentWords);
NumberOfBeta = atoi(currentWords[1].c_str());
cout<<"Number of beta electrons: " <<NumberOfBeta <<endl;
search(fin,"elec_tot_num",aline);
parsewords(aline.c_str(),currentWords);
NumberOfEls = atoi(currentWords[1].c_str());
cout<<"Number of electrons: " <<NumberOfEls <<endl;
search(fin,"spin_multiplicity",aline);
parsewords(aline.c_str(),currentWords);
SpinMultiplicity = atoi(currentWords[1].c_str());
cout<<"SPIN MULTIPLICITY: " <<SpinMultiplicity <<endl;
SpinRestricted=true;
search(fin,"nucl_num",aline);
parsewords(aline.c_str(),currentWords);
NumberOfAtoms = atoi(currentWords[1].c_str());
cout<<"NUMBER OF ATOMS: " <<NumberOfAtoms <<endl;
IonSystem.create(NumberOfAtoms);
GroupName.resize(NumberOfAtoms);
getGeometry(fin);
fin.seekg(pivot_begin);
getGaussianCenters(fin);
fin.seekg(pivot_begin);
MOtype = "Canonical";
readtype=0;
getMO(fin);
fin.close();
if(multideterminant)
{
fin.open(fname.c_str());
QP=true;
pivot_begin= fin.tellg();
search(fin,"BEGIN_DET",aline);
getwords(currentWords,fin);//Empty Line
getwords(currentWords,fin);//mo_num
getwords(currentWords,fin);//Number Of determinants
cout<<"Found "<<currentWords[1]<<" Determinants"<<endl;
ci_size=atoi(currentWords[1].c_str());
NFZC = 0;
NAC = numMO;
NEXT =0;
NTOT=NEXT+NAC;
cout<<"# core, #active, #external: " <<NFZC <<" " <<NAC <<" " <<NEXT <<endl;
//fin.seekg(pivot_begin);
getQPCI(fin);
}
}
void QPParser::getGeometry(std::istream& is)
{
//atomic numbers
vector<int> atomic_number,core;
vector<double> q,pos;
int natms=0;
tags.clear();
is.seekg(pivot_begin);
//read atomic info
bool notfound=true;
do
{
if(is.eof())
{
cerr<<"Could not find atomic coordinates. \n";
abort();
}
getwords(currentWords,is);
if(currentWords.size() < 4 )
continue;
if (currentWords.size() == 4)
if(currentWords[0] == "Atomic" &&
currentWords[1] == "coord" &&
currentWords[2] == "in" &&
currentWords[3] == "Bohr" )
{
notfound=false;
getwords(currentWords,is);
while(currentWords.size() != 0)
{
if(currentWords[0] == "BEGIN_BASIS_SET")
break;
natms++;
double z=atof(currentWords[1].c_str());
int zint = (int)z; // is this dangerous???
atomic_number.push_back(zint);
q.push_back(z); // if using ECPs, change below
tags.push_back(currentWords[0]);
pos.push_back(atof(currentWords[2].c_str()));
pos.push_back(atof(currentWords[3].c_str()));
pos.push_back(atof(currentWords[4].c_str()));
getwords(currentWords,is);
}
}
}
while(notfound);
if(natms != NumberOfAtoms)
{
cerr<<"Could not find atomic coordinates for all atoms. \n";
abort();
}
//ECP PART!!!
if(usingECP==true){
is.seekg(pivot_begin);
notfound=true;
while(notfound)
{
if(is.eof())
{
cerr<<"Problem looking for ECPs, this should not happen. Contact developers for help. \n";
abort();
}
getwords(currentWords,is);
// this should appear below the ECP section in the output file
// so use this to avoid going all the way to the bottom
if(currentWords.size() == 0 )
continue;
if( currentWords[0] == "BEGIN_PSEUDO")
{
core.resize(NumberOfAtoms);
// getwords(currentWords,is); // -------------
bool done=false;
while(!done)
{
if(is.eof())
{
cerr<<"2 Found ECPs, but problem looking ZCORE data.\n";
abort();
}
getwords(currentWords,is);
if(currentWords.size() == 0)
continue;
if(currentWords.size() == 1)
{
if(currentWords[0] == "END_PSEUDO")
done=true;
}
if(currentWords[0] == "PARAMETERS" &&
currentWords[1] == "FOR" )
{
//done=true;
std::vector<std::string>::iterator it,it0;
it = find(currentWords.begin(),currentWords.end(),"ZCORE");
it0 = find(currentWords.begin(),currentWords.end(),"ATOM");
if(it0 == currentWords.end())
{
cerr<<"Problem with ECP data. Didn't found ATOM tag\n";
cerr<<is <<endl;
abort();
}
it0++;
int nq0 = atoi(it0->c_str())-1;
if(it != currentWords.end())
{
it++;
cout<<"Avant nq0="<<nq0<<" core[nq0]="<<core[nq0]<<" q[nq0]="<<q[nq0]<<endl;
core[nq0] = atoi(it->c_str());
q[nq0] -= core[nq0];
cout<<"Apres nq0="<<nq0<<" core[nq0]="<<core[nq0]<<" q[nq0]="<<q[nq0]<<endl;
cout<<"1 Found ECP for atom " <<nq0 <<" with zcore " <<core[nq0] <<endl;
}
else
{
it = find(currentWords.begin(),currentWords.end(),"ATOM");
if(it == currentWords.end())
{
cerr<<"Problem with ECP data. Didn't found ATOM tag\n";
cerr<<"Atom: " <<nq0 <<endl;
abort();
}
std::vector<std::string>::iterator it2=it;
it2++;
int nq = atoi(it2->c_str());
if(nq != nq0+1)
{
cerr<<"Problem with ECP data. ID's don't agree\n";
cerr<<"Atom: " <<nq0 <<endl;
abort();
}
it = find(it2,currentWords.end(),"ATOM");
if(it == currentWords.end())
{
cerr<<"Problem with ECP data (2).\n";
cerr<<"Atom: " <<nq0 <<endl;
abort();
}
nq = atoi((it+1)->c_str());
core[nq0] = core[nq-1];
q[nq0] -= core[nq0];
cout<<"2 Found ECP for atom " <<nq0 <<" with zcore " <<core[nq0] <<endl;
}
}
}
notfound=false;
}
else
{
if(currentWords.size() < 2 )
continue;
if( currentWords[0] == "END_PSEUDO")
break;
}
}
}
SpeciesSet& species(IonSystem.getSpeciesSet());
for(int i=0, ii=0; i<NumberOfAtoms; i++)
{
IonSystem.R[i][0]=pos[ii++];
IonSystem.R[i][1]=pos[ii++];
IonSystem.R[i][2]=pos[ii++];
GroupName[i]=IonName[atomic_number[i]];
int speciesID = species.addSpecies(GroupName[i]);
IonSystem.GroupID[i]=speciesID;
species(AtomicNumberIndex,speciesID)=atomic_number[i];
species(IonChargeIndex,speciesID)=q[i];
}
}
void QPParser::getGaussianCenters(std::istream& is)
{
string Shell_temp;
gBound.resize(NumberOfAtoms+1);
int ng,nx;
string aline;
std::map<std::string,int> basisDataMap;
int nUniqAt=0;
for(int i=0; i<NumberOfAtoms; i++)
{
std::map<std::string,int>::iterator it(basisDataMap.find(tags[i]));
if(it == basisDataMap.end())
{
basisDataMap[tags[i]]=nUniqAt++;
}
}
vector<vector<double> > expo(nUniqAt),coef(nUniqAt),coef2(nUniqAt);
vector<int> nshll(nUniqAt,0); //use this to
vector<vector<int> > ncoeffpershell(nUniqAt);
vector<vector<std::string> > shID(nUniqAt);
std::map<std::string,int> gsMap;
gsMap[std::string("S")]=1;
gsMap[std::string("SP")]=2;
gsMap[std::string("P")]=3;
gsMap[std::string("D")]=4;
gsMap[std::string("F")]=5;
gsMap[std::string("G")]=6;
gsMap[std::string("H")]=7;
gsMap[std::string("I")]=8;
is.seekg(pivot_begin);
bool found=false;
while(!found)
{
if(is.eof())
{
cerr<<"Problem with basis set data.\n";
abort();
}
getwords(currentWords,is);
if(currentWords.size() >2)
continue;
if(currentWords[0] == "BEGIN_BASIS_SET")
found=true;
}
is.seekg(pivot_begin);
int currPos=-1;
lookFor(is,"BEGIN_BASIS_SET");
int NbCoeffperShell=0;
getwords(currentWords,is);
bool allbase=true;
while(allbase==true)
{
if(currentWords.empty()) getwords(currentWords,is);
if(currentWords.size() ==1)
{
std::map<std::string,int>::iterator it(basisDataMap.find(currentWords[0]));
if(it == basisDataMap.end())
{
cerr<<"Error in parser.\n";
abort();
}
currPos=it->second;
nshll[currPos]=0;
ncoeffpershell[currPos].clear();
ncoeffpershell[currPos].push_back(0);
shID[currPos].clear();
shID[currPos].push_back("NONE");
bool group=true;
do{
getwords(currentWords,is);
if(currentWords[0]=="S"||currentWords[0]=="P"||currentWords[0]=="D"||currentWords[0]=="F"||currentWords[0]=="G"||currentWords[0]=="H"||currentWords[0]=="I") //empty line after the shell
{
shID[currPos].push_back(currentWords[0]);
Shell_temp=currentWords[0];
NbCoeffperShell=atoi(currentWords[1].c_str());
ncoeffpershell[currPos][nshll[currPos]]=NbCoeffperShell;
for(int nbcoef=0;nbcoef<NbCoeffperShell;nbcoef++){
getwords(currentWords,is);
expo[currPos].push_back(atof(currentWords[1].c_str()));
coef[currPos].push_back(atof(currentWords[2].c_str()));
shID[currPos][nshll[currPos]] = Shell_temp;
cout << currPos << ":" <<expo[currPos].back() << " " << coef[currPos].back() << " "
<< ncoeffpershell[currPos][nshll[currPos]]
<< " " << shID[currPos][nshll[currPos]]<<endl;
}
nshll[currPos]++;
ncoeffpershell[currPos].push_back(0);
shID[currPos].push_back("NONE");
}
else{
if(currentWords[0] == "END_BASIS_SET")
{
ng=SizeOfBasisSet;
group=false;
allbase=false;
break;
}
else
{
break;
}
}
}while (group==true);
}
else
{
cerr<<"error in parser"<<endl;
abort();
}
}
gShell.clear();
gNumber.clear();
gExp.clear();
gC0.clear();
gC1.clear();
int gtot=0;
for(int i=0; i<NumberOfAtoms; i++)
{
std::map<std::string,int>::iterator it(basisDataMap.find(tags[i]));
if(it == basisDataMap.end())
{
cerr<<"Error in parser.\n";
abort();
}
gBound[i] = gtot;
int indx = it->second;
gtot+=nshll[indx];
for(int k=0; k<nshll[indx]; k++){
gShell.push_back(gsMap[shID[indx][k]]);
}
for(int k=0; k<nshll[indx]; k++)
gNumber.push_back(ncoeffpershell[indx][k]);
for(int k=0; k<expo[indx].size(); k++)
gExp.push_back(expo[indx][k]);
for(int k=0; k<coef[indx].size(); k++)
gC0.push_back(coef[indx][k]);
}
gBound[NumberOfAtoms] = gtot;
}
void QPParser::getMO(std::istream& is)
{
EigVal_alpha.resize(numMO);
EigVal_beta.resize(numMO);
EigVec.resize(2*SizeOfBasisSet*numMO);
std::string aline;
search(is,"BEGIN_MO");
getwords(currentWords,is); // empty line
vector<double> dummy(50);
Matrix<double> CartMat(numMO,SizeOfBasisSet);
streampos pivot;
pivot= is.tellg();
std::vector<std::string> CartLabels(SizeOfBasisSet);
getwords(currentWords,is);
for(int k=0; k<SizeOfBasisSet; k++)
{
getwords(currentWords,is);
CartLabels[k] = currentWords[0];
}
is.seekg(pivot);
getMO_single_set(is, CartMat, EigVal_alpha);
int cnt=0;
for(int i=0; i<numMO; i++)
for(int k=0; k<SizeOfBasisSet; k++)
EigVec[cnt++] = CartMat[i][k];
// beta states for now
if(!SpinRestricted)
{
search(is,"BEGIN_MO");
getwords(currentWords,is); // empty line
getMO_single_set(is, CartMat, EigVal_beta);
}
for(int i=0; i<numMO; i++)
for(int k=0; k<SizeOfBasisSet; k++)
EigVec[cnt++] = CartMat[i][k];
cout<<"Finished reading MO." <<endl;
}
void QPParser::getMO_single_set(std::istream& is, Matrix<double> &CartMat, std::vector<value_type>& EigVal)
{
int nq = numMO/4;
int rem = numMO%4;
int cnt=0;
for(int i=0; i<nq; i++)
{
getwords(currentWords,is);
for(int k=0; k<SizeOfBasisSet; k++)
{
getwords(currentWords,is);
if(currentWords.size() == 6)
{
CartMat[cnt][k] = atof(currentWords[2].c_str()) ;
CartMat[cnt+1][k] = atof(currentWords[3].c_str()) ;
CartMat[cnt+2][k] = atof(currentWords[4].c_str()) ;
CartMat[cnt+3][k] = atof(currentWords[5].c_str()) ;
}
else
{
cerr<<"Problem reading orbitals!!" <<endl;
abort();
}
}
getwords(currentWords,is);
cnt+=4;
}
if(rem > 0)
{
getwords(currentWords,is);
for(int k=0; k<SizeOfBasisSet; k++)
{
getwords(currentWords,is);
for(int i=0; i<rem; i++)
{
CartMat[cnt+i][k] = atof(currentWords[2+i].c_str()) ;
}
}
getwords(currentWords,is);
}
}
void QPParser::getQPCI(std::istream& is)
{
int Ci_indx;
double Ci_Coeff;
CIcoeff.clear();
CIalpha.clear();
CIbeta.clear();
CIcoeff.resize(ci_size);
CIalpha.resize(ci_size);
CIbeta.resize(ci_size);
for (int i=0; i<ci_size;i++)
{
getwords(currentWords,is); //Empty Line
if(currentWords[0] == "END_DET")
{
cout<<"Done reading determinants"<<endl;
break;
}
getwords(currentWords,is); //Coeff
CIcoeff[i]=atof(currentWords[0].c_str());
//Alpha spin
getwords(currentWords,is);
CIalpha[i]=currentWords[0];
//Beta spin
getwords(currentWords,is);
CIbeta[i]=currentWords[0];
}
ci_nea=ci_neb=0;
for(int i=0; i<CIalpha[0].size(); i++)
if(CIalpha[0].at(i) == '1')
ci_nea++;
for(int i=0; i<CIbeta[0].size(); i++)
if(CIbeta[0].at(i) == '1')
ci_neb++;
if(CIalpha[0].size() != CIbeta[0].size())
{
cerr<<"QMCPack can't handle different number of active orbitals in alpha and beta channels right now. Contact developers for help (Miguel).\n";
abort();
}
int ds=SpinMultiplicity-1;
int neb= (NumberOfEls-ds)/2;
int nea= NumberOfEls-NumberOfBeta;
ci_nca = nea-ci_nea;
ci_ncb = neb-ci_neb;
cout<<" Done reading CIs!!"<<endl;
ci_nstates = CIalpha[0].size();
}

42
src/QMCTools/QPParser.h Normal file
View File

@ -0,0 +1,42 @@
//#ifndef QMCPLUSPLUS_TOOLS_GAMESS_OUT_H
//#define QMCPLUSPLUS_TOOLS_GAMESS_OUT_H
#include "QMCTools/QMCGaussianParserBase.h"
#include <iostream>
#include <sstream>
#include <iomanip>
#include <vector>
#include "OhmmsPETE/TinyVector.h"
#include "OhmmsData/OhmmsElementBase.h"
class QPParser: public QMCGaussianParserBase,
public OhmmsAsciiParser
{
public:
QPParser();
QPParser(int argc, char** argv);
streampos pivot_begin;
vector<std::string> tags;
bool usingECP;
std::string MOtype;
int readtype, numAO;
int NFZC, NEXT, NTOT, NAC;
void parse(const std::string& fname);
void getGeometry(std::istream& is);
void getGaussianCenters(std::istream& is);
void getMO(std::istream& is);
void getMO_single_set(std::istream& is, Matrix<double> &CartMat, std::vector<value_type>& EigVal_alpha);
void getQPCI(std::istream& is);
};
//#endif

97
src/QMCTools/convert4qmc.cpp
View File

@ -2,6 +2,9 @@
#include "QMCTools/GaussianFCHKParser.h"
#include "QMCTools/GamesXmlParser.h"
#include "QMCTools/GamesAsciiParser.h"
#include "QMCTools/VSVBParser.h"
#include "QMCTools/QPParser.h"
#include "QMCTools/GamesFMOParser.h"
#include "QMCTools/BParser.h"
#include "Message/Communicate.h"
#include "Utilities/OhmmsInfo.h"
@ -13,8 +16,8 @@ int main(int argc, char **argv)
{
if(argc<2)
{
std::cout << "Usage: convert [-gaussian|-casino|-gamesxml|-gamessAscii] filename ";
std::cout << "[-nojastrow -hdf5 -psi_tag psi0 -ion_tag ion0 -gridtype log|log0|linear -first ri -last rf -size npts -ci file.out -threshold cimin -NaturalOrbitals NumToRead -add3BodyJ -prefix title]"
std::cout << "Usage: convert [-gaussian|-casino|-gamesxml|-gamessAscii|-gamessFMO |-VSVB| -QP] filename ";
std::cout << "[-nojastrow -hdf5 -psi_tag psi0 -ion_tag ion0 -gridtype log|log0|linear -first ri -last rf -size npts -ci file.out -threshold cimin -NaturalOrbitals NumToRead -add3BodyJ -prefix title -addCusp]"
<< std::endl;
std::cout << "Defaults : -gridtype log -first 1e-6 -last 100 -size 1001 -ci required -threshold 0.01 -prefix sample" << std::endl;
std::cout << "When the input format is missing, the extension of filename is used to determine the parser " << std::endl;
@ -36,7 +39,7 @@ int main(int argc, char **argv)
string ion_tag("ion0");
string prefix("sample");
bool usehdf5=false;
bool ci=false,zeroCI=false,orderByExcitation=false;
bool ci=false,zeroCI=false,orderByExcitation=false,VSVB=false, fmo=false,addCusp=false;
double thres=0.01;
int readNO=0; // if > 0, read Natural Orbitals from gamess output
int readGuess=0; // if > 0, read Initial Guess from gamess output
@ -58,6 +61,23 @@ int main(int argc, char **argv)
parser = new GamesAsciiParser(argc,argv);
in_file =argv[++iargc];
}
else if(a == "-QP")
{
parser = new QPParser(argc,argv);
in_file =argv[++iargc];
}
else if(a == "-VSVB")
{
parser = new VSVBParser(argc,argv);
in_file =argv[++iargc];
VSVB=true;
}
else if(a == "-gamessFMO")
{
parser = new GamesFMOParser(argc,argv);
in_file =argv[++iargc];
fmo=true;
}
else if(a == "-casino")
{
parser = new CasinoParser(argc,argv);
@ -89,6 +109,10 @@ int main(int argc, char **argv)
ci=true;
punch_file = argv[++iargc];
}
else if(a == "-addCusp" )
{
addCusp = true;
}
else if(a == "-threshold" )
{
thres = atof(argv[++iargc]);
@ -155,46 +179,33 @@ int main(int argc, char **argv)
exit(1);
}
}
parser->Title=prefix;
parser->UseHDF5=usehdf5;
parser->IonSystem.setName(ion_tag);
parser->multideterminant=ci;
parser->ci_threshold=thres;
parser->readNO=readNO;
parser->orderByExcitation=orderByExcitation;
parser->zeroCI = zeroCI;
parser->readGuess=readGuess;
parser->outputFile=punch_file;
parser->parse(in_file);
parser->dump(psi_tag, ion_tag);
OHMMS::Controller->finalize();
if(fmo)
{
parser->Title=prefix;
parser->UseHDF5=usehdf5;
parser->DoCusp=addCusp;
parser->parse(in_file);
}
else
{
parser->Title=prefix;
parser->DoCusp=addCusp;
parser->UseHDF5=usehdf5;
parser->IonSystem.setName(ion_tag);
parser->multideterminant=ci;
parser->ci_threshold=thres;
parser->readNO=readNO;
parser->orderByExcitation=orderByExcitation;
parser->zeroCI = zeroCI;
parser->readGuess=readGuess;
parser->outputFile=punch_file;
parser->VSVB=VSVB;
parser->parse(in_file);
parser->dump(psi_tag, ion_tag);
OHMMS::Controller->finalize();
}
return 0;
}
/*
int main(int argc, char **argv) {
char buffer[200];
std::string _txt;
std::string _data;
double _t;
std::cout.setf(std::ios::scientific, std::ios::floatfield);
std::cout.setf(std::ios::right,std::ios::adjustfield);
std::cout.precision(12);
while(std::cin.getline( buffer, sizeof ( buffer ) ) ){
std::istringstream stream(buffer);
if(isdigit(buffer[1])) {
while(stream>>_t) {
std::cout << std::setw(21) << _t ;
}
std::cout << std::endl;
} else {
if(stream>>_t) {
std::cout << "probably numbers " << buffer << std::endl;
} else {
std::cout << "probably statement " << buffer << std::endl;
}
}
}
}
*/