Homec4science
Diffusion lammps (lammps-icms) Browse

lammps/lib/atc/ElectronChargeDensity.h
db1e17197fcflammps-icms

Last Change
Raw File

ElectronChargeDensity.h
View Options

#ifndef ELECTRON_DENSITY_H
#define ELECTRON_DENSITY_H
#include <map>
#include <string>
#include <fstream>
#include "ATC_TypeDefs.h"
#include "Function.h"
const double tol = 1.e-8;
namespace ATC {
/**
* @class ElectronChargeDensity
* @brief Base class for models of extrinsic electric charges
*/
class ElectronChargeDensity
{
public:
ElectronChargeDensity() {};
virtual ~ElectronChargeDensity() {};
virtual bool electron_charge_density(const FIELD_MATS &fields,
DENS_MAT &flux) const { return false; };
virtual void D_electron_charge_density(const FieldName fieldName,
const FIELD_MATS &fields,
DENS_MAT &flux) const
{ throw ATC_Error("Charge density D_electron_charge_density unimplemented function");}
virtual void band_edge_potential(const FIELD_MATS &fields,
DENS_MAT &density) const
{ throw ATC_Error("Charge density band_edge_potential unimplemented function");}
};
//-----------------------------------------------------------------------
/**
* @class ElectronChargeDensityInterpolation
* @brief Class for models of electron charge density as a tabular function of electric potential
*/
class ElectronChargeDensityInterpolation : public ElectronChargeDensity
{
public:
ElectronChargeDensityInterpolation(std::fstream &matfile,std::map<std::string,double> & parameters);
virtual ~ElectronChargeDensityInterpolation() {};
virtual bool electron_charge_density(const FIELD_MATS &fields,
DENS_MAT &flux) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
const DENS_MAT & phi = phi_field->second;
int nNodes = phi.nRows();
flux.reset(nNodes,1,false);
for (int i = 0; i < nNodes; i++) { // a mapping of a vector
flux(i,0) = n_.f(phi(i,0));
}
flux *= -1.;
return true;
};
virtual void D_electron_charge_density(const FieldName field,
const FIELD_MATS &fields,
DENS_MAT &coef) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
const DENS_MAT & phi = phi_field->second;
int nNodes = phi.nRows();
coef.reset(nNodes,1,false);
for (int i = 0; i < nNodes; i++) {
coef(i,0) = n_.dfdt(phi(i,0));
coef(i,0) = n_.dfdt(phi(i,0));
}
coef *= -1.;
}
private:
InterpolationFunction n_;
};
//-----------------------------------------------------------------------
/**
* @class ElectronChargeDensityLinear
* @brief Class for models of electron charge density proportional to electric potential
*/
class ElectronChargeDensityLinear : public ElectronChargeDensity
{
public:
ElectronChargeDensityLinear(std::fstream &matfile,std::map<std::string,double> & parameters);
virtual ~ElectronChargeDensityLinear() {};
virtual bool electron_charge_density(const FIELD_MATS &fields,
DENS_MAT &flux) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
flux = phi_field->second;
flux *= -C_;
return true;
};
virtual void D_electron_charge_density(const FieldName field,
const FIELD_MATS &fields,
DENS_MAT &coef) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
const DENS_MAT & phi = phi_field->second;
int nNodes = phi.nRows();
coef.reset(nNodes,1,false);
coef = -C_;
}
private:
double C_;
};
//-----------------------------------------------------------------------
/**
* @class ElectronChargeDensityExponential
* @brief Class for models of electron charge density dependent on difference between electric potential and the Fermi level n = n_i exp ( (phi-E_i) / kB T)
*/
class ElectronChargeDensityExponential : public ElectronChargeDensity
{
public:
ElectronChargeDensityExponential(std::fstream &matfile,std::map<std::string,double> & parameters);
virtual ~ElectronChargeDensityExponential() {};
double n(const double phi, double T) const
{
return -intrinsicConcentration_*exp((phi-intrinsicEnergy_)/(kBeV_*T));
}
double dndphi(const double phi, double T) const
{
return n(phi,T)/(kBeV_*T);
}
virtual bool electron_charge_density(const FIELD_MATS &fields,
DENS_MAT &density) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
FIELD_MATS::const_iterator T_field = fields.find(TEMPERATURE);
double T = 300;
bool hasTref = (referenceTemperature_ > 0 );
const DENS_MAT & phi = phi_field->second;
int nNodes = phi.nRows();
density.resize(nNodes,1);
if (hasTref) {
T = referenceTemperature_;
for (int i = 0; i < nNodes; i++) {
density(i,0) = n(phi(i,0),T); }
}
else {
const DENS_MAT & temp = T_field->second;
for (int i = 0; i < nNodes; i++) {
density(i,0) = n(phi(i,0),temp(i,0)); }
}
density *= -1.;
return true;
};
virtual void D_electron_charge_density(const FieldName field,
const FIELD_MATS &fields,
DENS_MAT &coef) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
FIELD_MATS::const_iterator T_field = fields.find(TEMPERATURE);
double T = 300;
bool hasTref = (referenceTemperature_ > 0 );
const DENS_MAT & phi = phi_field->second;
int nNodes = phi.nRows();
coef.resize(nNodes,1);
if (hasTref) {
T = referenceTemperature_;
for (int i = 0; i < nNodes; i++) {
coef(i,0) = dndphi(phi(i,0),T); }
}
else {
const DENS_MAT & temp = T_field->second;
for (int i = 0; i < nNodes; i++) {
coef(i,0) = dndphi(phi(i,0),temp(i,0)); }
}
coef *= -1.;
};
protected:
double intrinsicConcentration_,intrinsicEnergy_;
double referenceTemperature_;
};
//-----------------------------------------------------------------------
/**
* @class ElectronChargeDensityFermiDirac
* @brief Class for models of electron charge density based on Fermi-Dirac statistics
*/
class ElectronChargeDensityFermiDirac : public ElectronChargeDensity
{
public:
ElectronChargeDensityFermiDirac(std::fstream &matfile,std::map<std::string,double> & parameters);
virtual ~ElectronChargeDensityFermiDirac() {};
double fermi_dirac(const double E, const double T) const
{
double f = 1.0;
if (T > 0) f = 1.0 / ( exp((E-Ef_)/kBeV_/T)+1.0 );
else if (E > Ef_) f = 0;
return f;
};
virtual bool electron_charge_density(const FIELD_MATS &fields,
DENS_MAT &density) const
{
// psi : the inhomogeneous solution
FIELD_MATS::const_iterator psi_field = fields.find(ELECTRON_WAVEFUNCTION);
const DENS_MAT & psi = psi_field->second;
FIELD_MATS::const_iterator psis_field = fields.find(ELECTRON_WAVEFUNCTIONS);
// if (psis_field==fields.end())
//throw ATC_Error("Wavefunctions not defined");
const DENS_MAT & psis = psis_field->second;
FIELD_MATS::const_iterator E_field = fields.find(ELECTRON_WAVEFUNCTION_ENERGIES);
const DENS_MAT & Es = E_field->second;
FIELD_MATS::const_iterator T_field = fields.find(ELECTRON_TEMPERATURE);
const DENS_MAT & Ts = T_field->second;
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
const DENS_MAT & phi = phi_field->second;
int nNodes = psi.nRows();
density.reset(nNodes,1);
double T = referenceTemperature_;
int count = 0;
for (int i = 0; i < nNodes; i++) {
if (!hasReferenceTemperature_) { T = Ts(i,0); }
int j = 0;
for (j = 0; j < psis.nCols(); j++) {
double E = Es(j,0); // Eigenvalue
double f = fermi_dirac(E,T);
if (f < tol) break;
else count++;
density(i,0) -= psis(i,j)*psis(i,j)*f; // < 0
}
if (donorIonization_) {
double E = -1.0* phi(i,0);// units [eV] E = - |e| phi
if ( E + Eb_ > Ef_+Ed_) density(i,0) += Nd_; // > 0
}
}
return true;
};
virtual void D_electron_charge_density(const FieldName fieldName,
const FIELD_MATS &fields,
DENS_MAT &coef) const
{
FIELD_MATS::const_iterator phi_field = fields.find(ELECTRIC_POTENTIAL);
const DENS_MAT & phi = phi_field->second;
int nNodes = phi.nRows();
coef.reset(nNodes,1,false);
}
virtual void band_edge_potential(const FIELD_MATS &fields,
DENS_MAT &density) const
{
FIELD_MATS::const_iterator p_field = fields.find(ELECTRIC_POTENTIAL);
const DENS_MAT & phi = p_field->second;
int nNodes = phi.nRows();
density.reset(nNodes,1,false);
density = Eb_;
};
protected:
double Ef_;
double referenceTemperature_;
double Ed_, Nd_;
double Eb_;
bool hasReferenceTemperature_, donorIonization_;
};
}
#endif
c4science · Help