SMILE  v2.5
Schwarzschild Modelling Interactive expLoratory Environment
Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
smile::CMassModel Class Reference

class defining a spherical mass model. More...

#include <massmodel.h>

Public Types

enum  MODELSTATUS {
  MS_OK, MS_ARRAY_NOT_EQ, MS_ARRAY_SHORT, MS_RADII_NOT_ASC,
  MS_MASS_NOT_ASC, MS_DF_UNDEFINED, MS_DF_NEGATIVE, MS_DIFCOEF_DIVERGE
}
 internal status or error code of the model More...
 
enum  FITTYPE { FT_MASS, FT_DF_ENERGY, FT_DF_PHASEVOL }
 modes of operation for fitMassModel More...
 

Public Member Functions

 CMassModel (const CMassModel &)
 copy constructor (makes deep copies of gsl vectors and splines)
 
 CMassModel (const std::vector< double > &initRadii, const std::vector< double > &initMass, const CMassModel *model_pot=NULL)
 create the model from an array of { r, m(r) } pairs and pre-computes all necessary quantities. More...
 
MODELSTATUS getStatus () const
 < return the status code set in the constructor
 
bool initDistrFuncEnergy (const std::vector< double > &initE, const std::vector< double > &initF)
 initialize spline for externally supplied distribution function in energy. More...
 
bool computeNewPotential ()
 reinitialize potential from the distribution function in energy assigned before
 
double mass (double rad) const
 return enclosed mass at given radius
 
double rad (double mass) const
 find radius for enclosed mass
 
double dens (double rad) const
 density at given radius
 
double pot (double rad) const
 potential at given radius
 
double drhodphi (double E) const
 d(rho)/d(phi) at given energy E=phi
 
double d2rhodphi2 (double E) const
 d^2(rho)/d(phi)^2 at given energy E=phi
 
double distrf2 (double E) const
 distribution function of (negative) energy, computed directly from integral with d^2(rho)/d(Phi)^2
 
double distrf (double E) const
 distribution function of (negative) energy, spline approximation to distrf2
 
double ded (double E) const
 differential energy distribution g(E)
 
double phasevol (double E) const
 phase volume for energy below E
 
double masse (double E) const
 mass of particles having energies <E
 
double trad (double E) const
 period of radial orbit at energy E
 
double tepi (double E) const
 epicyclic period at energy E (with circular radius rcirc(E))
 
double rmax (double E) const
 max radius at given energy
 
double rcirc (double E) const
 radius of circular orbit at energy E
 
double lcirc (double E) const
 ang.momentum at circular orbit
 
double veldisp (double rad) const
 velocity dispersion at given radius
 
double densproj (double rad) const
 projected density at given projected radius
 
double losvd (double rad) const
 line-of-sight velocity dispersion at a given projected radius
 
double samplev (double r) const
 sample a velocity from distribution function at a given radius
 
void DifCoefVel (double r, double v, double *dvpar, double *dv2par, double *dv2per, bool exact=false) const
 compute local (position-dependent) velocity diffusion coefficients
 
void DifCoefEL (double E, double *DEE=NULL, double *DE=NULL, double *DLL=NULL, double *DRRoverR=NULL, bool exact=false) const
 compute averaged diffusion coefs depending on E. More...
 
void getMassProfile (std::vector< double > *radii, std::vector< double > *inmass) const
 obtain the original mass profile
 
bool hasDensityCore (double *rho0_, double *drhodr_, double *alpha_) const
 tell whether the model has no BH and a constant-density core, and if yes, give its parameters
 
double getTotalMass () const
 return total mass of the model
 
double getPotentialAtCenter () const
 return potential at r=0 without the contribution of central point mass (if any)
 
double scaledE (double E) const
 return log-scaled energy variable (ranging from -infinity to infinity)
 
double unscaledE (double scaledE) const
 return energy for the given log-scaled energy variable
 
double scaledEder (double E) const
 return d(scaledE)/dE
 
double getMinEnergyOrbitAvg (double dcmult) const
 
void DifVexact (double E, double Phi, double *I0, double *I12, double *I32) const
 compute integrals for velocity diffusion coefficients using exact quadrature
 
void DifVapprox (double E, double Phi, double *I0, double *I12, double *I32) const
 compute integrals for velocity diffusion coefficients using interpolation
 
double computeSplineInScaledEnergy (double E, const gsl_spline *spl) const
 compute spline-interpolated quantity in internally scaled energy
 
void computeSplineDerivInScaledEnergy (double E, const gsl_spline *spl, double *value, double *deriv, double *deriv2) const
 compute value and optionally first and second derivative of a spline-interpolated quantity in internally scaled energy
 

Private Member Functions

CMassModeloperator= (const CMassModel &)
 assignment operator is forbidden, use copy constructor
 
void initSplineMass (const std::vector< double > &initRadii, const std::vector< double > &initMass)
 initialize M(r) spline
 
void initSplinePhi (const CMassModel *model_pot)
 compute potential and initialize Phi(r) spline, or use external potential if model_pot!=NULL More...
 
void initEddingtonDF ()
 compute distribution function by Eddington inversion and initialize relevant splines
 
gsl_spline * initSplineInScaledEnergy (const std::vector< double > &initE, const std::vector< double > &initVal)
 initialize spline in scaled energy variable mu=-log(1/(1/Phi(0)-1/E)). More...
 
void initSplineDifCoefE ()
 
