LORENE
Lorene::Param_elliptic Class Reference

This class contains the parameters needed to call the general elliptic solver. More...

#include <param_elliptic.h>

Public Member Functions

 Param_elliptic (const Scalar &)
 Standard constructor from a Scalar. More...
 
 ~Param_elliptic ()
 Destructor. More...
 
const Map_radialget_mp () const
 Returns the mapping. More...
 
double get_alpha (int) const
 
double get_beta (int) const
 
int get_type (int) const
 
void set_helmholtz_minus (int zone, double mas, Scalar &so)
 Set the operator to $\left(\Delta - m^2\right)$ in one domain (not in the nucleus). More...
 
void set_poisson_pseudo_1d (Scalar &so)
 Set the operator to $\left(\Delta - 2 \partial_r / r\right)$ everywhere but in the compactified domain. More...
 
void set_helmholtz_minus_pseudo_1d (int zone, double mas, Scalar &so)
 Set the operator to $\left(\Delta - 2 \partial_r / r - m^2\right)$ in one domain. More...
 
void set_helmholtz_plus (int zone, double mas, Scalar &so)
 Set the operator to $\left(\Delta + m^2\right)$ in one domain (only in the shells). More...
 
void set_poisson_2d (const Scalar &, bool indic=false)
 Set everything to do a 2d-Poisson, with or without l=0 (not put by default...) More...
 
void set_helmholtz_minus_2d (int zone, double mas, const Scalar &)
 Set the 2D Helmholtz operator (with minus sign) More...
 
void set_sec_order_r2 (int zone, double a, double b, double c)
 Set the operator to $a r^2 \partial^2/\partial r^2 + b r \partial /\partial r + c$ in one domain (only in the shells). More...
 
void set_sec_order (int zone, double a, double b, double c)
 Set the operator to $a \partial^2/\partial r^2 + b \partial /\partial r + c$ in one domain (only in the shells). More...
 
void set_ope_vorton (int zone, Scalar &so)
 Set the operator to $\Delta - 2\partial /\partial r$ in one domain (not implemented in the nucleus). More...
 
void set_poisson_vect_r (int zone, bool only_l_zero=false)
 Sets the operator to $\Delta + \frac{2}{r} \frac{\partial}{\partial r} + \frac{2 - l(l+1)}{r^2} $ in all domains, for $ l \not= 0 $; and to $\frac{\partial^2}{\partial r^2} + \frac{2}{r} \frac{\partial}{\partial r} - \frac{2}{r^2} $ in all domains otherwise. More...
 
void set_poisson_vect_eta (int zone)
 Sets the operator to be a regular elliptic operator to solve for the $\eta $ component of the vector Poisson equation. More...
 
void set_poisson_tens_rr (int zone)
 Sets the operator to $\Delta + \frac{4}{r} \frac{\partial}{\partial r} + \frac{6 - l(l+1)}{r^2} $ in all domains, for $l \geq 2$ only. More...
 
void inc_l_quant (int zone)
 Increases the quantum number l in the domain zone . More...
 
void set_variable_F (const Scalar &)
 Changes the variable function F. More...
 
void set_variable_G (const Scalar &)
 Changes the variable function G. More...
 
double F_plus (int zone, int k, int j) const
 Returns the value of F, for a given angular point, at the outer boundary of the domain zone ;. More...
 
double F_minus (int zone, int k, int j) const
 Returns the value of F, for a given angular point, at the inner boundary of the domain zone ;. More...
 
double dF_plus (int zone, int k, int j) const
 Returns the value of the radial derivative of F, for a given angular point, at the outer boundary of the domain zone ;. More...
 
double dF_minus (int zone, int k, int j) const
 Returns the value of the radial derivative of F, for a given angular point, at the inner boundary of the domain zone ;. More...
 
double G_plus (int zone) const
 Returns the value of G, for a given angular point, at the outer boundary of the domain zone ;. More...
 
double G_minus (int zone) const
 Returns the value of G, for a given angular point, at the inner boundary of the domain zone ;. More...
 
double dG_plus (int zone) const
 Returns the value of the radial derivative of G, for a given angular point, at the outer boundary of the domain zone ;. More...
 
double dG_minus (int zone) const
 Returns the value of the radial derivative of G, for a given angular point, at the inner boundary of the domain zone ;. More...
 

