EshdfFile Class Reference

Collaboration diagram for EshdfFile:

Public Member Functions
	EshdfFile (const std::string &hdfFileName)
	~EshdfFile ()
void	writeQEBoilerPlate (const XmlNode &qeXml)
void	writeQboxBoilerPlate (const XmlNode &qboxSample)
void	writeQESupercell (const XmlNode &qeXml)
void	writeQboxSupercell (const XmlNode &qboxSample)
void	writeQEAtoms (const XmlNode &qeXml)
void	writeQboxAtoms (const XmlNode &qboxSample)
void	writeQEElectrons (const XmlNode &qeXml, const std::string &dir_name)
void	writeQboxElectrons (const XmlNode &qboxSample)

Private Types
using	momap_t = std::map< std::vector< int >, std::complex< double > >
using	mopair_t = std::pair< std::vector< int >, std::complex< double > >

Private Member Functions
int	wrapped (int i, int size) const
int	getIntsOnly (const std::string &str) const
void	writeApplication (const std::string &appName, int major, int minor, int sub)
void	writeCreator ()
void	writeFormat ()
void	writeVersion ()
void	readInEigFcn (const XmlNode &nd, FftContainer &cont)
void	handleSpinGroup (const XmlNode *nd, double &nocc, FftContainer &cont)
double	getOccupation (const XmlNode *nd) const
void	handleDensity (const XmlNode &qeXml, const std::string &dir_name, int spinpol)
std::vector< double >	getPtvs (const XmlNode &qeXml)
void	processKPts (const XmlNode &band_structure_xml, const std::vector< double > &ptvs, std::vector< std::vector< double >> &eigenvals, std::vector< std::vector< double >> &occupations, std::vector< KPoint > &kpts, std::vector< double > &weights, std::vector< int > &ngvecs)
void	getNumElectrons (std::vector< std::vector< double >> &occupations, std::vector< double > &weights, int &nup, int &ndn, int spinpol, int ncol)
void	handleKpt (int kpt_num, const std::string &dir_name, KPoint &kpt, const std::vector< double > &eigenvalues, double weight, int spinpol, int noncol, const momap_t &moref)
void	readKptGvecs (int kpt_num, const std::string &dir_name, int spinpol, momap_t &morefmap)
qmcplusplus::hdf_archive	openHdfFileForRead (const std::string &fname)
	EshdfFile (const EshdfFile &f)
EshdfFile &	operator= (const EshdfFile &f)

Private Attributes
qmcplusplus::hdf_archive	outfile_

Detailed Description

helper class meant to control access ot the eshdf file that is being created and provide helper functions so that the proper xml data from QBOX or quantum espresso can be handed to it and the data can be processed and written to the eshdf file

Definition at line 32 of file WriteEshdf.h.

Member Typedef Documentation

momap_t

using momap_t = std::map<std::vector<int>, std::complex<double> >

private

Definition at line 34 of file WriteEshdf.h.

mopair_t

using mopair_t = std::pair<std::vector<int>, std::complex<double> >

private

Definition at line 35 of file WriteEshdf.h.

Constructor & Destructor Documentation

EshdfFile() [1/2]

EshdfFile ( const EshdfFile &  f )

private

EshdfFile() [2/2]

EshdfFile

(

const std::string &

hdfFileName

)

constructor that sets up the hdf_archive data structure that will hold the eshdf file to be created. Needs to be give the name for that file (hdfFileName)

Definition at line 79 of file WriteEshdf.cpp.

{ outfile_.create(hdfFileName); }

~EshdfFile()

~EshdfFile

(

)

Definition at line 81 of file WriteEshdf.cpp.

{ outfile_.close(); }

Member Function Documentation

getIntsOnly()

int getIntsOnly ( const std::string &  str ) const

private

helper function that takes a string and returns an int consisting of only the numerals in the string

Definition at line 261 of file WriteEshdf.cpp.

{
  stringstream numerals;
  int temp;
  for (int i = 0; i < str.size(); ++i)
  {
    string blah = str.substr(i, 1);
    if (stringstream(blah) >> temp)
    {
      numerals << str[i];
    }
  }
  int result;
  numerals >> result;
  return result;
}

getNumElectrons()

void getNumElectrons ( std::vector< std::vector< double >> &  occupations, std::vector< double > &  weights, int &  nup, int &  ndn, int  spinpol, int  ncol  )

private

will use the occupations and weights to figure out the total number of up and down electrons (nup,ndn). Also needs a flag to see if the calculation was spin_polarized (0 for no, 1 for yes) and another to see if noncolinear spins were used (0 for no, 1 for yes)

Definition at line 608 of file WriteEshdf.cpp.

{
  nup = 0;
  ndn = 0;

  double nup_flt = 0.0;
  double ndn_flt = 0.0;

  for (int i = 0; i < weights.size(); i++)
  {
    if (noncol == 1)
    {
      double sum = 0.0;
      for (int j = 0; j < occupations[i].size(); j++)
      {
        sum += occupations[i][j];
      }
      nup_flt += sum * weights[i];
      ndn_flt = -1.0;
    }
    else if (spinpol == 0)
    {
      double sum = 0.0;
      for (int j = 0; j < occupations[i].size(); j++)
      {
        sum += occupations[i][j];
      }
      nup_flt += sum * weights[i] / 2.0;
      ndn_flt += sum * weights[i] / 2.0;
    }
    else
    {
      double sum_up = 0.0;
      double sum_dn = 0.0;
      for (int j = 0; j < occupations[i].size() / 2; j++)
      {
        sum_up += occupations[i][j];
        sum_dn += occupations[i][j + occupations[i].size() / 2];
      }
      nup_flt += sum_up * weights[i];
      ndn_flt += sum_dn * weights[i];
    }
  }
  nup = static_cast<int>(round(nup_flt));
  ndn = static_cast<int>(round(ndn_flt));
}

getOccupation()

double getOccupation ( const XmlNode *  nd ) const

private

to be handed a slater_determinant tag from the qbox sample file. Will read in the occupations from the density_matrix subnode add the total which will be returned as output

Definition at line 96 of file WriteEshdf.cpp.

References XmlNode::getChild(), and XmlNode::getValue().

{
  double result = 0.0;
  vector<double> occupancies;
  const XmlNode& densMat = nd->getChild("density_matrix");
  densMat.getValue(occupancies);
  for (int i = 0; i < occupancies.size(); i++)
  {
    result += occupancies[i];
  }
  return result;
}

getPtvs()

vector< double > getPtvs ( const XmlNode &  qeXml )

private

to be handed the top level xml node from data-file-schema, reads the reciprocal lattice vectors from the xml, calculates the primitive translation vectors and puts them sequentially in a 9 element long std::vector

Definition at line 660 of file WriteEshdf.cpp.

References qmcplusplus::det(), XmlNode::getChild(), and XmlNode::getValue().

{
  const XmlNode& reciprocal_lattice_xml = qeXml.getChild("output").getChild("basis_set").getChild("reciprocal_lattice");
  vector<double> rlv; // reciprocal_lattice_vectors
  reciprocal_lattice_xml.getChild("b1").getValue(rlv);
  reciprocal_lattice_xml.getChild("b2").getValue(rlv);
  reciprocal_lattice_xml.getChild("b3").getValue(rlv);

  const double det = rlv[0] * (rlv[4] * rlv[8] - rlv[5] * rlv[7]) - rlv[1] * (rlv[3] * rlv[8] - rlv[5] * rlv[6]) +
      rlv[2] * (rlv[3] * rlv[7] - rlv[4] * rlv[6]);
  const double invdet = 1.0 / det;
  vector<double> ptv;
  ptv.push_back(invdet * (rlv[4] * rlv[8] - rlv[5] * rlv[7]));
  ptv.push_back(invdet * (rlv[5] * rlv[6] - rlv[3] * rlv[8]));
  ptv.push_back(invdet * (rlv[3] * rlv[7] - rlv[4] * rlv[6]));
  ptv.push_back(invdet * (rlv[2] * rlv[7] - rlv[1] * rlv[8]));
  ptv.push_back(invdet * (rlv[0] * rlv[8] - rlv[2] * rlv[6]));
  ptv.push_back(invdet * (rlv[1] * rlv[6] - rlv[0] * rlv[7]));
  ptv.push_back(invdet * (rlv[1] * rlv[5] - rlv[2] * rlv[4]));
  ptv.push_back(invdet * (rlv[2] * rlv[3] - rlv[0] * rlv[5]));
  ptv.push_back(invdet * (rlv[0] * rlv[4] - rlv[1] * rlv[3]));
  return ptv;
}

handleDensity()

void handleDensity ( const XmlNode &  qeXml, const std::string &  dir_name, int  spinpol  )

private

to be handed the top level xml node from data-file-schema, also a string naming the directory where the charge-density.hdf5 file is located and an integer spinpol which is 0 for a non spin polarized calculation and 1 for one that is

Definition at line 492 of file WriteEshdf.cpp.

References qmcplusplus::dims, XmlNode::getAttribute(), XmlNode::getChild(), hdf_archive::getShape(), XmlNode::getValue(), hdf_archive::read(), and hdf_archive::readSlabReshaped().

{
  const XmlNode& output_xml = qeXml.getChild("output");
  // read in the grid to put in the charge density
  const XmlNode& basis_set_xml = output_xml.getChild("basis_set");
  const XmlNode& fft_grid      = basis_set_xml.getChild("fft_grid");
  int nr1;
  int nr2;
  int nr3;
  fft_grid.getAttribute("nr1", nr1);
  fft_grid.getAttribute("nr2", nr2);
  fft_grid.getAttribute("nr3", nr3);
  // also figure out how many gvectors there will be
  const XmlNode& ngm_xml = basis_set_xml.getChild("ngm");
  int num_dens_gvecs;
  ngm_xml.getValue(num_dens_gvecs);

  // now open the hdf file and read the necessary quantities
  const string dens_fname = dir_name + "charge-density.hdf5";
  hdf_archive densfile    = openHdfFileForRead(dens_fname);

  vector<int> readDims;
  densfile.getShape<int>("MillerIndices", readDims);
  vector<int> gvecs;
  array<int, 2> shape{readDims[0], readDims[1]};
  densfile.readSlabReshaped(gvecs, shape, "MillerIndices");

  vector<double> dens;
  densfile.read(dens, "rhotot_g");

  vector<double> diffdens;
  if (spinpol == 1)
  {
    densfile.read(diffdens, "rhodiff_g");
  }

  // now need to write everything out
  outfile_.push("density");
  array<int, 2> dims{num_dens_gvecs, 3};
  outfile_.writeSlabReshaped(gvecs, dims, "gvectors");

  vector<int> grid{nr1, nr2, nr3};
  outfile_.write(grid, "mesh");
  outfile_.write(num_dens_gvecs, "number_of_gvectors");

  array<int, 2> dims_dens{num_dens_gvecs, 2};
  if (spinpol == 0)
  {
    outfile_.push("spin_0");
    outfile_.writeSlabReshaped(dens, dims_dens, "density_g");
    outfile_.pop();
  }
  else
  {
    // do spin up
    vector<double> working(2 * num_dens_gvecs);
    for (int i = 0; i < num_dens_gvecs * 2; i++)
    {
      working[i] = (dens[i] + diffdens[i]) / 2.0;
    }
    outfile_.push("spin_0");
    outfile_.writeSlabReshaped(working, dims_dens, "density_g");
    outfile_.pop();

    // do spin down
    for (int i = 0; i < num_dens_gvecs * 2; i++)
    {
      working[i] = (dens[i] - diffdens[i]) / 2.0;
    }
    outfile_.push("spin_1");
    outfile_.writeSlabReshaped(working, dims_dens, "density_g");
    outfile_.pop();
  }
  outfile_.pop(); // get out of density group
}

handleKpt()

void handleKpt ( int  kpt_num, const std::string &  dir_name, KPoint &  kpt, const std::vector< double > &  eigenvalues, double  weight, int  spinpol, int  noncol, const momap_t &  moref  )

private

helper function to read in the data from the wfc(kpt_num).hdf5 files located in directory dir_name. Needs to know the KPoint, its eigenvalues, weight, the flags for spinolarized and noncolinear as well as the momap_t that handles the global mapping from g_vectors at individual k_points to global ones

Definition at line 713 of file WriteEshdf.cpp.

