source: src/Element/periodentafel.cpp@ adbeca

/*
 * Project: MoleCuilder
 * Description: creates and alters molecular systems
 * Copyright (C)  2010-2012 University of Bonn. All rights reserved.
 * 
 *
 *   This file is part of MoleCuilder.
 *
 *    MoleCuilder is free software: you can redistribute it and/or modify
 *    it under the terms of the GNU General Public License as published by
 *    the Free Software Foundation, either version 2 of the License, or
 *    (at your option) any later version.
 *
 *    MoleCuilder is distributed in the hope that it will be useful,
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *    GNU General Public License for more details.
 *
 *    You should have received a copy of the GNU General Public License
 *    along with MoleCuilder.  If not, see <http://www.gnu.org/licenses/>.
 */

/** \file periodentafel.cpp
 *
 * Function implementations for the class periodentafel.
 *
 */

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

//#include "CodePatterns/MemDebug.hpp"

#include <cstring>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <sstream>

#include "CodePatterns/Assert.hpp"
#include "CodePatterns/Log.hpp"
#include "element.hpp"
#include "elements_db.hpp"
#include "Helpers/defs.hpp"
#include "ion.hpp"
#include "periodentafel.hpp"

using namespace std;

/************************************* Functions for class periodentafel ***************************/

/** constructor for class periodentafel
 * Initialises start and end of list and resets periodentafel::checkliste to false.
 */
periodentafel::periodentafel(const bool DoLoad)
{
  if (DoLoad) {
    ScanPeriodentafel();
  }
};

/** destructor for class periodentafel
 * Removes every element and afterwards deletes start and end of list.
 * TODO: Handle when elements have changed and store databases then
 */
periodentafel::~periodentafel()
{
  CleanupPeriodtable();
};

/** Adds element to period table list
 * \param *pointer element to be added
 * \return iterator to added element
 */
periodentafel::iterator periodentafel::AddElement(element * pointer)
{
  atomicNumber_t Z = pointer->getAtomicNumber();
  ASSERT(!elements.count(Z), "Element is already present.");
  if (pointer->getAtomicNumber() < 1 && pointer->getAtomicNumber() >= MAX_ELEMENTS)
    ELOG(0, "Invalid Z number!");
  pair<iterator,bool> res = elements.insert(pair<atomicNumber_t,element*>(Z,pointer));
  return res.first;
};

/** Removes element from list.
 * \param *pointer element to be removed
 */
size_t periodentafel::RemoveElement(const element * pointer)
{
  return RemoveElement(pointer->getAtomicNumber());
};

/** Removes element from list.
 * \param Z element to be removed
 */
size_t periodentafel::RemoveElement(atomicNumber_t Z)
{
  return elements.erase(Z);
};

/** Removes every element from the period table.
 */
void periodentafel::CleanupPeriodtable()
{
  for(iterator iter=elements.begin();iter!=elements.end();++iter){
    delete(*iter).second;
  }
  elements.clear();
};

/** Finds an element by its atomic number.
 * If element is not yet in list, returns NULL.
 * \param Z atomic number
 * \return pointer to element or NULL if not found
 */
const element * periodentafel::FindElement(atomicNumber_t Z) const
{
  const_iterator res = elements.find(Z);
  return res!=elements.end()?((*res).second):0;
};

/** Returns the desired ion to a specific element.
 *
 * If the respective element is not in the list, we return
 * NULL.
 * \return pointer to an element or NULL if not found
 */
