SMILE  v2.5
Schwarzschild Modelling Interactive expLoratory Environment
Public Member Functions | List of all members
smile::CPotential Class Referenceabstract

base class for potential-density model including potential derivatives (=forces) More...

#include <potential.h>

Inheritance diagram for smile::CPotential:
Inheritance graph
[legend]
Collaboration diagram for smile::CPotential:
Collaboration graph
[legend]

Public Member Functions

virtual CPotentialclone () const =0
 Return a pointer to a copy of this instance of potential. More...
 
virtual double Phi (double X, double Y, double Z, double t=0) const =0
 Return potential at a given spatial point (possibly a time-varying one). More...
 
virtual void Force (const double xyz[N_DIM], const double t, double *force, double *forceDeriv=NULL) const =0
 Compute forces and, optionally, force derivatives at a given point. More...
 
- Public Member Functions inherited from smile::CDensity
virtual POTENTIALTYPE PotentialType () const =0
 enumerable potential type
 
virtual const char * PotentialName () const =0
 string representation of potential type
 
virtual SYMMETRYTYPE symmetry () const =0
 returns symmetry type of this potential
 
virtual double Rho (double X, double Y, double Z, double t=0) const =0
 returns density at given coordinates, this should obviously be overriden in derivative classes
 
virtual double Mass (const double r) const
 returns mass inside given radius (approximately! not necessary to integrate density over sphere, just a rough estimate used e.g. in choosing radial nodes of Schwarzschild grid)
 
virtual double totalMass () const
 returns estimated M(r=infinity) or -1 if mass is infinite
 
double getRadiusByMass (const double m) const
 solves for Mass(r)=m
 
void getRadiiByMass (const vectord &masses, vectord *radii) const
 solves for Mass(r)=m for an array of sorted values of m (more efficient than doing it one-by-one)
 
bool checkMassMonotonic () const
 safety measure: check (roughly) that mass is increasing with radius
 
bool checkDensityNonzero () const
 another safety measure: check that density doesn't drop to zero along any of three axes (important to assess spherical-harmonic approximation quality)
 
virtual double getGamma () const
 returns inner density slope estimate (only used in BSE potential expansion for the automatic selection of shape parameter Alpha)
 

Additional Inherited Members

- Public Types inherited from smile::CDensity
enum  POTENTIALTYPE {
  PT_UNKNOWN, PT_DIRECT, PT_COMPOSITE, PT_COEFS,
  PT_NB, PT_BSE, PT_BSECOMPACT, PT_SPLINE,
  PT_CYLSPLINE, PT_LOG, PT_HARMONIC, PT_SCALEFREE,
  PT_SCALEFREESH, PT_SPHERICAL, PT_DEHNEN, PT_MIYAMOTONAGAI,
  PT_FERRERS, PT_PLUMMER, PT_ISOCHRONE, PT_PERFECTELLIPSOID,
  PT_NFW, PT_SERSIC, PT_EXPDISK, PT_ELLIPSOIDAL,
  PT_MGE
}
 list of all existing types of density or density/potential models, each of them implemented in its own class More...
 
enum  SYMMETRYTYPE {
  ST_NONE = 0, ST_REFLECTION = 1, ST_PLANESYM = 2, ST_ZROTSYM = 4,
  ST_SPHSYM = 8, ST_TRIAXIAL = ST_REFLECTION | ST_PLANESYM, ST_AXISYMMETRIC = ST_TRIAXIAL | ST_ZROTSYM, ST_SPHERICAL = ST_AXISYMMETRIC | ST_SPHSYM,
  ST_DEFAULT = ST_TRIAXIAL
}
 Type of symmetry. More...
 

Detailed Description

base class for potential-density model including potential derivatives (=forces)

Member Function Documentation

virtual CPotential* smile::CPotential::clone ( ) const
pure virtual
virtual void smile::CPotential::Force ( const double  xyz[N_DIM],
const double  t,
double *  force,
double *  forceDeriv = NULL 
) const
pure virtual

Compute forces and, optionally, force derivatives at a given point.

Parameters
[in]xyz- coordinates of the point to compute forces (array of 3 numbers)
[in]t- time to compute forces (matters only if potential is time-dependent)
[out]force- computed values of -d Phi/d x_i; output array must exist and contain N_DIM values.
[out]forceDeriv- if not NULL, then also compute the second derivatives of potential (which is a symmetric matrix of size N_DIM^2, thus the array must be of size N_DIM*(N_DIM+1)/2 ): first 3 values contain $ -d^2 \Phi/d x_i^2 $, second three contain mixed derivatives $ -d^2 \Phi/d x_1 d x_2, -d^2 \Phi/d x_2 d x_3, -d^2 \Phi/d x_3 d x_1 $.

Implemented in smile::CPotentialNB, smile::CPotentialScaleFree, smile::CPotentialCylSpline, smile::CPotentialSH, smile::CPotentialDehnen, smile::CPotentialFerrers, smile::CPotentialMiyamotoNagai, smile::CPotentialSpherical, smile::CPotentialHarmonic, smile::CPotentialLog, smile::CPotentialComposite, and smile::CPotentialAxisymmetrized.

virtual double smile::CPotential::Phi ( double  X,
double  Y,
double  Z,
double  t = 0 
) const
pure virtual

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