References hdf_archive::close(), qmcplusplus::dims, hdf_archive::getShape(), KPoint::kx, KPoint::ky, KPoint::kz, hdf_archive::open(), hdf_archive::readSlabReshaped(), and hdf_archive::readSlabSelection().

{
  stringstream ss;
  ss << "kpoint_" << kpt_num;
  outfile_.push(ss.str());
  outfile_.write(weight, "weight");

  vector<double> kpt_vector{kpt.kx, kpt.ky, kpt.kz};
  outfile_.write(kpt_vector, "reduced_k");

  //reread gvecs for this kpts to create a map to coeffs
  //then will merge with the maximal gvec map in moref
  //this will fill unused gvecs with zeros
  string fname;
  if (spinpol == 1)
  {
    stringstream ss;
    ss << dir_name << "wfcup" << (kpt_num + 1) << ".hdf5";
    fname = ss.str();
  }
  else
  {
    stringstream ss;
    ss << dir_name << "wfc" << (kpt_num + 1) << ".hdf5";
    fname = ss.str();
  }
  hdf_archive wfc_file = openHdfFileForRead(fname);
  vector<int> readShape;
  wfc_file.getShape<int>("MillerIndices", readShape);
  vector<int> gvecs;
  array<int, 2> shape{readShape[0], readShape[1]};
  wfc_file.readSlabReshaped(gvecs, shape, "MillerIndices");
  wfc_file.close();

  vector<vector<int>> gvs;
  for (int i = 0; i < readShape[0]; i++)
  {
    vector<int> g;
    for (int d = 0; d < readShape[1]; d++)
      g.push_back(gvecs[i * readShape[1] + d]);
    gvs.push_back(g);
  }

  // now do the non-gvectors related stuff
  hdf_archive spin_0_file;
  hdf_archive spin_1_file;

  if (spinpol == 1)
  {
    stringstream ss0;
    ss0 << dir_name << "wfcup" << (kpt_num + 1) << ".hdf5";
    spin_0_file.open(ss0.str());
    stringstream ss1;
    ss1 << dir_name << "wfcdw" << (kpt_num + 1) << ".hdf5";
    spin_1_file.open(ss1.str());
  }
  else
  {
    stringstream ss0;
    ss0 << dir_name << "wfc" << (kpt_num + 1) << ".hdf5";
    spin_0_file.open(ss0.str());
  }

  // set up storage space for the coefficients to be read from the file
  vector<double> upcoefs;
  vector<double> dncoefs;
  int ng = gvs.size();
  if (spinpol == 0)
  {
    if (noncol == 0)
      upcoefs.resize(ng * 2);
    if (noncol == 1)
    {
      upcoefs.resize(ng * 4);
      dncoefs.resize(ng * 2);
    }
  }
  else
  {
    upcoefs.resize(ng * 2);
    dncoefs.resize(ng * 2);
  }

  vector<momap_t> states_up;
  vector<momap_t> states_dn;
  int states_to_loop = eigenvalues.size();
  if (spinpol == 1)
    states_to_loop /= 2;
  for (int state = 0; state < states_to_loop; state++)
  {
    // set the elements of coefs to 0
    for (int i = 0; i < upcoefs.size(); i++)
      upcoefs[i] = 0.0;
    for (int i = 0; i < dncoefs.size(); i++)
      dncoefs[i] = 0.0;

    // do what is necessary to read only this band's coefficients
    array<int, 2> read_from{state, -1};
    if (spinpol == 0)
    {
      spin_0_file.readSlabSelection(upcoefs, read_from, "evc");
    }
    else
    {
      spin_0_file.readSlabSelection(upcoefs, read_from, "evc");
      spin_1_file.readSlabSelection(dncoefs, read_from, "evc");
    }

    if (spinpol == 0)
    {
      if (noncol == 0)
      {
        momap_t moup;
        for (int i = 0; i < ng; i++)
        {
          complex<double> c(upcoefs[i * 2], upcoefs[i * 2 + 1]);
          mopair_t p(gvs[i], c);
          moup.insert(p);
        }
        //now fill in rest of gvecs with zeros
        moup.insert(moref.begin(), moref.end());
        states_up.push_back(moup);
      }
      else
      {
        //dn part of spinor in second half of upcoefs
        for (int i = 0; i < ng * 2; i++)
          dncoefs[i] = upcoefs[i + 2 * ng];
        momap_t moup, modn;
        for (int i = 0; i < ng; i++)
        {
          complex<double> c(upcoefs[i * 2], upcoefs[i * 2 + 1]);
          mopair_t p(gvs[i], c);
          moup.insert(p);
          c = complex<double>(dncoefs[i * 2], dncoefs[i * 2 + 1]);
          p = mopair_t(gvs[i], c);
          modn.insert(p);
        }
        //now fill in rest of gvecs with zeros
        moup.insert(moref.begin(), moref.end());
        modn.insert(moref.begin(), moref.end());
        states_up.push_back(moup);
        states_dn.push_back(modn);
      }
    }
    else
    {
      momap_t moup, modn;
      for (int i = 0; i < ng; i++)
      {
        complex<double> c(upcoefs[i * 2], upcoefs[i * 2 + 1]);
        mopair_t p(gvs[i], c);
        moup.insert(p);
        c = complex<double>(dncoefs[i * 2], dncoefs[i * 2 + 1]);
        p = mopair_t(gvs[i], c);
        modn.insert(p);
      }
      moup.insert(moref.begin(), moref.end());
      modn.insert(moref.begin(), moref.end());
      states_up.push_back(moup);
      states_dn.push_back(modn);
    }
  }
  spin_0_file.close();
  spin_1_file.close();

  //now write to eshdf
  if (kpt_num == 0)
  {
    vector<int> allgvs;
    int nallgvs = moref.size();
    int dim     = moref.begin()->first.size();
    array<int, 2> shape{nallgvs, dim};
    for (auto& v : moref)
    {
      for (int d = 0; d < dim; d++)
        allgvs.push_back(v.first[d]);
    }
    outfile_.writeSlabReshaped(allgvs, shape, "gvectors");
    outfile_.write(nallgvs, "number_of_gvectors");
  }

  for (int state = 0; state < states_up.size(); state++)
  {
    stringstream ss;
    ss << "state_" << state;
    outfile_.push("spin_0");
    outfile_.push(ss.str());
    vector<double> c;
    for (auto& v : states_up[state])
    {
      c.push_back(v.second.real());
      c.push_back(v.second.imag());
    }
    array<int, 2> dims{static_cast<int>(states_up[state].size()), 2};
    outfile_.writeSlabReshaped(c, dims, "psi_g");
    outfile_.pop();
    outfile_.pop();
    if (noncol == 1 || spinpol == 1)
    {
      outfile_.push("spin_1");
      outfile_.push(ss.str());
      vector<double> c;
      for (auto& v : states_dn[state])
      {
        c.push_back(v.second.real());
        c.push_back(v.second.imag());
      }
      array<int, 2> dims{static_cast<int>(states_dn[state].size()), 2};
      outfile_.writeSlabReshaped(c, dims, "psi_g");
      outfile_.pop();
      outfile_.pop();
    }
  }

  // now all the states are written, so write out eigenvalues and number of states
  vector<double> eigval = eigenvalues;
  if (spinpol == 0)
  {
    outfile_.push("spin_0");
    int int_eig_sz = static_cast<int>(eigenvalues.size());
    outfile_.write(int_eig_sz, "number_of_states");
    outfile_.write(eigval, "eigenvalues");
    outfile_.pop();
    if (noncol == 1)
    {
      outfile_.push("spin_1");
      outfile_.write(int_eig_sz, "number_of_states");
      outfile_.write(eigval, "eigenvalues");
      outfile_.pop();
    }
  }
  else // spin polarized case
  {
    int totstates = eigenvalues.size();
    int upstates  = totstates / 2;
    int dnstates  = totstates / 2;

    vector<double> upeig;
    vector<double> dneig;
    for (int i = 0; i < upstates; i++)
    {
      upeig.push_back(eigenvalues[i]);
      dneig.push_back(eigenvalues[i + upstates]);
    }

    outfile_.push("spin_0");
    outfile_.write(upstates, "number_of_states");
    outfile_.write(upeig, "eigenvalues");
    outfile_.pop();

    outfile_.push("spin_1");
    outfile_.write(dnstates, "number_of_states");
    outfile_.write(dneig, "eigenvalues");
    outfile_.pop();
  }

  outfile_.pop(); // get out of the kpoint_ kpt_num group
}

handleSpinGroup()

void handleSpinGroup ( const XmlNode *  nd, double &  nocc, FftContainer &  cont  )

private

to be handed a slater_determinant tag of the qbox sample file. Will read in the occupations and call readInEigFcn for each state, also needs an FftContainer to pass on to readInEigFcn

Definition at line 110 of file WriteEshdf.cpp.

References FftContainer::fullSize, XmlNode::getChild(), XmlNode::getName(), XmlNode::getNumChildren(), and FftContainer::kspace.

{
  nocc             = getOccupation(nd);
  int stateCounter = 0;
  vector<double> eigvals;

  for (int chIdx = 0; chIdx < nd->getNumChildren(); chIdx++)
  {
    if (nd->getChild(chIdx).getName() == "grid_function")
    {
      //cout << "Working on state " << stateCounter << endl;
      stringstream statess;
      statess << "state_" << stateCounter;
      outfile_.push(statess.str());

      // HACK_HACK_HACK!!!
      // QBOX does not write out the eigenvalues for the states in the
      // sample file, so just make sure they are treated in ascending order
      eigvals.push_back(-5000.0 + stateCounter);

      const XmlNode& eigFcnNode = nd->getChild(chIdx);
      readInEigFcn(eigFcnNode, cont);
      // write eigfcn to proper place
      array<int, 2> psig_dims{cont.fullSize, 2};

      vector<double> temp;
      for (int i = 0; i < cont.fullSize; i++)
      {
        temp.push_back(cont.kspace[i][0]);
        temp.push_back(cont.kspace[i][1]);
      }

      outfile_.writeSlabReshaped(temp, psig_dims, "psi_g");
      stateCounter++;
      outfile_.pop();
    }
  }

  outfile_.write(stateCounter, "number_of_states");
  outfile_.write(eigvals, "eigenvalues");
}

openHdfFileForRead()

hdf_archive openHdfFileForRead ( const std::string &  fname )

private

will open the hdf file fname and will make an hdf_archive

Definition at line 477 of file WriteEshdf.cpp.

References hdf_archive::open().

{
  hdf_archive file;
  bool result = file.open(fname);
  if (!result)
  {
    cout << "could not find " << fname << ", make sure espresso is" << endl;
    cout << "compiled with hdf support and you do not explicitly set" << endl;
    cout << "wf_collect=.false. in your input" << endl;
    exit(1);
  }
  return file;
}

operator=()

EshdfFile& operator= ( const EshdfFile &  f )

private

processKPts()

void processKPts ( const XmlNode &  band_structure_xml, const std::vector< double > &  ptvs, std::vector< std::vector< double >> &  eigenvals, std::vector< std::vector< double >> &  occupations, std::vector< KPoint > &  kpts, std::vector< double > &  weights, std::vector< int > &  ngvecs  )

private

to be handed a band_structure xml node from data-file-schema as well as the 9 element vector of ptvs. Then will populate the eigenvalues per kpoint, occupations per kpoint, kpoints, and the weights and number of g-vectors for each kpoint

Definition at line 568 of file WriteEshdf.cpp.

References XmlNode::getAttribute(), XmlNode::getChild(), XmlNode::getName(), XmlNode::getNumChildren(), XmlNode::getValue(), KPoint::kx, KPoint::ky, and KPoint::kz.

{
  for (int i = 0; i < band_structure_xml.getNumChildren(); i++)
  {
    const XmlNode& child_xml = band_structure_xml.getChild(i);
    if (child_xml.getName() == "ks_energies")
    {
      vector<double> eigs;
      vector<double> occs;
      vector<double> kpt_vec;
      KPoint kpt;
      double weight;
      int ngvec;

      const XmlNode& k_point_xml = child_xml.getChild("k_point");
      k_point_xml.getAttribute("weight", weight);
      k_point_xml.getValue(kpt_vec);
      kpt.kx = kpt_vec[0] * ptvs[0] + kpt_vec[1] * ptvs[1] + kpt_vec[2] * ptvs[2];
      kpt.ky = kpt_vec[0] * ptvs[3] + kpt_vec[1] * ptvs[4] + kpt_vec[2] * ptvs[5];
      kpt.kz = kpt_vec[0] * ptvs[6] + kpt_vec[1] * ptvs[7] + kpt_vec[2] * ptvs[8];

      child_xml.getChild("npw").getValue(ngvec);
      child_xml.getChild("eigenvalues").getValue(eigs);
      child_xml.getChild("occupations").getValue(occs);

      eigenvals.push_back(eigs);
      occupations.push_back(occs);
      kpts.push_back(kpt);
      weights.push_back(weight);
      ngvecs.push_back(ngvec);
    }
  }
}

readInEigFcn()

void readInEigFcn ( const XmlNode &  nd, FftContainer &  cont  )

private

to be handed a grid_function tag of the qbox sample file. Will then read in the eigenfunciton (spo) in real-space and use the FftContainer cont to transform them to k-space

Definition at line 152 of file WriteEshdf.cpp.

References FftContainer::executeFFT(), FftContainer::fixKsNorm(), FftContainer::fullSize, XmlNode::getAttribute(), FftContainer::getNx(), FftContainer::getNy(), FftContainer::getNz(), FftContainer::getQboxIndex(), XmlNode::getValue(), FftContainer::rspace, and qmcplusplus::sqrt().

{
  const string type     = nd.getAttribute("type");
  const string encoding = nd.getAttribute("encoding");

  if (encoding != "text")
  {
    cerr << "Don't yet know how to handle encoding of wavefunction values other than text" << endl;
    exit(1);
  }

  vector<double> values;
  nd.getValue(values);
  const double fixnorm = 1 / std::sqrt(static_cast<double>(cont.fullSize));
  if (type == "complex")
  {
    int index = 0;
    for (int ix = 0; ix < cont.getNx(); ix++)
    {
      for (int iy = 0; iy < cont.getNy(); iy++)
      {
        for (int iz = 0; iz < cont.getNz(); iz++)
        {
          const int qbx         = cont.getQboxIndex(ix, iy, iz);
          cont.rspace[index][0] = values[2 * qbx] * fixnorm;
          cont.rspace[index][1] = values[2 * qbx + 1] * fixnorm;
          index++;
        }
      }
    }
  }
  else if (type == "double")
  {
    int index = 0;
    for (int ix = 0; ix < cont.getNx(); ix++)
    {
      for (int iy = 0; iy < cont.getNy(); iy++)
      {
        for (int iz = 0; iz < cont.getNz(); iz++)
        {
          const int qbx         = cont.getQboxIndex(ix, iy, iz);
          cont.rspace[index][0] = values[qbx] * fixnorm;
          cont.rspace[index][1] = 0.0;
          index++;
        }
      }
    }
  }
  //cout << "in readInEigFcn, before fft, real space L2 norm = " << cont.getL2NormRS() << endl;
  cont.executeFFT();
  cont.fixKsNorm(fixnorm);
  //cout << "in readInEigFcn, after fft, k space L2 norm = " << cont.getL2NormKS() << endl;
}

readKptGvecs()

void readKptGvecs ( int  kpt_num, const std::string &  dir_name, int  spinpol, momap_t &  morefmap  )

private

handles reading in the g-vectors for a particular k-point so that they can be inserted into a momap_t. Needs to know the number of the kpoint (kpt_num), the directory where it is located (dir_name) and whether it is spin polarized (spinpol)

Definition at line 684 of file WriteEshdf.cpp.

References hdf_archive::close(), hdf_archive::getShape(), and hdf_archive::readSlabReshaped().

{
  stringstream ss;
  string fname;
  if (spinpol == 1)
    ss << dir_name << "wfcup" << (kpt_num + 1) << ".hdf5";
  else
    ss << dir_name << "wfc" << (kpt_num + 1) << ".hdf5";
  fname                = ss.str();
  hdf_archive wfc_file = openHdfFileForRead(fname);
  vector<int> readShape;
  wfc_file.getShape<int>("MillerIndices", readShape);
  vector<int> gvecs;
  array<int, 2> shape{readShape[0], readShape[1]};
  wfc_file.readSlabReshaped(gvecs, shape, "MillerIndices");
  wfc_file.close();

  for (int i = 0; i < readShape[0]; i++)
  {
    vector<int> gv;
    for (int d = 0; d < readShape[1]; d++)
      gv.push_back(gvecs[i * readShape[1] + d]);
    mopair_t p(gv, std::complex<double>(0, 0));
    morefmap.insert(p);
  }
}

wrapped()

int wrapped ( int  i, int  size  ) const

private

helper indexing function to go from the interval [0,size) to (-size/2,size/2)

Definition at line 84 of file WriteEshdf.cpp.

{
  if (i < size / 2)
  {
    return i;
  }
  else
  {
    return wrapped(i - size, size);
  }
}

writeApplication()

void writeApplication ( const std::string &  appName, int  major, int  minor, int  sub  )

private

writes to the application field of outfile_ (code -> appName), (version -> major,minor,sub)

Definition at line 206 of file WriteEshdf.cpp.

References qmcplusplus::hdf::version.

{
  vector<int> version{major, minor, sub};

  outfile_.push("application");
  string str = appName;
  outfile_.write(str, "code");
  outfile_.write(version, "version");
  outfile_.pop();
}

writeCreator()

void writeCreator ( )

private

writes to the creator field of outfile_ (program_name -> convertpw4qmc), (version -> 0,1,0)

Definition at line 224 of file WriteEshdf.cpp.

References qmcplusplus::hdf::version.

{
  vector<int> version{0, 1, 0};
  outfile_.push("creator");
  string tmp_str = "convertpw4qmc";
  outfile_.write(tmp_str, "program_name");
  outfile_.write(version, "version");
  outfile_.pop();
}

writeFormat()

void writeFormat ( )

private

writes to the format field of outfile_ (format -> ES-HDF)

Definition at line 234 of file WriteEshdf.cpp.

{
  string tmp_str = "ES-HDF";
  outfile_.write(tmp_str, "format");
}

writeQboxAtoms()