const element * periodentafel::FindElement(atomicNumber_t Z, const int ionization)
{
  // if not ionization given, fall back to other function
  if (ionization == 0) {
    return FindElement(Z);
  }
  // element present?
  const_iterator elementiter = elements.find(Z);
  if (elementiter == elements.end())
    return NULL;
  const element & element_base = *(elementiter->second);
        
  // element has already got ions?
  IonsPerElement::iterator setiter =
      ions.find(Z);
  if (setiter != ions.end()) {
    // yes, found ion list
    ionSet::const_iterator res = setiter->second.find(ionization);
        
        if (res != setiter->second.end()) {
          // ion present already
          element * const _ion = res->second;
          return _ion;
        } else {
          // ion not present yet
          ion * const _ion = new ion(element_base, ionization);
          // insert ion
          setiter->second.insert( std::make_pair( ionization, _ion) );
          return _ion;
        }
  } else {
    // no ions yet, create map
    ion * const _ion = new ion(element_base, ionization);
    std::pair<IonsPerElement::iterator,bool> inserter =
        ions.insert( std::make_pair(Z, ionSet()) );
    // insert ion
    ASSERT( inserter.second,
        "periodentafel::FindElement() - could not insert new ionSet to element.");
    inserter.first->second.insert( std::make_pair( ionization, _ion) );
    return _ion;
  }
  return NULL;
}

/** Finds an element by its atomic number.
 * If element is not yet in list, datas are asked and stored in database.
 * \param shorthand chemical symbol of the element, e.g. H for hydrogene
 * \return pointer to element
 */
const element * periodentafel::FindElement(const string &shorthand) const
{
  element *res = 0;
  for(const_iterator iter=elements.begin();iter!=elements.end();++iter) {
    if((*iter).second->getSymbol() == shorthand){
      res = (*iter).second;
      break;
    }
  }
  return res;
};

/** Asks for element number and returns pointer to element
 * \return desired element or NULL
 */
const element * periodentafel::AskElement() const
{
  const element * walker = NULL;
  int Z;
  do {
    std::cout << "Atomic number Z: ";
    std::cin >> Z;
    walker = this->FindElement(Z);  // give type
  } while (walker == NULL);
  return walker;
};

/** Asks for element and if not found, presents mask to enter info.
 * \return pointer to either present or newly created element
 */
const element * periodentafel::EnterElement()
{
  atomicNumber_t Z = 0;
  std::cout << "Atomic number: " << Z;
  cin >> Z;
  const element *res = FindElement(Z);
  if (!res) {
    // TODO: make this using the constructor
    std::cout << "Element not found in database, please enter." << std::endl;
    element *tmp = new element;
    tmp->Z = Z;
    std::cout << "Mass: ";
    cin >> tmp->mass;
    std::cout << "Name [max 64 chars]: ";
    cin >> tmp->name;
    std::cout << "Short form [max 3 chars]: ";
    cin >> tmp->symbol;
    AddElement(tmp);
    return tmp;
  }
  return res;
};


/******************** Access to iterators ****************************/
periodentafel::const_iterator periodentafel::begin() const{
  return elements.begin();
}

periodentafel::const_iterator periodentafel::end() const{
  return elements.end();
}

periodentafel::reverse_iterator periodentafel::rbegin() const{
  return reverse_iterator(elements.end());
}

periodentafel::reverse_iterator periodentafel::rend() const{
  return reverse_iterator(elements.begin());
}

/** Prints element data to \a *out.
 * \param *out outstream
 */
void periodentafel::OutputElement(ostream * const out, const element *elem) const
{
  *out << elem->getName() << "\t";
  *out << elem->getSymbol() << "\t";
  *out << elem->getAtomicNumber() << "\t";
  *out << elem->getMass() << "\t";
  *out << elem->getCovalentRadius() << "\t";
  *out << elem->getVanDerWaalsRadius() << std::endl;
  //*out << elem->getSymbol() << "\t"  << fixed << setprecision(11) << showpoint << elem->getMass() << "g/mol\t" << elem->getName() << "\t" << elem->getSymbol() << "\t" << endl;
};


/** Prints period table to given stream.
 * \param output stream
 */
bool periodentafel::Output(ostream * const output) const
{
  if (output != NULL) {
    for(elementSet::const_iterator iter = elements.begin(); iter != elements.end(); ++iter) {
      OutputElement(output, iter->second);
    }
    return true;
  }
  return false;
}

/** Scan periodentafel contents from internal databases.
 *
 */
void periodentafel::ScanPeriodentafel()
{
  {
    stringstream input(elementsDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadElementsDatabase(input);
    ASSERT(status,  "General element initialization failed");
  }
  {
    stringstream input(ElectronegativitiesDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadElectronegativityDatabase(input);
    ASSERT(status, "Electronegativities entry of element initialization failed");
  }
  {
    stringstream input(valenceDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadValenceDatabase(input);
    ASSERT(status, "Valence entry of element initialization failed");
  }
  {
    stringstream input(orbitalsDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadOrbitalsDatabase(input);
    ASSERT(status, "Orbitals entry of element initialization failed");
  }
  {
    stringstream input(HbondangleDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadHBondAngleDatabase(input);
    ASSERT(status, "HBond angle entry of element initialization failed");
  }
  {
    stringstream input(HbonddistanceDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadHBondLengthsDatabase(input);
    ASSERT(status, "HBond distance entry of element initialization failed");
  }
  {
    stringstream input(ColorDB,ios_base::in);
#ifndef NDEBUG
    bool status =
#endif
        LoadColorDatabase(input);
    ASSERT(status, "color entry of element initialization failed");
  }
}

/** Loads element list from file.
 * \param *path to to standard file names
 */
bool periodentafel::LoadPeriodentafel(const char *path)
{
  ifstream input;
  bool status = true;
  bool otherstatus = true;
  char filename[MAXSTRINGSIZE];

  // fill elements DB
  if (strlen(path)+1+strlen(STANDARDELEMENTSDB) > MAXSTRINGSIZE-1)
    ELOG(2, "Generated path '" << path << "' will be too long.");
  filename[0] = '\0';
  strncat(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDELEMENTSDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as elements database.");
  status = status && LoadElementsDatabase(input);
  input.close();
  input.clear();

  // fill valence DB per element
  strncpy(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDELECTRONEGATIVITYDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as electronegativity database.");
  otherstatus = otherstatus && LoadElectronegativityDatabase(input);
  input.close();
  input.clear();

  // fill valence DB per element
  strncpy(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDVALENCEDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as valence database.");
  otherstatus = otherstatus && LoadValenceDatabase(input);
  input.close();
  input.clear();

  // fill orbitals DB per element
  strncpy(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDORBITALDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as orbitals database.");
  otherstatus = otherstatus && LoadOrbitalsDatabase(input);
  input.close();
  input.clear();

  // fill H-BondAngle DB per element
  strncpy(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDHBONDANGLEDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as H bond angle database.");
  otherstatus = otherstatus && LoadHBondAngleDatabase(input);
  input.close();
  input.clear();

  // fill H-BondDistance DB per element
  strncpy(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDHBONDDISTANCEDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as H bond length database.");
  otherstatus = otherstatus && LoadHBondLengthsDatabase(input);
  input.close();
  input.clear();

  // fill color DB per element
  strncpy(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDCOLORDB, MAXSTRINGSIZE-strlen(filename));
  input.open(filename);
  if (!input.fail())
    LOG(0, "Using " << filename << " as color database.");
  otherstatus = otherstatus && LoadColorDatabase(input);
  input.close();
  input.clear();

  if (!otherstatus){
    ELOG(2, "Something went wrong while parsing the other databases!");
  }

  return status;
};

/** load the element info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadElementsDatabase(istream &input)
{
  bool status = true;
  string header1tmp,header2tmp;
//  std::stringstream parsedelements;
  // first parse into a map, so we can revert to old status in case something goes wront
  map<atomicNumber_t,element*> parsedElements;
  if (!input.fail()) {
    getline(input,header1tmp);
    getline(input,header2tmp); // skip first two header lines
    //cout << "First header: " << header1tmp << endl;
    //cout << "Second header: " << header2tmp << endl;
//    parsedelements <<  "Parsed elements:");
    while (!input.eof()) {
      element *neues = new element;
      input >> neues->name;
      //input >> ws;
      input >> neues->symbol;
      //input >> ws;
      input >> neues->period;
      //input >> ws;
      input >> neues->group;
      //input >> ws;
      input >> neues->block;
      //input >> ws;
      input >> neues->Z;
      //input >> ws;
      input >> neues->mass;
      //input >> ws;
      input >> neues->CovalentRadius;
      //input >> ws;
      input >> neues->VanDerWaalsRadius;
      //input >> ws;
      input >> ws;
      //neues->Output((ofstream *)&cout);
      if ((neues->getAtomicNumber() > 0) && (neues->getAtomicNumber() < MAX_ELEMENTS)) {
        parsedElements[neues->Z] = neues;
//        parsedelements << " " << *neues);
      } else {
        ELOG(2, "Detected empty line or invalid element in elements db, discarding.");
//        parsedelements << " <?>");
        delete(neues);
      }
      // when the input is in failed state, we most likely just read garbage
      if(input.fail()) {
        ELOG(2, "Error parsing elements db.");
        status = false;
        break;
      }
    }
  } else {
    ELOG(1, "Could not open the database.");
    status = false;
  }
  //LOG(0, parsedElements.str());

  if (!parsedElements.size())
    status = false;

  if(status){
    for(map<atomicNumber_t,element*>::iterator iter=parsedElements.begin();
                                               iter!=parsedElements.end();
                                               ++iter){
      if (elements.count(iter->first)) {
        // if element already present, replace the old one
        // pointer to old element might still be in use, so we have to replace into the old element
        *(elements[iter->first])=*iter->second;
        delete(iter->second);
      }
      else {
        // no such element in periodentafel... we can just insert
        elements[iter->first] = iter->second;
      }
    }
    // all went well.. we now copy the header
    strncpy(header1,header1tmp.c_str(),MAXSTRINGSIZE);
    header1[MAXSTRINGSIZE-1]=0;
    strncpy(header2,header2tmp.c_str(),MAXSTRINGSIZE);
    header2[MAXSTRINGSIZE-1]=0;
  }

  return status;
}

/** load the electronegativity info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadElectronegativityDatabase(std::istream &input)
{
  char dummy[MAXSTRINGSIZE];
  if (!input.fail()) {
    input.getline(dummy, MAXSTRINGSIZE);
    while (!input.eof()) {
      atomicNumber_t Z;
      input >> Z;
      ASSERT(elements.count(Z), "Element not present");
      input >> ws;
      input >> elements[Z]->Electronegativity;
      input >> ws;
      //LOG(1, "INFO: Element " << Z << " has " << FindElement(Z)->Electronegativity << " valence electrons.");
    }
    return true;
  } else
    return false;
}

/** load the valence info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadValenceDatabase(istream &input)
{
  char dummy[MAXSTRINGSIZE];
  if (!input.fail()) {
    input.getline(dummy, MAXSTRINGSIZE);
    while (!input.eof()) {
      atomicNumber_t Z;
      input >> Z;
      ASSERT(elements.count(Z), "Element not present");
      input >> ws;
      input >> elements[Z]->Valence;
      input >> ws;
      //LOG(3, "INFO: Element " << Z << " has " << FindElement(Z)->Valence << " valence electrons.");
    }
    return true;
  } else
                return false;
}

/** load the orbitals info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadOrbitalsDatabase(istream &input)
{
  char dummy[MAXSTRINGSIZE];
  if (!input.fail()) {
    input.getline(dummy, MAXSTRINGSIZE);
    while (!input.eof()) {
      atomicNumber_t Z;
      input >> Z;
      ASSERT(elements.count(Z), "Element not present");
      input >> ws;
      input >> elements[Z]->NoValenceOrbitals;
      input >> ws;
      //LOG(3, "Element " << Z << " has " << FindElement(Z)->NoValenceOrbitals << " number of singly occupied valence orbitals.");
    }
    return true;
  } else
    return false;
}

/** load the hbond angles info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadHBondAngleDatabase(istream &input)
{
  char dummy[MAXSTRINGSIZE];
  if (!input.fail()) {
    input.getline(dummy, MAXSTRINGSIZE);
    while (!input.eof()) {
      atomicNumber_t Z;
      input >> Z;
      ASSERT(elements.count(Z), "Element not present");
      input >> ws;
      input >> elements[Z]->HBondAngle[0];
      input >> elements[Z]->HBondAngle[1];
      input >> elements[Z]->HBondAngle[2];
      input >> ws;
      //LOG(3, "Element " << (int)tmp << " has " << FindElement((int)tmp)->HBondAngle[0] << ", " << FindElement((int)tmp)->HBondAngle[1] << ", " << FindElement((int)tmp)->HBondAngle[2] << " degrees bond angle for one, two, three connected hydrogens.");
    }
    return true;
  } else
                return false;
}

/** load the hbond lengths info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadHBondLengthsDatabase(istream &input)
{
  char dummy[MAXSTRINGSIZE];
  if (!input.fail()) {
    input.getline(dummy, MAXSTRINGSIZE);
    while (!input.eof()) {
      atomicNumber_t Z;
      input >> Z;
      ASSERT(elements.count(Z), "Element not present");
      input >> ws;
      input >> elements[Z]->HBondDistance[0];
      input >> elements[Z]->HBondDistance[1];
      input >> elements[Z]->HBondDistance[2];
      input >> ws;
      //LOG(3, "Element " << (int)tmp << " has " << FindElement((int)tmp)->HBondDistance[0] << " Angstrom typical distance to hydrogen.");
    }
    return true;
  } else
                return false;
}

/** load the color info.
 * \param *input stream to parse from
 * \return true - parsing successful, false - something went wrong
 */
bool periodentafel::LoadColorDatabase(istream &input)
{
  char dummy[MAXSTRINGSIZE];
  if (!input.fail()) {
    input.getline(dummy, MAXSTRINGSIZE);
    while (!input.eof()) {
      atomicNumber_t Z;
      input >> Z;
      ASSERT(elements.count(Z), "Element not present");
      input >> ws;
      input >> dummy;
      {
        int tmpcolor;   // char here will only parse a single char (i.e. only "2" out of "255")
        for (int i=0;i<3;++i) {
          input >> ws;
          input >> tmpcolor;
          elements[Z]->color[i] = (unsigned char)tmpcolor;
        }
      }
      input >> ws;
//      {
//        const element * tmp = FindElement(Z);
//        LOG(0, "Element " << tmp->getName() << " has ("
//            << (int)tmp->color[0] << "," << (int)tmp->color[1] << "," << (int)tmp->color[2]
//            << ") colors.");
//      }
    }
    return true;
  } else
    return false;
}

/** Stores element list to file.
 */
bool periodentafel::StorePeriodentafel(const char *path) const
{
  bool result = true;
  ofstream f;
  char filename[MAXSTRINGSIZE];

  if (strlen(path)+1+strlen(STANDARDELEMENTSDB) > MAXSTRINGSIZE-1)
    ELOG(2, "Generated path '" << path << "' will be too long.");
  filename[0] = '\0';
  strncat(filename, path, MAXSTRINGSIZE);
  strncat(filename, "/", MAXSTRINGSIZE-strlen(filename));
  strncat(filename, STANDARDELEMENTSDB, MAXSTRINGSIZE-strlen(filename));
  f.open(filename);
  if (f != NULL) {
    f << header1 << endl;
    f << header2 << endl;
    for(const_iterator iter=elements.begin();iter!=elements.end();++iter){
      OutputElement(&f, iter->second);
    }
    f.close();
    return true;
  } else
    return result;
};

/** Comparison operator for periodentafel.
 *
 * @param other other instance to compare to
 * @return true when both contain same elements
 */
bool periodentafel::operator==(const periodentafel &other) const
{
  // there are only pointers in the elementSet, hence we have to compare ourselves
  if (elements.size() != other.elements.size()) return false;
  const_iterator iter = elements.begin();
  const_iterator otheriter = other.elements.begin();
  for (;(iter != elements.end()) && (otheriter != other.elements.end());
      ++iter, ++otheriter) {
    bool status = true;
    status = status && (iter->first == otheriter->first);
    status = status && (*(iter->second) == *(otheriter->second));
    if (!status) {
      std::cout << *(iter->second) << " not equal to " << *(otheriter->second) << "." << std::endl;
      return false;
    }
//    else
//      std::cout << (iter->second)->getName() << " are equal to " << (otheriter->second)->getName() << "." << std::endl;
  }
  if (strncmp(header1, other.header1, MAXSTRINGSIZE) != 0) return false;
  if (strncmp(header2, other.header2, MAXSTRINGSIZE) != 0) return false;
  return true;
}