gaussbas.h

//
// gaussbas.h
//
// Copyright (C) 1996 Limit Point Systems, Inc.
//
// Author: Curtis Janssen <cljanss@limitpt.com>
// Maintainer: LPS
//
// This file is part of the SC Toolkit.
//
// The SC Toolkit is free software; you can redistribute it and/or modify
// it under the terms of the GNU Library General Public License as published by
// the Free Software Foundation; either version 2, or (at your option)
// any later version.
//
// The SC Toolkit 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 Library General Public License for more details.
//
// You should have received a copy of the GNU Library General Public License
// along with the SC Toolkit; see the file COPYING.LIB.  If not, write to
// the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
//
// The U.S. Government is granted a limited license as per AL 91-7.
//
 
#ifndef _chemistry_qc_basis_gaussbas_h
#define _chemistry_qc_basis_gaussbas_h
 
#include <vector>
#include <iostream>
 
#include <util/state/state.h>
#include <util/keyval/keyval.h>
#include <math/scmat/matrix.h>
#include <math/scmat/vector3.h>
#include <chemistry/molecule/molecule.h>
#include <chemistry/qc/basis/gaussshell.h>
#include <math/mmisc/grid.h>
 
#ifdef MPQC_NEW_FEATURES
#include "mpqc/range.hpp"
#endif
 