void writeQboxAtoms

(

const XmlNode &

qboxSample

)

Definition at line 412 of file WriteEshdf.cpp.

References qmcplusplus::dims, XmlNode::getAttribute(), XmlNode::getChild(), XmlNode::getName(), XmlNode::getNumChildren(), and XmlNode::getValue().

Referenced by main().

{
  const XmlNode& atomset = qboxSample.getChild("atomset");

  //make group
  outfile_.push("atoms");

  map<string, int> SpeciesNameToInt;
  //go through each species, extract:
  //   atomic_number, mass, name, pseudopotential and valence_charge
  //   then write to a file, also set up mapping between species name and number
  int speciesNum = 0;
  for (int i = 0; i < atomset.getNumChildren(); i++)
  {
    if (atomset.getChild(i).getName() == "species")
    {
      const XmlNode& species = atomset.getChild(i);

      string spName            = species.getAttribute("name");
      SpeciesNameToInt[spName] = speciesNum;

      int atomic_number;
      species.getChild("atomic_number").getValue(atomic_number);
      double mass;
      species.getChild("mass").getValue(mass);
      string name = species.getChild("symbol").getValue();
      int val_charge;
      species.getChild("norm_conserving_pseudopotential").getChild("valence_charge").getValue(val_charge);

      stringstream gname;
      gname << "species_" << speciesNum;
      outfile_.push(gname.str());
      outfile_.write(atomic_number, "atomic_number");
      outfile_.write(mass, "mass");
      outfile_.write(val_charge, "valence_charge");
      outfile_.write(name, "name");
      string tmp_str = "unknown";
      outfile_.write(tmp_str, "pseudopotential");
      speciesNum++;
      outfile_.pop();
    }
  }
  outfile_.write(speciesNum, "number_of_species");

  // go through atoms and extract their position and type
  std::vector<int> species_ids;
  std::vector<double> positions;
  int at_num = 0;
  for (int i = 0; i < atomset.getNumChildren(); i++)
  {
    if (atomset.getChild(i).getName() == "atom")
    {
      const XmlNode& atNode = atomset.getChild(i);
      species_ids.push_back(SpeciesNameToInt[atNode.getAttribute("species")]);
      atNode.getChild("position").getValue(positions);
      at_num++;
    }
  }
  array<int, 2> dims{at_num, 3};
  outfile_.writeSlabReshaped(positions, dims, "positions");
  outfile_.write(species_ids, "species_ids");
  outfile_.write(at_num, "number_of_atoms");
  outfile_.pop();
}

writeQboxBoilerPlate()

void writeQboxBoilerPlate

(

const XmlNode &

qboxSample

)

Definition at line 240 of file WriteEshdf.cpp.

References XmlNode::getChild(), and XmlNode::getValue().

Referenced by main().

{
  const string appName       = "qbox";
  const XmlNode& description = qboxSample.getChild("description");
  string desString           = description.getValue();
  int versionStart           = desString.find("qbox-");
  versionStart += 5;
  string versionStr = desString.substr(versionStart);

  const int firstDotIdx  = versionStr.find_first_of('.');
  const int secondDotIdx = versionStr.find_last_of('.');
  const int major        = stoi(versionStr.substr(0, firstDotIdx));
  const int minor        = stoi(versionStr.substr(firstDotIdx + 1, secondDotIdx - firstDotIdx - 1));
  const int sub          = stoi(versionStr.substr(secondDotIdx + 1));

  writeApplication(appName, major, minor, sub);
  writeVersion();
  writeCreator();
  writeFormat();
}

writeQboxElectrons()

void writeQboxElectrons

(

const XmlNode &

qboxSample

)

Definition at line 1044 of file WriteEshdf.cpp.

References qmcplusplus::dims, qmcplusplus::floor(), XmlNode::getAttribute(), XmlNode::getAttributeIndex(), XmlNode::getChild(), XmlNode::getName(), XmlNode::getNumChildren(), KPoint::kx, KPoint::ky, and KPoint::kz.

Referenced by main().

{
  const XmlNode& wfnNode = qboxSample.getChild("wavefunction");
  int nspin, nel;
  wfnNode.getAttribute("nspin", nspin);
  wfnNode.getAttribute("nel", nel);
  const XmlNode& gridNode = wfnNode.getChild("grid");
  int nx, ny, nz;
  gridNode.getAttribute("nx", nx);
  gridNode.getAttribute("ny", ny);
  gridNode.getAttribute("nz", nz);

  FftContainer fftCont(nx, ny, nz);

  vector<KPoint> kpts;
  map<KPoint, const XmlNode*> kptToUpNode;
  map<KPoint, const XmlNode*> kptToDnNode;

  for (int i = 0; i < wfnNode.getNumChildren(); i++)
  {
    if (wfnNode.getChild(i).getName() == "slater_determinant")
    {
      const XmlNode& sdNode = wfnNode.getChild(i);
      int spinIdx           = sdNode.getAttributeIndex("spin");
      string species;
      if (spinIdx >= 0)
      {
        species = sdNode.getAttribute(spinIdx);
      }
      else
      {
        species = "up";
      }

      KPoint kpt;
      stringstream ss;
      ss.str(sdNode.getAttribute("kpoint"));
      ss >> kpt.kx;
      ss >> kpt.ky;
      ss >> kpt.kz;

      bool newKpt = true;
      for (int j = 0; j < kpts.size(); j++)
      {
        if (kpts[j] == kpt)
        {
          newKpt = false;
        }
      }
      if (newKpt)
      {
        kpts.push_back(kpt);
      }

      if (species == "up")
      {
        kptToUpNode[kpt] = std::addressof(sdNode);
      }
      else
      {
        kptToDnNode[kpt] = std::addressof(sdNode);
      }
    }
  }
  outfile_.push("electrons");
  if (kpts.size() > 1)
  {
    std::cout << "Warning: Due to limitations of the current tool, extreme care" << std::endl;
    std::cout << "is required if tiling to a supercell from qbox calculations with" << std::endl;
    std::cout << "multiple k-points.  Specifically spo eigenvalues are not properly" << std::endl;
    std::cout << "included, so improper choice of orbitals may result." << std::endl;
  }
  const int int_kpts_sz = static_cast<int>(kpts.size());
  outfile_.write(int_kpts_sz, "number_of_kpoints");
  outfile_.write(nspin, "number_of_spins");

  double avgNup = 0.0;
  double avgNdn = 0.0;
  // go through kpt by kpt and write
  for (int i = 0; i < kpts.size(); i++)
  {
    stringstream kptElemName;
    kptElemName << "kpoint_" << i;
    outfile_.push(kptElemName.str());

    vector<double> kvector;
    kvector.push_back(kpts[i].kx);
    kvector.push_back(kpts[i].ky);
    kvector.push_back(kpts[i].kz);
    //writeNumsToHDF("numsym", 1, kpt_group);
    //writeNumsToHDF("symgroup", 1, kpt_group);
    outfile_.write(kvector, "reduced_k");
    const double dbl_weight = 1.0 / static_cast<double>(kpts.size());
    outfile_.write(dbl_weight, "weight");

    if (i == 0)
    {
      // figure out the order of the g-vectors
      // write them to an integer array and then put it in the element
      vector<int> gvectors;
      for (int ix = 0; ix < nx; ix++)
      {
        for (int iy = 0; iy < ny; iy++)
        {
          for (int iz = 0; iz < nz; iz++)
          {
            gvectors.push_back(wrapped(ix, nx));
            gvectors.push_back(wrapped(iy, ny));
            gvectors.push_back(wrapped(iz, nz));
          }
        }
      }

      array<int, 2> dims{nx * ny * nz, 3};
      outfile_.writeSlabReshaped(gvectors, dims, "gvectors");
      const int int_tot_num = nx * ny * nz;
      outfile_.write(int_tot_num, "number_of_gvectors");
    }


    // here is where we will read in both the occupations for the
    // kpoint (species depenent if needed)
    // and also the wavefunction (in real space)
    // fourier transform it and write it into the appropriate hdf element
    double nup = 0;
    double ndn = 0;

    outfile_.push("spin_0");
    const XmlNode* upnd = kptToUpNode[kpts[i]];
    handleSpinGroup(upnd, nup, fftCont);
    outfile_.pop();

    if (nspin == 2)
    {
      outfile_.push("spin_1");
      const XmlNode* dnnd = kptToDnNode[kpts[i]];
      handleSpinGroup(dnnd, ndn, fftCont);
      outfile_.pop();
    }
    else
    {
      ndn = nup;
    }

    avgNup += nup / static_cast<double>(kpts.size());
    avgNdn += ndn / static_cast<double>(kpts.size());
  }
  outfile_.pop();
  vector<int> nels;
  nels.push_back(static_cast<int>(std::floor(avgNup + 0.1)));
  nels.push_back(static_cast<int>(std::floor(avgNdn + 0.1)));
  outfile_.write(nels, "number_of_electrons");
  outfile_.pop();
}

writeQboxSupercell()

void writeQboxSupercell

(

const XmlNode &

qboxSample

)

Definition at line 304 of file WriteEshdf.cpp.

References qmcplusplus::dims, XmlNode::getAttribute(), and XmlNode::getChild().

Referenced by main().

{
  // grab the primitive translation vectors from the atomset tag's attributes and put the entries in the vector ptvs
  const XmlNode& atomset   = qboxSample.getChild("atomset");
  const XmlNode& unit_cell = atomset.getChild("unit_cell");

  stringstream ss;
  ss << unit_cell.getAttribute("a") << "  " << unit_cell.getAttribute("b") << "  " << unit_cell.getAttribute("c");

  vector<double> ptvs;
  double temp;
  while (ss >> temp)
  {
    ptvs.push_back(temp);
  }

  // write the ptvs to the supercell group of the hdf file
  array<int, 2> dims{3, 3};
  outfile_.push("supercell");
  outfile_.writeSlabReshaped(ptvs, dims, "primitive_vectors");
  outfile_.pop();
}

writeQEAtoms()

void writeQEAtoms

(

const XmlNode &

qeXml

)

functions to write the atoms section of the eshdf file needs to be handed the top level xml node (data-file-schema or qbox.sample)

Definition at line 352 of file WriteEshdf.cpp.

References qmcplusplus::dims, XmlNode::getAttribute(), XmlNode::getChild(), XmlNode::getName(), XmlNode::getNumChildren(), and XmlNode::getValue().

Referenced by main().

{
  const XmlNode& atomic_species_xml = qeXml.getChild("output").getChild("atomic_species");

  //make group
  outfile_.push("atoms");

  map<string, int> species_name_to_int;
  //go through each species, extract:
  //   name, mass
  //   in future, it would be good to extract atomic_number, and valence_charge
  //   then write to a file, also set up mapping between species name and number

  int species_num = 0;
  for (int i = 0; i < atomic_species_xml.getNumChildren(); i++)
  {
    if (atomic_species_xml.getChild(i).getName() == "species")
    {
      const XmlNode& species_xml   = atomic_species_xml.getChild(i);
      string sp_name               = species_xml.getAttribute("name");
      species_name_to_int[sp_name] = species_num;
      double mass;
      species_xml.getChild("mass").getValue(mass);

      stringstream gname;
      gname << "species_" << species_num;
      outfile_.push(gname.str());
      outfile_.write(sp_name, "name");
      outfile_.write(mass, "mass");
      species_num++;
      outfile_.pop();
    }
  }
  outfile_.write(species_num, "number_of_species");

  const XmlNode& atomic_positions_xml =
      qeXml.getChild("output").getChild("atomic_structure").getChild("atomic_positions");
  // go through atoms and extract their position and type
  std::vector<int> species_ids;
  std::vector<double> positions;
  int at_num = 0;
  for (int i = 0; i < atomic_positions_xml.getNumChildren(); i++)
  {
    if (atomic_positions_xml.getChild(i).getName() == "atom")
    {
      const XmlNode& at_node_xml = atomic_positions_xml.getChild(i);
      species_ids.push_back(species_name_to_int[at_node_xml.getAttribute("name")]);
      at_node_xml.getValue(positions); // this will append the three numbers to the verctor
      at_num++;
    }
  }
  array<int, 2> dims{at_num, 3};
  outfile_.writeSlabReshaped(positions, dims, "positions");
  outfile_.write(species_ids, "species_ids");
  outfile_.write(at_num, "number_of_atoms");

  outfile_.pop();
}