void initSplineDifCoefV ()
 

Private Attributes

MODELSTATUS errorcode
 if not MS_OK, signal an error in the constructor
 
double potentialcenter
 central potential not accounting for SMBH
 
double totalmass
 total mass at infinity
 
double masscenter
 BH mass if it exists;.
 
double totalmasscenter
 total central mass (equal to masscenter if no external potential is used, or taken from the external potential), used for scaling the energy
 
double rin
 
double rout
 inner and outer grid radii
 
double gammain
 
double gammaout
 power-law slopes of density at r->0 and r->infinity
 
double coefmrout
 extrapolation for M(r) at large radii: M=totalmass-coefmrout*r^(3-gammaout)
 
double potentialin
 total potential(including BH) at r=rin
 
bool densconst
 precautions for constant-density cores:
 
double alpha
 
double rho0
 
double coefdrhodr
 
double coefdrhodphi
 coefficients for extrapolation to r->0 for constant-density core
 
std::vector< std::vector
< double > > 
splineJ1
 
std::vector< std::vector
< double > > 
splineJ3
 
gsl_spline * spl_loginmass
 
gsl_spline * spl_logrmass
 
gsl_spline * spl_logpotential
 
gsl_spline * spl_logdenspot
 
gsl_spline * spl_logdistrf
 
gsl_spline * spl_logphasevol
 
gsl_spline * spl_logtrad
 
gsl_spline * spl_logrcirc
 
gsl_spline * spl_logfint
 
gsl_spline * spl_logfgint
 
gsl_spline * spl_logfhint
 

Detailed Description

class defining a spherical mass model.

These models are constructed from an array of { r, m(r) } pairs; the dependence of mass, potential, distribution function, and other quantities on radius or energy is pre-computed in the constructor and stored as a spline function of (logarithmically scaled) radius or energy, with a careful extrapolation beyond the definition region of the spline. The isotropic distribution function is computed by the Eddington inversion formula. These models may be used independently of the main SMILE machinery, or for generating initial conditions for the orbit library in COrbitLibrary::prepareInitialConditions

Member Enumeration Documentation

modes of operation for fitMassModel

Enumerator
FT_MASS 

use locations of particles to compute M(r)

FT_DF_ENERGY 

use energies of particles to compute f(E)

FT_DF_PHASEVOL 

use energies of particles to compute f(h) as function of phase volume

internal status or error code of the model

Enumerator
MS_OK 

everything is fine

MS_ARRAY_NOT_EQ 

input arrays are of different lengths

MS_ARRAY_SHORT 

input arrays are too short (less than 5 radial points)

MS_RADII_NOT_ASC 

radii are not in ascending order

MS_MASS_NOT_ASC 

mass is not increasing with radius

MS_DF_UNDEFINED 

could not initialize Eddington distribution function at all energies

MS_DF_NEGATIVE 

distribution function is negative somewhere (not a critical error)

MS_DIFCOEF_DIVERGE 

initialization of diffusion coefficients interpolation table failed

Constructor & Destructor Documentation

smile::CMassModel::CMassModel ( const std::vector< double > &  initRadii,
const std::vector< double > &  initMass,
const CMassModel model_pot = NULL 
)

create the model from an array of { r, m(r) } pairs and pre-computes all necessary quantities.

If another mass model in model_pot is supplied, it is used for the "external" potential, while mass and density profile are taken from the input array. Consequently, quantities such as distribution function and velocity dispersion reflect the "external" potential but "internal" density.

Returns
error code (MS_OK if no error)

Member Function Documentation

void smile::CMassModel::DifCoefEL ( double  E,
double *  DEE = NULL,
double *  DE = NULL,
double *  DLL = NULL,
double *  DRRoverR = NULL,
bool  exact = false 
) const

compute averaged diffusion coefs depending on E.

These are averaged over entire phase volume at given energy. The coefficients are:

Parameters
[out]DEE= <Delta E^2>,
[out]DE= <Delta e>=""> – diffusion and drift coefs in energy;
[out]DLL= <Delta L^2> – diff.coef. in angular momentum;
[out]DRRoverR= <Delta R^2>/(2*R) – limiting slope of diff.coef. in R=L^2/Lcirc^2 at R->0, used for loss-cone problems. If a pointer is NULL, this coef is not computed.
[in]Eis the energy at which the coefficients are computed.
[in]exact– if true, compute using exact quadrature, otherwise using spline interpolation (much faster).
bool smile::CMassModel::initDistrFuncEnergy ( const std::vector< double > &  initE,
const std::vector< double > &  initF 
)

initialize spline for externally supplied distribution function in energy.

return true on success, false if input arrays are unequal or initE is not sorted in increasing order, or values in initE are above zero or below minimum of potential, or DF supplied is non-positive

gsl_spline * smile::CMassModel::initSplineInScaledEnergy ( const std::vector< double > &  initE,
const std::vector< double > &  initVal 
)
private

initialize spline in scaled energy variable mu=-log(1/(1/Phi(0)-1/E)).

input arrays must be either of equal size or the second array be longer by two elements, in which case a clamped spline will be created (with given endpoint derivatives provided in the last two elements of the initVal array); derivatives are with respect to scaled energy.

void smile::CMassModel::initSplinePhi ( const CMassModel model_pot)
private

compute potential and initialize Phi(r) spline, or use external potential if model_pot!=NULL

!!! doesn't work properly for external potential!!!?


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