Slater_Orbital.h

#ifndef QNANO_NEW_SLATER_ORBITAL_DEFINED_H
#define QNANO_NEW_SLATER_ORBITAL_DEFINED_H

#include "orbitals/Orbital_Index.h"
#include "tools/Integrator.h"
#include "orbitals/Spherical_Harmonics.h"

class Slater_Orbital: public Integrable_Function{
public:
  //Wave function: N*r^a*exp(-b*r/a_Bohr)*Y(angular)
  //Note: Bohr radius: a_Bohr = 0.0529177 nm
  size_t angular;
  double a, b, N;

  static double f(double *x, size_t dim, void *params){
    const Slater_Orbital *cl = (const Slater_Orbital*)params;
    double r=sqrt(x[0]*x[0]+x[1]*x[1]+x[2]*x[2]);
    return cl->N*pow(r,cl->a)*exp(-cl->b*r/0.0529177)*Spherical_Harmonics::Y(cl->angular,x);
  }
  
  void normalize_Riemann(double lim=3., size_t Ndiscr=300);
  void normalize_Simpson(double lim=3., size_t Ndiscr=300);
  void normalize(double lim=3., size_t Ndiscr=300);

  void print(std::ostream &os=std::cout)const;

  Slater_Orbital(): Integrable_Function(f, 3, this), a(0), b(1.), N(-1.) {
  }
  Slater_Orbital(size_t angular_, double a_, double b_, double N_=-1.)
    : Integrable_Function(f, 3, this), angular(angular_%10), a(a_), b(b_), N(N_){
  }

};

#endif