Protected Attributes

int type_map
 Type of mapping either MAP_AFF or MAP_LOG. More...
 
const Map_afmp_af
 The mapping, if affine. More...
 
const Map_logmp_log
 The mapping if log type. More...
 
Ope_elementary ** operateurs
 Array on the elementary operators. More...
 
Scalar var_F
 Additive variable change function. More...
 
Scalar var_G
 Multiplicative variable change that must be sphericaly symetric ! More...
 
Itbl done_F
 Stores what has been computed for F. More...
 
Itbl done_G
 Stores what has been computed for G. More...
 
Tbl val_F_plus
 Values of F at the outer boundaries of the various domains. More...
 
Tbl val_F_minus
 Values of F at the inner boundaries of the various domains. More...
 
Tbl val_dF_plus
 Values of the derivative of F at the outer boundaries of the various domains. More...
 
Tbl val_dF_minus
 Values of the derivative of F at the inner boundaries of the various domains. More...
 
Tbl val_G_plus
 Values of G at the outer boundaries of the various domains. More...
 
Tbl val_G_minus
 Values of G at the inner boundaries of the various domains. More...
 
Tbl val_dG_plus
 Values of the derivative of G at the outer boundaries of the various domains. More...
 
Tbl val_dG_minus
 Values of the derivative of G at the inner boundaries of the various domains. More...
 

Private Member Functions

void compute_val_F (int, int, int) const
 Computes the various values of F. More...
 
void compute_val_G (int) const
 Computes the various values of G. More...
 

Friends

Mtbl_cf elliptic_solver (const Param_elliptic &, const Mtbl_cf &)
 
Mtbl_cf elliptic_solver_boundary (const Param_elliptic &, const Mtbl_cf &, const Mtbl_cf &bound, double fact_dir, double fact_neu)
 
Mtbl_cf elliptic_solver_no_zec (const Param_elliptic &, const Mtbl_cf &, double)
 
Mtbl_cf elliptic_solver_only_zec (const Param_elliptic &, const Mtbl_cf &, double)
 
Mtbl_cf elliptic_solver_sin_zec (const Param_elliptic &, const Mtbl_cf &, double *, double *)
 
Mtbl_cf elliptic_solver_fixe_der_zero (double, const Param_elliptic &, const Mtbl_cf &)
 
void Map_af::sol_elliptic (Param_elliptic &, const Scalar &, Scalar &) const
 
void Map_af::sol_elliptic_boundary (Param_elliptic &, const Scalar &, Scalar &, const Mtbl_cf &, double, double) const
 
void Map_af::sol_elliptic_boundary (Param_elliptic &, const Scalar &, Scalar &, const Scalar &, double, double) const
 
void Map_af::sol_elliptic_no_zec (Param_elliptic &, const Scalar &, Scalar &, double) const
 
void Map_af::sol_elliptic_only_zec (Param_elliptic &, const Scalar &, Scalar &, double) const
 
void Map_af::sol_elliptic_sin_zec (Param_elliptic &, const Scalar &, Scalar &, double *, double *) const
 
void Map_af::sol_elliptic_fixe_der_zero (double, Param_elliptic &, const Scalar &, Scalar &) const
 
void Map_af::sol_elliptic_2d (Param_elliptic &, const Scalar &, Scalar &) const
 
void Map_af::sol_elliptic_pseudo_1d (Param_elliptic &, const Scalar &, Scalar &) const
 
void Map_log::sol_elliptic (Param_elliptic &, const Scalar &, Scalar &) const
 
void Map_log::sol_elliptic_boundary (Param_elliptic &, const Scalar &, Scalar &, const Mtbl_cf &, double, double) const
 
void Map_log::sol_elliptic_boundary (Param_elliptic &, const Scalar &, Scalar &, const Scalar &, double, double) const
 
void Map_log::sol_elliptic_no_zec (Param_elliptic &, const Scalar &, Scalar &, double) const
 

Detailed Description

This class contains the parameters needed to call the general elliptic solver.

For every domain and every spherical harmonics, it contains the appropriate operator of type Ope_elementary and the appropriate variable given by a Change_var .()

This class is only defined on an affine mapping Map_af or a logarithmic one Map_log

Definition at line 135 of file param_elliptic.h.

Constructor & Destructor Documentation

◆ Param_elliptic()

Lorene::Param_elliptic::Param_elliptic ( const Scalar so)

◆ ~Param_elliptic()

Lorene::Param_elliptic::~Param_elliptic ( )

Member Function Documentation

◆ compute_val_F()

void Lorene::Param_elliptic::compute_val_F ( int  zone,
int  k,
int  j 
) const
private

◆ compute_val_G()

void Lorene::Param_elliptic::compute_val_G ( int  zone) const
private

◆ dF_minus()

double Lorene::Param_elliptic::dF_minus ( int  zone,
int  k,
int  j 
) const

Returns the value of the radial derivative of F, for a given angular point, at the inner boundary of the domain zone ;.

Definition at line 82 of file param_elliptic_val_lim.C.

References compute_val_F(), done_F, and val_dF_minus.

◆ dF_plus()

double Lorene::Param_elliptic::dF_plus ( int  zone,
int  k,
int  j 
) const

Returns the value of the radial derivative of F, for a given angular point, at the outer boundary of the domain zone ;.

Definition at line 74 of file param_elliptic_val_lim.C.

References compute_val_F(), done_F, and val_dF_plus.

◆ dG_minus()

double Lorene::Param_elliptic::dG_minus ( int  zone) const

Returns the value of the radial derivative of G, for a given angular point, at the inner boundary of the domain zone ;.

Definition at line 115 of file param_elliptic_val_lim.C.

References compute_val_G(), done_G, and val_dG_minus.

◆ dG_plus()

double Lorene::Param_elliptic::dG_plus ( int  zone) const

Returns the value of the radial derivative of G, for a given angular point, at the outer boundary of the domain zone ;.

Definition at line 107 of file param_elliptic_val_lim.C.

References compute_val_G(), done_G, and val_dG_plus.

◆ F_minus()

double Lorene::Param_elliptic::F_minus ( int  zone,
int  k,
int  j 
) const

Returns the value of F, for a given angular point, at the inner boundary of the domain zone ;.

Definition at line 66 of file param_elliptic_val_lim.C.

References compute_val_F(), done_F, and val_F_minus.

◆ F_plus()

double Lorene::Param_elliptic::F_plus ( int  zone,
int  k,
int  j 
) const

Returns the value of F, for a given angular point, at the outer boundary of the domain zone ;.

Definition at line 58 of file param_elliptic_val_lim.C.

References compute_val_F(), done_F, and val_F_plus.

◆ G_minus()

double Lorene::Param_elliptic::G_minus ( int  zone) const

Returns the value of G, for a given angular point, at the inner boundary of the domain zone ;.

Definition at line 99 of file param_elliptic_val_lim.C.

References compute_val_G(), done_G, and val_G_minus.

◆ G_plus()

double Lorene::Param_elliptic::G_plus ( int  zone) const

Returns the value of G, for a given angular point, at the outer boundary of the domain zone ;.

Definition at line 91 of file param_elliptic_val_lim.C.

References compute_val_G(), done_G, and val_G_plus.

◆ get_mp()

const Map_radial & Lorene::Param_elliptic::get_mp ( ) const

Returns the mapping.

Definition at line 124 of file param_elliptic.C.

References mp_af, mp_log, and type_map.

◆ inc_l_quant()

void Lorene::Param_elliptic::inc_l_quant ( int  zone)

Increases the quantum number l in the domain zone .

Definition at line 329 of file param_elliptic.C.

References Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Ope_elementary::inc_l_quant(), and operateurs.

◆ set_helmholtz_minus()

void Lorene::Param_elliptic::set_helmholtz_minus ( int  zone,
double  mas,
Scalar so 
)

Set the operator to $\left(\Delta - m^2\right)$ in one domain (not in the nucleus).

Parameters
zone[input] : the domain.
mas[input] : the masse $m$.
so[input] : the source (used only to get the right basis).

Definition at line 349 of file param_elliptic.C.

References Lorene::Valeur::base, Lorene::Valeur::coef(), Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_va(), Lorene::Map_log::get_type(), Lorene::Base_val::give_quant_numbers(), mp_log, operateurs, Lorene::Scalar::set_spectral_va(), type_map, and Lorene::Valeur::ylm().

◆ set_helmholtz_minus_2d()

void Lorene::Param_elliptic::set_helmholtz_minus_2d ( int  zone,
double  mas,
const Scalar source 
)

◆ set_helmholtz_minus_pseudo_1d()

void Lorene::Param_elliptic::set_helmholtz_minus_pseudo_1d ( int  zone,
double  mas,
Scalar so 
)

Set the operator to $\left(\Delta - 2 \partial_r / r - m^2\right)$ in one domain.

Parameters
zone[input] : the domain.
mas[input] : the masse $m$.
so[input] : the source (used only to get the right basis).

Definition at line 101 of file param_elliptic_pseudo_1d.C.

References Lorene::Valeur::base, Lorene::Scalar::check_dzpuis(), Lorene::Scalar::get_dzpuis(), Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_va(), Lorene::Base_val::give_quant_numbers(), operateurs, and type_map.

◆ set_helmholtz_plus()

void Lorene::Param_elliptic::set_helmholtz_plus ( int  zone,
double  mas,
Scalar so 
)

◆ set_ope_vorton()

void Lorene::Param_elliptic::set_ope_vorton ( int  zone,
Scalar so 
)

◆ set_poisson_2d()

void Lorene::Param_elliptic::set_poisson_2d ( const Scalar source,
bool  indic = false 
)

◆ set_poisson_pseudo_1d()

void Lorene::Param_elliptic::set_poisson_pseudo_1d ( Scalar so)

Set the operator to $\left(\Delta - 2 \partial_r / r\right)$ everywhere but in the compactified domain.

Parameters
so[input] : the source (used only to get the right basis).

Definition at line 65 of file param_elliptic_pseudo_1d.C.

References Lorene::Valeur::base, Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_va(), Lorene::Base_val::give_quant_numbers(), operateurs, and type_map.

◆ set_poisson_tens_rr()

void Lorene::Param_elliptic::set_poisson_tens_rr ( int  zone)

◆ set_poisson_vect_eta()

void Lorene::Param_elliptic::set_poisson_vect_eta ( int  zone)

Sets the operator to be a regular elliptic operator to solve for the $\eta $ component of the vector Poisson equation.

The operator is $\frac{\partial^2}{\partial r^2} + \frac{2}{r} \frac{\partial}{\partial r} - \frac{l(l-1)}{r^2} $ (Poisson with the decrease of l by one unit) in all domains but the CED, for $ l \not= 0 $; it is not defined for l = 0. In the CED, the operator is also the Laplace one, but with l increased by one unit: $\frac{\partial^2}{\partial r^2} + \frac{2}{r} \frac{\partial}{\partial r} - \frac{(l+1)(l+2)}{r^2} $. This is intended to solve the equation for $ \eta $ arising in the decomposition of the vector Poisson equation.

Parameters
zone[input] : the domain.

Definition at line 539 of file param_elliptic.C.

References Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_nzone(), and Lorene::Mg3d::get_type_r().

◆ set_poisson_vect_r()

void Lorene::Param_elliptic::set_poisson_vect_r ( int  zone,
bool  only_l_zero = false 
)

Sets the operator to $\Delta + \frac{2}{r} \frac{\partial}{\partial r} + \frac{2 - l(l+1)}{r^2} $ in all domains, for $ l \not= 0 $; and to $\frac{\partial^2}{\partial r^2} + \frac{2}{r} \frac{\partial}{\partial r} - \frac{2}{r^2} $ in all domains otherwise.

Parameters
zone[input] : the domain.
only_l_zero[input] : the operator in built only for l=0

Definition at line 493 of file param_elliptic.C.

References Lorene::Ope_elementary::get_base_r(), Lorene::Ope_poisson::get_dzpuis(), Lorene::Ope_poisson::get_lquant(), Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), operateurs, and type_map.

◆ set_sec_order()

void Lorene::Param_elliptic::set_sec_order ( int  zone,
double  a,
double  b,
double  c 
)

Set the operator to $a \partial^2/\partial r^2 + b \partial /\partial r + c$ in one domain (only in the shells).

Parameters
zone[input] : the domain.
a[input] : the parameter $a$.
b[input] : the parameter $b$.
c[input] : the parameter $c$.

Definition at line 675 of file param_elliptic.C.

References Lorene::Ope_elementary::get_base_r(), Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Map_log::get_type(), mp_log, operateurs, R_CHEBI, and type_map.

◆ set_sec_order_r2()

void Lorene::Param_elliptic::set_sec_order_r2 ( int  zone,
double  a,
double  b,
double  c 
)

Set the operator to $a r^2 \partial^2/\partial r^2 + b r \partial /\partial r + c$ in one domain (only in the shells).

Parameters
zone[input] : the domain.
a[input] : the parameter $a$.
b[input] : the parameter $b$.
c[input] : the parameter $c$.

Definition at line 608 of file param_elliptic.C.

References Lorene::Ope_elementary::get_base_r(), Lorene::Map::get_mg(), get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Map_log::get_type(), mp_log, operateurs, R_CHEBI, and type_map.

◆ set_variable_F()

void Lorene::Param_elliptic::set_variable_F ( const Scalar so)

Changes the variable function F.

Definition at line 780 of file param_elliptic.C.

References Lorene::Scalar::get_etat().

◆ set_variable_G()

void Lorene::Param_elliptic::set_variable_G ( const Scalar so)

Changes the variable function G.

Definition at line 789 of file param_elliptic.C.

References Lorene::Scalar::get_etat().

Member Data Documentation

◆ done_F

Itbl Lorene::Param_elliptic::done_F
mutableprotected

Stores what has been computed for F.

Definition at line 147 of file param_elliptic.h.

◆ done_G

Itbl Lorene::Param_elliptic::done_G
mutableprotected

Stores what has been computed for G.

Definition at line 148 of file param_elliptic.h.

◆ mp_af

const Map_af* Lorene::Param_elliptic::mp_af
protected

The mapping, if affine.

Definition at line 139 of file param_elliptic.h.

◆ mp_log

const Map_log* Lorene::Param_elliptic::mp_log
protected

The mapping if log type.

Definition at line 140 of file param_elliptic.h.

◆ operateurs

Ope_elementary** Lorene::Param_elliptic::operateurs
protected

Array on the elementary operators.

Definition at line 142 of file param_elliptic.h.

◆ type_map

int Lorene::Param_elliptic::type_map
protected

Type of mapping either MAP_AFF or MAP_LOG.

Definition at line 138 of file param_elliptic.h.

◆ val_dF_minus

Tbl Lorene::Param_elliptic::val_dF_minus
mutableprotected

Values of the derivative of F at the inner boundaries of the various domains.

Definition at line 152 of file param_elliptic.h.

◆ val_dF_plus

Tbl Lorene::Param_elliptic::val_dF_plus
mutableprotected

Values of the derivative of F at the outer boundaries of the various domains.

Definition at line 151 of file param_elliptic.h.

◆ val_dG_minus

Tbl Lorene::Param_elliptic::val_dG_minus
mutableprotected

Values of the derivative of G at the inner boundaries of the various domains.

Definition at line 156 of file param_elliptic.h.

◆ val_dG_plus

Tbl Lorene::Param_elliptic::val_dG_plus
mutableprotected

Values of the derivative of G at the outer boundaries of the various domains.

Definition at line 155 of file param_elliptic.h.

◆ val_F_minus

Tbl Lorene::Param_elliptic::val_F_minus
mutableprotected

Values of F at the inner boundaries of the various domains.

Definition at line 150 of file param_elliptic.h.

◆ val_F_plus

Tbl Lorene::Param_elliptic::val_F_plus
mutableprotected

Values of F at the outer boundaries of the various domains.

Definition at line 149 of file param_elliptic.h.

◆ val_G_minus

Tbl Lorene::Param_elliptic::val_G_minus
mutableprotected

Values of G at the inner boundaries of the various domains.

Definition at line 154 of file param_elliptic.h.

◆ val_G_plus

Tbl Lorene::Param_elliptic::val_G_plus
mutableprotected

Values of G at the outer boundaries of the various domains.

Definition at line 153 of file param_elliptic.h.

◆ var_F

Scalar Lorene::Param_elliptic::var_F
protected

Additive variable change function.

Definition at line 144 of file param_elliptic.h.

◆ var_G

Scalar Lorene::Param_elliptic::var_G
protected

Multiplicative variable change that must be sphericaly symetric !

Definition at line 145 of file param_elliptic.h.


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