writeQEBoilerPlate()

void writeQEBoilerPlate

(

const XmlNode &

qeXml

)

functions to get the boilerplate written for qbox or quantum espresso calls writeApplication, writeVersion, writeCreator and writeFormat nees to be handed the top level xml node (data-file-schema or qbox.sample)

Definition at line 279 of file WriteEshdf.cpp.

References XmlNode::getAttribute(), and XmlNode::getChild().

Referenced by main().

{
  const string appName       = "espresso";
  const XmlNode& creatorNode = qeXml.getChild("general_info").getChild("creator");
  const string versionStr    = creatorNode.getAttribute("VERSION");
  int minor                  = 0;
  int sub                    = 0;
  const int firstDotIdx      = versionStr.find_first_of('.');
  const int secondDotIdx     = versionStr.find_last_of('.');
  const int major            = getIntsOnly(versionStr.substr(0, firstDotIdx));
  if (firstDotIdx == secondDotIdx) // this means on subersion is provided
  {
    minor = getIntsOnly(versionStr.substr(firstDotIdx + 1));
  }
  else
  {
    minor = getIntsOnly(versionStr.substr(firstDotIdx + 1, secondDotIdx - firstDotIdx - 1));
    sub   = getIntsOnly(versionStr.substr(secondDotIdx + 1));
  }
  writeApplication(appName, major, minor, sub);
  writeVersion();
  writeCreator();
  writeFormat();
}

writeQEElectrons()

void writeQEElectrons	(	const XmlNode &	qeXml,
		const std::string &	dir_name
	)

functions to write the electrons section of the eshdf file need to be handed the top level xml node (data-file-schema or qbox.sample) in the case of espresso, also needs to be handed the directory name where the wfc.hdf5 and charge-density.hdf5 files are located

Definition at line 980 of file WriteEshdf.cpp.

References XmlNode::getChild(), and XmlNode::getValue().

Referenced by main().

{
  // make electrons group in hdf file
  outfile_.push("electrons");
  const XmlNode& output_xml = qeXml.getChild("output");

  // need to figure out if this is spin polarized, or if it is noncolinear
  int spinpol                       = 0;
  int noncol                        = 0;
  const XmlNode& band_structure_xml = output_xml.getChild("band_structure");

  const XmlNode& lsda_xml     = band_structure_xml.getChild("lsda");
  string lsda_string_bool     = lsda_xml.getValue();
  const XmlNode& noncolin_xml = band_structure_xml.getChild("noncolin");
  string noncolin_string_bool = noncolin_xml.getValue();
  if (lsda_string_bool == "true")
    spinpol = 1;
  if (noncolin_string_bool == "true")
    noncol = 1;

  // scrape xml file and associated hdf for density and write out
  handleDensity(qeXml, dir_name, spinpol);

  // read in information about kpts from the xml
  vector<double> ptv = getPtvs(qeXml);
  vector<vector<double>> eigenvals;
  vector<vector<double>> occupations;
  vector<KPoint> kpts;
  vector<double> weights;
  vector<int> ngvecs;
  processKPts(band_structure_xml, ptv, eigenvals, occupations, kpts, weights, ngvecs);

  // write number of kpoints, number of spins and spinors if appropriate
  int int_kpt_sz = static_cast<int>(kpts.size());
  outfile_.write(int_kpt_sz, "number_of_kpoints");
  outfile_.write(noncol, "has_spinors");

  if (noncol == 0)
  {
    int nspins = 1;
    if (spinpol == 1)
      nspins = 2;
    outfile_.write(nspins, "number_of_spins");
  }

  // figure out how many electrons of each spin and write to file
  int nup;
  int ndn;
  getNumElectrons(occupations, weights, nup, ndn, spinpol, noncol);
  vector<int> nels;
  nels.push_back(nup);
  nels.push_back(ndn);
  outfile_.write(nels, "number_of_electrons");

  //find maximal set of gvecs
  momap_t moref;
  for (int i = 0; i < kpts.size(); i++)
    readKptGvecs(i, dir_name, spinpol, moref);
  //moref now has maximal set of gvecs
  for (int i = 0; i < kpts.size(); i++)
    handleKpt(i, dir_name, kpts[i], eigenvals[i], weights[i], spinpol, noncol, moref);
}

writeQESupercell()

void writeQESupercell

(

const XmlNode &

qeXml

)

functions to write the supercell and primitive_vectors sections of the eshdf file. Needs to be handed the top level xml node (data-file-schema or qbox.sample)

Definition at line 327 of file WriteEshdf.cpp.

References qmcplusplus::dims, XmlNode::getChild(), and XmlNode::getValue().

Referenced by main().

{
  // grab the primitive translation vectors the output tag's atomic structure part
  const XmlNode& cell = qeXml.getChild("output").getChild("atomic_structure").getChild("cell");
  const XmlNode& a    = cell.getChild("a1");
  const XmlNode& b    = cell.getChild("a2");
  const XmlNode& c    = cell.getChild("a3");

  stringstream ss;
  ss << a.getValue() << "  " << b.getValue() << "  " << c.getValue();

  vector<double> ptvs;
  double temp;
  while (ss >> temp)
  {
    ptvs.push_back(temp);
  }

  // write the ptvs to the supercell group of the hdf file
  array<int, 2> dims{3, 3};
  outfile_.push("supercell");
  outfile_.writeSlabReshaped(ptvs, dims, "primitive_vectors");
  outfile_.pop();
}

writeVersion()

void writeVersion ( )

private

writes to the version field of outfile_ (version -> 2,1,0)

Definition at line 218 of file WriteEshdf.cpp.

References qmcplusplus::hdf::version.

{
  vector<int> version{2, 1, 0};
  outfile_.write(version, "version");
}

Member Data Documentation

outfile_

qmcplusplus::hdf_archive outfile_

private

Definition at line 38 of file WriteEshdf.h.

---

The documentation for this class was generated from the following files:

• /home/pk7/projects/qmc/for_cron_doxygen/qmcpack/src/QMCTools/WriteEshdf.h
• /home/pk7/projects/qmc/for_cron_doxygen/qmcpack/src/QMCTools/WriteEshdf.cpp