namespace sc {
 
class BasisFileSet;
class Integral;
 
class CartesianIter;
class SphericalTransformIter;
 
 
class GaussianBasisSet: virtual public SavableState
#ifdef MPQC_NEW_FEATURES
, virtual public DescribedXMLWritable
#endif
{
  public:
 
    class Shell : public GaussianShell {
      public:
        Shell(const GaussianBasisSet* basis, unsigned int center, const GaussianShell& shell);
        Shell(const Shell& other) :
          GaussianShell(static_cast<const GaussianShell&>(other)),
          basis_(other.basis_),
          center_(other.center_) {}
        ~Shell() {}
 
        Shell& operator=(const Shell& other) {
          static_cast<GaussianShell&>(*this) = static_cast<const GaussianShell&>(other);
          basis_ = other.basis_;
          center_ = other.center_;
          return *this;
        }
 
        bool equiv(const Shell& s) const;
 
        const GaussianBasisSet* basis() const { return basis_; }
        unsigned int center() const { return center_; }
 
#ifdef MPQC_NEW_FEATURES
        virtual boost::property_tree::ptree& write_xml(boost::property_tree::ptree& parent, const XMLWriter& writer);
#endif
 
      private:
        friend class GaussianBasisSet;
 
        static ClassDesc class_desc_;
 
        const GaussianBasisSet* basis_;
        unsigned int center_;
    };
 
 
  private:
    friend class UnionBasisSet;
 
    // primary data
    std::string name_;    // non-empty if keyword "name" was provided
    std::string label_;   // same as name_ if name_ non-empty, else something else
 
  protected: Ref<Molecule> molecule_;
  private: std::vector<Shell> shells_;
 
    // computes secondary data from the primary data
    void init2(int skip_ghosts=0,bool include_q=0);
 
    // secondary data (expensive or convenient to have precomputed)
    std::vector<int> shell_to_function_;
    std::vector<int> function_to_shell_;
    std::vector<int> shell_to_primitive_;
    std::vector<int> center_to_shell_;
    std::vector<unsigned int> center_to_nshell_;
    std::vector<unsigned int> center_to_nfunction_;
    unsigned int nbasis_;
    unsigned int nprim_;
    bool has_pure_;
    // old-school SCMatrix stuff, may go away
    Ref<SCMatrixKit> matrixkit_;
    Ref<SCMatrixKit> so_matrixkit_;
    RefSCDimension basisdim_;
 
    // Counts shells in this basis for this chemical element
    int count_shells_(Ref<KeyVal>& keyval, const char* elemname, const char* sbasisname, BasisFileSet& bases,
              int havepure, int pure, bool missing_ok);
    // constructs shells on @c atom
    void get_shells_(unsigned int atom, Ref<KeyVal>& keyval, const char* elemname, const char* sbasisname, BasisFileSet& bases,
             int havepure, int pure, bool missing_ok);
    // Counts shells in an even-tempered primitive basis
    int count_even_temp_shells_(Ref<KeyVal>& keyval, const char* elemname, const char* sbasisname,
                                int havepure, int pure);
    // Constructs an even-tempered primitive basis
    void get_even_temp_shells_(unsigned int atom, Ref<KeyVal>& keyval, const char* elemname, const char* sbasisname,
                               int havepure, int pure);
    // Constructs basis set specified as an array of shells
    void recursively_get_shell(unsigned int atom,Ref<KeyVal>&,
                               const char*,const char*,BasisFileSet&,
                               int,int,int,bool missing_ok);
 
    // finishes up the KeyVal constructor
    void init(Ref<Molecule>&,Ref<KeyVal>&,
              BasisFileSet&,
              int have_userkeyval,
              int pure);
 
    void validate_point_group() const;
 
  protected:
    /* This CTOR leaves it up to the derived class to completely
        initialize the basis set. */
    GaussianBasisSet();
    virtual void set_matrixkit(const Ref<SCMatrixKit>&);
 
    void init(std::string name,
              std::string label,
              const Ref<Molecule> &molecule,
              const std::vector<Shell>& shell);
 
  public:
    class ValueData {
      protected:
        CartesianIter **civec_;
        SphericalTransformIter **sivec_;
        int maxam_;
      public:
        ValueData(const Ref<GaussianBasisSet> &, const Ref<Integral> &);
        ~ValueData();
        CartesianIter **civec() { return civec_; }
        SphericalTransformIter **sivec() { return sivec_; }
    };
 
 
    GaussianBasisSet(const Ref<KeyVal>&);
 
    GaussianBasisSet(const Ref<Molecule>& molecule,
                     const std::vector<GaussianShell>& shells,
                     const std::vector<unsigned int>& shell_to_center,
                     std::string name = std::string(),
                     std::string label = std::string());
 
    GaussianBasisSet(const GaussianBasisSet&);
 
    static Ref<GaussianBasisSet> unit();
 
 
 
    GaussianBasisSet(StateIn& si);
    void save_data_state(StateOut& so);
 
#ifdef MPQC_NEW_FEATURES
    virtual boost::property_tree::ptree& write_xml(boost::property_tree::ptree& parent, const XMLWriter& writer);
#endif
 
 
    virtual ~GaussianBasisSet();
 
    Ref<GaussianBasisSet> operator+(const Ref<GaussianBasisSet>& B);
    GaussianBasisSet& operator=(const GaussianBasisSet& A);
 
    const std::string& name() const { return name_; }
    const std::string& label() const { if (!name_.empty()) return name_; else return label_; }
 
    Ref<Molecule> molecule() const { return molecule_; }
    Ref<SCMatrixKit> matrixkit() { return matrixkit_; }
    Ref<SCMatrixKit> so_matrixkit() { return so_matrixkit_; }
    RefSCDimension basisdim() { return basisdim_; }
 
    unsigned int ncenter() const;
    unsigned int nshell() const { return shells_.size(); }
    int nshell_on_center(int icenter) const;
    int shell_on_center(int icenter, int shell) const;
 
 
    int shell_to_center(int ishell) const { return shells_[ishell].center(); }
    unsigned int nbasis() const { return nbasis_; }
    unsigned int nbasis_on_center(int icenter) const;
    unsigned int nprimitive() const { return nprim_; }
    bool has_pure() const { return has_pure_; }
 
    unsigned int max_nfunction_in_shell() const;
    unsigned int max_ncartesian_in_shell(int aminc=0) const;
    unsigned int max_nprimitive_in_shell() const;
    unsigned int max_angular_momentum() const;
    unsigned int max_ncontraction() const;
    unsigned int max_am_for_contraction(int con) const;
    unsigned int max_cartesian() const;
 
    int shell_to_function(int i) const { return shell_to_function_[i]; }
    int function_to_shell(int i) const;
 
#ifdef MPQC_NEW_FEATURES
    mpqc::range range(int s) const;
#endif
 
    const Shell& operator()(int i) const { return shells_[i]; }
    Shell& operator()(int i) { return shells_[i]; }
    const Shell& operator[](int i) const { return shells_[i]; }
    Shell& operator[](int i) { return shells_[i]; }
    const Shell& shell(int i) const { return shells_[i]; }
    Shell& shell(int i) { return shells_[i]; }
 
    const Shell& operator()(int icenter,int ishell) const;
    Shell& operator()(int icenter,int ishell);
    const Shell& shell(int i,int j) const { return operator()(i,j); }
    Shell& shell(int i,int j) { return operator()(i,j); }
 
    int find(int C, const GaussianShell& S) const;
 
    double r(int icenter,int xyz) const;
 
    int values(const SCVector3& r, ValueData *, double* basis_values) const;
    int grad_values(const SCVector3& r, ValueData *,
                    double*g_values,double* basis_values=0) const;
    int hessian_values(const SCVector3& r, ValueData *, double *h_values,
                       double*g_values=0,double* basis_values=0) const;
    int shell_values(const SCVector3& r, int sh,
                     ValueData *, double* basis_values) const;
    int grad_shell_values(const SCVector3& r, int sh,
                          ValueData *,
                          double*g_values, double* basis_values=0) const;
    int hessian_shell_values(const SCVector3& r, int sh,
                       ValueData *, double *h_values,
                       double*g_values=0,double* basis_values=0) const;
 
    int equiv(const Ref<GaussianBasisSet> &b);
 
    void print_brief(std::ostream& =ExEnv::out0()) const;
    virtual void print(std::ostream& =ExEnv::out0()) const;
};
 
Ref<GaussianBasisSet>
operator+(const Ref<GaussianBasisSet>& A, const Ref<GaussianBasisSet>& B);
 
int
ishell_on_center(int icenter, const Ref<GaussianBasisSet>& bs,
                 const GaussianShell& A);
 
Ref<GaussianBasisSet>
operator+(const Ref<GaussianBasisSet>& A, const Ref<GaussianBasisSet>& B);
 
std::vector<unsigned int> operator<<(const GaussianBasisSet& B, const GaussianBasisSet& A);
std::vector<int> map(const GaussianBasisSet& B, const GaussianBasisSet& A);
 
class GaussianBasisSetMap : public RefCount {
  public:
    GaussianBasisSetMap(const Ref<GaussianBasisSet>& Source, const Ref<GaussianBasisSet>& Target);
    ~GaussianBasisSetMap();
 
    const Ref<GaussianBasisSet>& source() const { return source_; }
    const Ref<GaussianBasisSet>& target() const { return target_; }
 
    int map_shell(int s) const;
    int map_function(int f) const;
    int nfblock() const;
    int fblock_to_function(int b) const;
    int fblock_size(int b) const;
 
  private:
    Ref<GaussianBasisSet> source_;
    Ref<GaussianBasisSet> target_;
 
    // maps
    std::vector<int> smap_;
    std::vector<int> fmap_;
    std::vector<int> fblock_to_function_;
    std::vector<int> fblock_size_;
 
};
 
class WriteBasisGrid : public WriteVectorGrid {
  private:
    struct BasisFunctionMap : public DimensionMap {
      std::vector<int> map;
      int operator()(int o) const { return map[o]; }
      std::size_t ndim() const { return map.size(); }
    };
    Ref<GaussianBasisSet> basis_;
    Ref<Integral> integral_;
 
  protected:
    BasisFunctionMap bfmap_;
 
    void label(char* buffer);
    Ref<Molecule> get_molecule() { return basis_->molecule(); }
    void calculate_values(const std::vector<SCVector3>& Points, std::vector<double>& Values);
    const DimensionMap& dimension_map() const { return bfmap_; }
    std::size_t ndim() const { return bfmap_.ndim(); }
    void initialize();
 
    static Ref<KeyVal> process_keyval_for_base_class(const Ref<KeyVal>& kv);
 
  public:
 
    WriteBasisGrid(const Ref<KeyVal>& keyval);
    WriteBasisGrid(
        const Ref<GaussianBasisSet>& basis,
        const Ref<Grid>& grid,
        std::string gridformat,
        std::string gridfile
    );
 
    virtual ~WriteBasisGrid();
 
#ifdef MPQC_NEW_FEATURES
    virtual boost::property_tree::ptree& write_xml(boost::property_tree::ptree& parent, const XMLWriter& writer);
#endif
 
};
 
 
// end of addtogroup ChemistryBasisGaussian
 
} // end namespace sc
 
 
#endif
 
// Local Variables:
// mode: c++
// c-file-style: "CLJ"
// End: