Tensor field of valence 0 (or component of a tensorial field). More...
#include <scalar.h>
Public Member Functions | |
Scalar (const Map &mpi) | |
Constructor from mapping. | |
Scalar (const Tensor &a) | |
Constructor from a Tensor (must be of valence 0). | |
Scalar (const Scalar &a) | |
Copy constructor. | |
Scalar (const Cmp &a) | |
Constructor by conversion of a Cmp. | |
Scalar (const Map &, const Mg3d &, FILE *) | |
Constructor from a file (see sauve(FILE*) ). | |
virtual | ~Scalar () |
Destructor. | |
virtual void | set_etat_nondef () |
Sets the logical state to ETATNONDEF (undefined). | |
virtual void | set_etat_zero () |
Sets the logical state to ETATZERO (zero). | |
virtual void | set_etat_qcq () |
Sets the logical state to ETATQCQ (ordinary state). | |
void | set_etat_one () |
Sets the logical state to ETATUN (one). | |
virtual void | allocate_all () |
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elements, down to the double arrays of the Tbl s. | |
void | annule_hard () |
Sets the Scalar to zero in a hard way. | |
int | get_etat () const |
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary). | |
int | get_dzpuis () const |
Returns dzpuis . | |
bool | dz_nonzero () const |
Returns true if the last domain is compactified and *this is not zero in this domain. | |
bool | check_dzpuis (int dzi) const |
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is not equal to dzi , otherwise return true. | |
void | operator= (const Scalar &a) |
Assignment to another Scalar defined on the same mapping. | |
virtual void | operator= (const Tensor &a) |
Assignment to a Tensor (of valence 0). | |
void | operator= (const Cmp &a) |
Assignment to a Cmp . | |
void | operator= (const Valeur &a) |
Assignment to a Valeur . | |
void | operator= (const Mtbl &a) |
Assignment to a Mtbl . | |
void | operator= (double) |
Assignment to a double . | |
void | operator= (int) |
Assignment to an int . | |
operator Cmp () const | |
Conversion to a Cmp . | |
const Valeur & | get_spectral_va () const |
Returns va (read only version). | |
Valeur & | set_spectral_va () |
Returns va (read/write version). | |
Tbl & | set_domain (int l) |
Read/write of the value in a given domain. | |
const Tbl & | domain (int l) const |
Read-only of the value in a given domain. | |
double | val_grid_point (int l, int k, int j, int i) const |
Returns the value of the field at a specified grid point. | |
double | val_point (double r, double theta, double phi) const |
Computes the value of the field at an arbitrary point , by means of the spectral expansion. | |
double & | set_grid_point (int l, int k, int j, int i) |
Setting the value of the field at a given grid point. | |
virtual void | annule (int l_min, int l_max) |
Sets the Scalar to zero in several domains. | |
void | set_inner_boundary (int l, double x) |
Sets the value of the Scalar at the inner boundary of a given domain. | |
void | set_outer_boundary (int l, double x) |
Sets the value of the Scalar at the outer boundary of a given domain. | |
Tbl | multipole_spectrum () const |
Gives the spectrum in terms of multipolar modes l . | |
Tbl | tbl_out_bound (int l_dom, bool leave_ylm=false) |
Returns the Tbl containing the values of angular coefficients at the outer boundary. | |
Tbl | tbl_in_bound (int n, bool leave_ylm=false) |
Returns the Tbl containing the values of angular coefficients at the inner boundary. | |
Scalar | scalar_out_bound (int n, bool leave_ylm=false) |
Returns the Scalar containing the values of angular coefficients at the outer boundary. | |
const Scalar & | dsdr () const |
Returns of *this . | |
const Scalar & | srdsdt () const |
Returns of *this . | |
const Scalar & | srstdsdp () const |
Returns of *this . | |
const Scalar & | dsdt () const |
Returns of *this . | |
const Scalar & | dsdradial () const |
Returns of *this if the mapping is affine (class Map_af ) and of *this if the mapping is logarithmic (class Map_log ). | |
const Scalar & | dsdrho () const |
Returns of *this . | |
const Scalar & | stdsdp () const |
Returns of *this . | |
const Scalar & | dsdx () const |
Returns of *this , where . | |
const Scalar & | dsdy () const |
Returns of *this , where . | |
const Scalar & | dsdz () const |
Returns of *this , where . | |
const Scalar & | deriv (int i) const |
Returns of *this , where . | |
const Vector & | derive_cov (const Metric &gam) const |
Returns the gradient (1-form = covariant vector) of *this . | |
const Vector & | derive_con (const Metric &gam) const |
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of the gradient (cf. | |
Scalar | derive_lie (const Vector &v) const |
Computes the derivative of this along a vector field v . | |
const Scalar & | laplacian (int ced_mult_r=4) const |
Returns the Laplacian of *this . | |
const Scalar & | lapang () const |
Returns the angular Laplacian of *this , where . | |
void | div_r () |
Division by r everywhere; dzpuis is not changed. | |
void | div_r_dzpuis (int ced_mult_r) |
Division by r everywhere but with the output flag dzpuis set to ced_mult_r . | |
void | div_r_ced () |
Division by r in the compactified external domain (CED), the dzpuis flag is not changed. | |
void | mult_r () |
Multiplication by r everywhere; dzpuis is not changed. | |
void | mult_r_dzpuis (int ced_mult_r) |
Multiplication by r everywhere but with the output flag dzpuis set to ced_mult_r . | |
void | mult_r_ced () |
Multiplication by r in the compactified external domain (CED), the dzpuis flag is not changed. | |
void | mult_rsint () |
Multiplication by everywhere; dzpuis is not changed. | |
void | mult_rsint_dzpuis (int ced_mult_r) |
Multiplication by but with the output flag dzpuis set to ced_mult_r . | |
void | div_rsint () |
Division by everywhere; dzpuis is not changed. | |
void | div_rsint_dzpuis (int ced_mult_r) |
Division by but with the output flag dzpuis set to ced_mult_r . | |
void | mult_cost () |
Multiplication by . | |
void | div_cost () |
Division by . | |
void | mult_sint () |
Multiplication by . | |
void | div_sint () |
Division by . | |
void | div_tant () |
Division by . | |
Scalar | primr (bool null_infty=true) const |
Computes the radial primitive which vanishes for . | |
double | integrale () const |
Computes the integral over all space of *this . | |
const Tbl & | integrale_domains () const |
Computes the integral in each domain of *this . | |
virtual void | dec_dzpuis (int dec=1) |
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the compactified external domain (CED). | |
virtual void | inc_dzpuis (int inc=1) |
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the compactified external domain (CED). | |
virtual void | change_triad (const Base_vect &new_triad) |
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly. | |
void | filtre (int n) |
Sets the n lasts coefficients in r to 0 in the external domain. | |
void | filtre_r (int *nn) |
Sets the n lasts coefficients in r to 0 in all domains. | |
void | filtre_r (int n, int nzone) |
Sets the n last coefficients in r to 0 in the domain nzone . | |
virtual void | exponential_filter_r (int lzmin, int lzmax, int p, double alpha=-16.) |
Applies an exponential filter to the spectral coefficients in the radial direction. | |
void | sarra_filter_r (int lzmin, int lzmax, double p, double alpha=-1E-16) |
Applies an exponential filter to the spectral coefficients in the radial direction. | |
void | exp_filter_r_all_domains (Scalar &ss, int p, double alpha=-16.) |
Applies an exponential filter in radial direction in all domains. | |
void | sarra_filter_r_all_domains (double p, double alpha=1E-16) |
Applies an exponential filter in radial direction in all domains for the case where p is a double (see Scalar:sarra_filter_r ). | |
virtual void | exponential_filter_ylm (int lzmin, int lzmax, int p, double alpha=-16.) |
Applies an exponential filter to the spectral coefficients in the angular directions. | |
void | annule_l (int l_min, int l_max, bool ylm_output=false) |
Sets all the multipolar components between l_min and l_max to zero. | |
void | filtre_phi (int n, int zone) |
Sets the n lasts coefficients in to 0 in the domain zone . | |
void | filtre_tp (int nn, int nz1, int nz2) |
Sets the n lasts coefficients in to 0 in the domain nz1 to nz2 when expressed in spherical harmonics. | |
void | fixe_decroissance (int puis) |
Substracts all the components behaving like in the external domain, with n strictly lower than puis , so that *this decreases at least like at infinity. | |
void | smooth_decay (int k, int n) |
Performs a matching of the last non-compactified shell with a decaying function where is the spherical harmonic index and n is some specifiable parameter. | |
void | raccord (int n) |
Performs the matching of the nucleus with respect to the first shell. | |
void | raccord_c1_zec (int puis, int nbre, int lmax) |
Performs the matching of the external domain with respect to the last shell using function like with for each spherical harmonics with . | |
void | raccord_externe (int puis, int nbre, int lmax) |
Matching of the external domain with the outermost shell. | |
void | match_tau (Param &par_bc, Param *par_mat=0x0) |
Method for matching accross domains and imposing boundary condition. | |
void | match_tau_shell (Param &par_bc, Param *par_mat=0x0) |
Method for matching accross domains and imposing boundary condition. | |
void | match_collocation (Param &par_bc, Param *par_mat=0x0) |
Method for matching accross domains and imposing boundary condition. | |
virtual void | sauve (FILE *) const |
Save in a file. | |
virtual void | spectral_display (const char *comment=0x0, double threshold=1.e-7, int precision=4, ostream &ostr=cout) const |
Displays the spectral coefficients and the associated basis functions. | |
void | visu_section (const char section_type, double aa, double umin, double umax, double vmin, double vmax, const char *title=0x0, const char *filename=0x0, bool start_dx=true, int nu=200, int nv=200) const |
3D visualization via a plane section. | |
void | visu_section (const Tbl &plane, double umin, double umax, double vmin, double vmax, const char *title=0x0, const char *filename=0x0, bool start_dx=true, int nu=200, int nv=200) const |
3D visualization via a plane section. | |
void | visu_section_anim (const char section_type, double aa, double umin, double umax, double vmin, double vmax, int jtime, double ttime, int jgraph=1, const char *title=0x0, const char *filename_root=0x0, bool start_dx=false, int nu=200, int nv=200) const |
3D visualization via time evolving plane section (animation). | |
void | visu_box (double xmin, double xmax, double ymin, double ymax, double zmin, double zmax, const char *title0=0x0, const char *filename0=0x0, bool start_dx=true, int nx=40, int ny=40, int nz=40) const |
3D visualization (volume rendering) via OpenDX. | |
void | operator+= (const Scalar &) |
+= Scalar | |
void | operator-= (const Scalar &) |
-= Scalar | |
void | operator*= (const Scalar &) |
*= Scalar | |
virtual void | std_spectral_base () |
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. | |
virtual void | std_spectral_base_odd () |
Sets the spectral bases of the Valeur va to the standard odd ones for a scalar field. | |
void | set_spectral_base (const Base_val &) |
Sets the spectral bases of the Valeur va . | |
const Base_val & | get_spectral_base () const |
Returns the spectral bases of the Valeur va . | |
void | set_dzpuis (int) |
Modifies the dzpuis flag. | |
Valeur ** | asymptot (int n, const int flag=0) const |
Asymptotic expansion at r = infinity. | |
Scalar | poisson () const |
Solves the scalar Poisson equation with *this as a source. | |
void | poisson (Param &par, Scalar &uu) const |
Solves the scalar Poisson equation with *this as a source (version with parameters to control the resolution). | |
Scalar | poisson_tau () const |
Solves the scalar Poisson equation with *this as a source using a real Tau method The source of the equation is represented by the Scalar *this . | |
void | poisson_tau (Param &par, Scalar &uu) const |
Solves the scalar Poisson equation with *this as a source using a real Tau method (version with parameters to control the resolution) The source of the equation is represented by the Scalar *this . | |
Scalar | poisson_dirichlet (const Valeur &limite, int num) const |
Is identicall to Scalar::poisson() . | |
Scalar | poisson_neumann (const Valeur &, int) const |
Idem as Scalar::poisson_dirichlet , the boundary condition being on the radial derivative of the solution. | |
Scalar | poisson_dir_neu (const Valeur &limite, int num, double fact_dir, double fact_neu) const |
Is identicall to Scalar::poisson() . | |
Scalar | poisson_frontiere_double (const Valeur &, const Valeur &, int) const |
Idem as Scalar::poisson_dirichlet , the boundary condition being on both the function and its radial derivative. | |
void | poisson_regular (int k_div, int nzet, double unsgam1, Param &par, Scalar &uu, Scalar &uu_regu, Scalar &uu_div, Tensor &duu_div, Scalar &source_regu, Scalar &source_div) const |
Solves the scalar Poisson equation with *this as a source (version with parameters to control the resolution). | |
Tbl | test_poisson (const Scalar &uu, ostream &ostr, bool detail=false) const |
Checks if a Poisson equation with *this as a source has been correctly solved. | |
Scalar | poisson_angu (double lambda=0) const |
Solves the (generalized) angular Poisson equation with *this as source. | |
Scalar | avance_dalembert (Param &par, const Scalar &fJm1, const Scalar &source) const |
Performs one time-step integration (from ) of the scalar d'Alembert equation with *this being the value of the function f at time J . | |
Scalar | sol_elliptic (Param_elliptic ¶ms) const |
Resolution of a general elliptic equation, putting zero at infinity. | |
Scalar | sol_elliptic_boundary (Param_elliptic ¶ms, const Mtbl_cf &bound, double fact_dir, double fact_neu) const |
Resolution of a general elliptic equation, putting zero at infinity and with inner boundary conditions. | |
Scalar | sol_elliptic_boundary (Param_elliptic ¶ms, const Scalar &bound, double fact_dir, double fact_neu) const |
Resolution of general elliptic equation, with inner boundary conditions as Scalars on mono-domain angulare grids. | |
Scalar | sol_elliptic_2d (Param_elliptic &) const |
Solves the scalar 2-dimensional elliptic equation with *this as a source. | |
Scalar | sol_elliptic_pseudo_1d (Param_elliptic &) const |
Solves a pseudo-1d elliptic equation with *this as a source. | |
Scalar | sol_elliptic_no_zec (Param_elliptic ¶ms, double val=0) const |
Resolution of a general elliptic equation, putting a given value at the outermost shell and not solving in the compactified domain. | |
Scalar | sol_elliptic_only_zec (Param_elliptic ¶ms, double val) const |
Resolution of a general elliptic equation solving in the compactified domain and putting a given value at the inner boundary. | |
Scalar | sol_elliptic_sin_zec (Param_elliptic ¶ms, double *coefs, double *phases) const |
General elliptic solver. | |
Scalar | sol_elliptic_fixe_der_zero (double val, Param_elliptic ¶ms) const |
Resolution of a general elliptic equation fixing the dericative at the origin and relaxing one continuity condition. | |
Scalar | sol_divergence (int n) const |
Resolution of a divergence-like equation. | |
void | import (const Scalar &ci) |
Assignment to another Scalar defined on a different mapping. | |
void | import_symy (const Scalar &ci) |
Assignment to another Scalar defined on a different mapping. | |
void | import_asymy (const Scalar &ci) |
Assignment to another Scalar defined on a different mapping. | |
void | import (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping. | |
void | import_symy (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping. | |
void | import_asymy (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping. | |
void | set_triad (const Base_vect &new_triad) |
Assigns a new vectorial basis (triad) of decomposition. | |
Scalar & | set (const Itbl &ind) |
Returns the value of a component (read/write version). | |
Scalar & | set (int i1, int i2) |
Returns the value of a component for a tensor of valence 2 (read/write version). | |
Scalar & | set (int i1, int i2, int i3) |
Returns the value of a component for a tensor of valence 3 (read/write version). | |
Scalar & | set (int i1, int i2, int i3, int i4) |
Returns the value of a component for a tensor of valence 4 (read/write version). | |
void | annule_domain (int l) |
Sets the Tensor to zero in a given domain. | |
void | annule_extern_cn (int l_0, int deg) |
Performs a smooth (C^n) transition in a given domain to zero. | |
const Tensor & | divergence (const Metric &gam) const |
Computes the divergence of this with respect to some metric . | |
Tensor | up (int ind, const Metric &gam) const |
Computes a new tensor by raising an index of *this . | |
Tensor | down (int ind, const Metric &gam) const |
Computes a new tensor by lowering an index of *this . | |
Tensor | up_down (const Metric &gam) const |
Computes a new tensor by raising or lowering all the indices of *this . | |
Tensor | trace (int ind1, int ind2) const |
Trace on two different type indices. | |
Tensor | trace (int ind1, int ind2, const Metric &gam) const |
Trace with respect to a given metric. | |
Scalar | trace () const |
Trace on two different type indices for a valence 2 tensor. | |
Scalar | trace (const Metric &gam) const |
Trace with respect to a given metric for a valence 2 tensor. | |
virtual int | position (const Itbl &ind) const |
Returns the position in the array cmp of a component given by its indices. | |
virtual Itbl | indices (int pos) const |
Returns the indices of a component given by its position in the array cmp . | |
const Map & | get_mp () const |
Returns the mapping. | |
const Base_vect * | get_triad () const |
Returns the vectorial basis (triad) on which the components are defined. | |
int | get_valence () const |
Returns the valence. | |
int | get_n_comp () const |
Returns the number of stored components. | |
int | get_index_type (int i) const |
Gives the type (covariant or contravariant) of the index number i . | |
Itbl | get_index_type () const |
Returns the types of all the indices. | |
int & | set_index_type (int i) |
Sets the type of the index number i . | |
Itbl & | set_index_type () |
Sets the types of all the indices. | |
const Scalar & | operator() (const Itbl &ind) const |
Returns the value of a component (read-only version). | |
const Scalar & | operator() (int i1, int i2) const |
Returns the value of a component for a tensor of valence 2 (read-only version). | |
const Scalar & | operator() (int i1, int i2, int i3) const |
Returns the value of a component for a tensor of valence 3 (read-only version). | |
const Scalar & | operator() (int i1, int i2, int i3, int i4) const |
Returns the value of a component for a tensor of valence 4 (read-only version). | |
Protected Member Functions | |
void | del_t () |
Logical destructor. | |
virtual void | del_deriv () const |
Logical destructor of the derivatives. | |
void | set_der_0x0 () const |
Sets the pointers for derivatives to 0x0. | |
void | import_gal (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings do not have a particular relative orientation. | |
void | import_align (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Cartesian axis. | |
void | import_anti (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned Cartesian axis (i.e. | |
void | import_align_symy (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Cartesian axis. | |
void | import_anti_symy (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned Cartesian axis (i.e. | |
void | import_align_asymy (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Cartesian axis. | |
void | import_anti_asymy (int nzet, const Scalar &ci) |
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned Cartesian axis (i.e. | |
virtual void | del_derive_met (int) const |
Logical destructor of the derivatives depending on the i-th element of met_depend . | |
void | set_der_met_0x0 (int) const |
Sets all the i-th components of met_depend , p_derive_cov , etc. | |
void | set_dependance (const Metric &) const |
To be used to describe the fact that the derivatives members have been calculated with met . | |
int | get_place_met (const Metric &) const |
Returns the position of the pointer on metre in the array met_depend . | |
void | compute_derive_lie (const Vector &v, Tensor &resu) const |
Computes the Lie derivative of this with respect to some vector field v (protected method; the public interface is method derive_lie ). | |
Protected Attributes | |
int | etat |
The logical state ETATNONDEF (undefined), ETATZERO (null), ETATUN (one), or ETATQCQ (ordinary). | |
int | dzpuis |
Power of r by which the quantity represented by this must be divided in the compactified external domain (CED) in order to get the correct physical values. | |
Valeur | va |
The numerical value of the Scalar . | |
Scalar * | p_dsdr |
Pointer on of *this (0x0 if not up to date). | |
Scalar * | p_srdsdt |
Pointer on of *this (0x0 if not up to date). | |
Scalar * | p_srstdsdp |
Pointer on of *this (0x0 if not up to date). | |
Scalar * | p_dsdt |
Pointer on of *this (0x0 if not up to date). | |
Scalar * | p_stdsdp |
Pointer on of *this (0x0 if not up to date). | |
Scalar * | p_dsdx |
Pointer on of *this , where (0x0 if not up to date). | |
Scalar * | p_dsdy |
Pointer on of *this , where (0x0 if not up to date). | |
Scalar * | p_dsdz |
Pointer on of *this , where (0x0 if not up to date). | |
Scalar * | p_lap |
Pointer on the Laplacian of *this (0x0 if not up to date). | |
Scalar * | p_lapang |
Pointer on the Laplacian of *this (0x0 if not up to date). | |
Scalar * | p_dsdradial |
Pointer on of *this . | |
Scalar * | p_dsdrho |
Pointer on of *this . | |
int | ind_lap |
Power of r by which the last computed Laplacian has been multiplied in the compactified external domain. | |
Tbl * | p_integ |
Pointer on the space integral of *this (values in each domain) (0x0 if not up to date). | |
const Map *const | mp |
Mapping on which the numerical values at the grid points are defined. | |
int | valence |
Valence of the tensor (0 = scalar, 1 = vector, etc...). | |
const Base_vect * | triad |
Vectorial basis (triad) with respect to which the tensor components are defined. | |
Itbl | type_indice |
1D array of integers (class Itbl ) of size valence containing the type of each index: COV for a covariant one and CON for a contravariant one. | |
int | n_comp |
Number of stored components, depending on the symmetry. | |
Scalar ** | cmp |
Array of size n_comp of pointers onto the components. | |
const Metric * | met_depend [N_MET_MAX] |
Array on the Metric 's which were used to compute derived quantities, like p_derive_cov , etc. | |
Tensor * | p_derive_cov [N_MET_MAX] |
Array of pointers on the covariant derivatives of this with respect to various metrics. | |
Tensor * | p_derive_con [N_MET_MAX] |
Array of pointers on the contravariant derivatives of this with respect to various metrics. | |
Tensor * | p_divergence [N_MET_MAX] |
Array of pointers on the divergence of this with respect to various metrics. | |
Friends | |
ostream & | operator<< (ostream &, const Scalar &) |
Display. | |
Scalar | operator- (const Scalar &) |
Scalar | operator+ (const Scalar &, const Scalar &) |
Scalar | operator+ (const Scalar &, const Mtbl &) |
Scalar | operator+ (const Scalar &, double) |
Scalar | operator- (const Scalar &, const Scalar &) |
Scalar | operator- (const Scalar &, const Mtbl &) |
Scalar | operator- (const Scalar &, double) |
Scalar | operator* (const Scalar &, const Scalar &) |
Scalar | operator% (const Scalar &, const Scalar &) |
Scalar | operator| (const Scalar &, const Scalar &) |
Scalar | operator* (const Mtbl &, const Scalar &) |
Scalar | operator* (double, const Scalar &) |
Scalar | operator/ (const Scalar &, const Scalar &) |
Scalar | operator/ (const Scalar &, const Mtbl &) |
Scalar | operator/ (const Mtbl &, const Scalar &) |
Scalar | operator/ (const Scalar &, double) |
Scalar | operator/ (double, const Scalar &) |
Scalar | sin (const Scalar &) |
Scalar | cos (const Scalar &) |
Scalar | tan (const Scalar &) |
Scalar | asin (const Scalar &) |
Scalar | acos (const Scalar &) |
Scalar | atan (const Scalar &) |
Scalar | exp (const Scalar &) |
Scalar | Heaviside (const Scalar &) |
Scalar | log (const Scalar &) |
Scalar | log10 (const Scalar &) |
Scalar | sqrt (const Scalar &) |
Scalar | racine_cubique (const Scalar &) |
Scalar | pow (const Scalar &, int) |
Scalar | pow (const Scalar &, double) |
Scalar | abs (const Scalar &) |
double | totalmax (const Scalar &) |
double | totalmin (const Scalar &) |
Tbl | max (const Scalar &) |
Tbl | min (const Scalar &) |
Tbl | norme (const Scalar &) |
Tbl | diffrel (const Scalar &a, const Scalar &b) |
Tbl | diffrelmax (const Scalar &a, const Scalar &b) |
class | Scalar |
class | Vector |
class | Sym_tensor |
class | Tensor_sym |
class | Metric |
Scalar | operator+ (const Tensor &, const Scalar &) |
Scalar | operator+ (const Scalar &, const Tensor &) |
Scalar | operator- (const Tensor &, const Scalar &) |
Scalar | operator- (const Scalar &, const Tensor &) |
Tensor_sym | operator* (const Tensor &, const Tensor_sym &) |
Tensor_sym | operator* (const Tensor_sym &, const Tensor &) |
Tensor_sym | operator* (const Tensor_sym &, const Tensor_sym &) |
Tensorial product of two symmetric tensors. |
Tensor field of valence 0 (or component of a tensorial field).
()
Definition at line 377 of file scalar.h.
Scalar::Scalar | ( | const Map & | mpi | ) | [explicit] |
Constructor from mapping.
Definition at line 197 of file scalar.C.
References Tensor::cmp, and set_der_0x0().
Scalar::Scalar | ( | const Tensor & | a | ) |
Constructor from a Tensor (must be of valence 0).
Definition at line 208 of file scalar.C.
References Tensor::cmp, set_der_0x0(), and Tensor::valence.
Scalar::Scalar | ( | const Scalar & | a | ) |
Copy constructor.
Definition at line 221 of file scalar.C.
References Tensor::cmp, and set_der_0x0().
Scalar::Scalar | ( | const Cmp & | a | ) |
Constructor by conversion of a Cmp.
Definition at line 231 of file scalar.C.
References Tensor::cmp, and set_der_0x0().
Constructor from a file (see sauve(FILE*)
).
Definition at line 241 of file scalar.C.
References Tensor::cmp, dzpuis, etat, fread_be(), Map::get_mg(), and set_der_0x0().
Scalar::~Scalar | ( | ) | [virtual] |
void Scalar::allocate_all | ( | ) | [virtual] |
Sets the logical state to ETATQCQ
(ordinary state) and performs the memory allocation of all the elements, down to the double
arrays of the Tbl
s.
This function performs in fact recursive calls to set_etat_qcq()
on each element of the chain Scalar
-> Valeur
-> Mtbl
-> Tbl
.
Reimplemented from Tensor.
Definition at line 358 of file scalar.C.
References Valeur::c, Mtbl::get_nzone(), Valeur::set_etat_c_qcq(), Tbl::set_etat_qcq(), Mtbl::set_etat_qcq(), set_etat_qcq(), Mtbl::t, and va.
void Scalar::annule | ( | int | l_min, | |
int | l_max | |||
) | [virtual] |
Sets the Scalar
to zero in several domains.
l_min | [input] The Scalar will be set (logically) to zero in the domains whose indices are in the range [l_min,l_max] . | |
l_max | [input] see the comments for l_min . |
Note that annule
(0,va.mg->get_nzone()-1) is equivalent to set_etat_zero()
.
Reimplemented from Tensor.
Definition at line 384 of file scalar.C.
References annule(), Valeur::annule(), etat, Mg3d::get_nzone(), Valeur::mg, p_dsdr, p_dsdradial, p_dsdt, p_dsdx, p_dsdy, p_dsdz, p_integ, p_lap, p_lapang, p_srdsdt, p_srstdsdp, p_stdsdp, set_etat_zero(), and va.
void Tensor::annule_domain | ( | int | l | ) | [inherited] |
Sets the Tensor
to zero in a given domain.
l | [input] Index of the domain in which the Tensor will be set (logically) to zero. |
Definition at line 662 of file tensor.C.
References Tensor::annule().
void Tensor::annule_extern_cn | ( | int | l_0, | |
int | deg | |||
) | [inherited] |
Performs a smooth (C^n) transition in a given domain to zero.
l_0 | [input] in the domain of index l0 the tensor is multiplied by the right polynomial (of degree 2n+1), to ensure continuty of the function and its n first derivative at both ends of this domain. The tensor is unchanged in the domains l < l_0 and set to zero in domains l > l_0. | |
deg | [input] the degree n of smoothness of the transition. |
Definition at line 686 of file tensor.C.
References allocate_all(), annule(), Itbl::annule_hard(), Tensor::cmp, Tensor::del_deriv(), Map::get_mg(), Mg3d::get_nr(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, Tensor::n_comp, pow(), Map::r, Tbl::set(), Itbl::set(), set_domain(), Tbl::set_etat_qcq(), std_spectral_base(), and Map::val_r().
void Scalar::annule_hard | ( | ) |
Sets the Scalar
to zero in a hard way.
1/ Sets the logical state to ETATQCQ
, i.e. to an ordinary state. 2/ Fills the Valeur
va
with zeros. NB: this function must be used for debugging purposes only. For other operations, the functions set_etat_zero()
or annule(int,int)
must be perferred.
Definition at line 373 of file scalar.C.
References Valeur::annule_hard(), del_deriv(), etat, and va.
void Scalar::annule_l | ( | int | l_min, | |
int | l_max, | |||
bool | ylm_output = false | |||
) |
Sets all the multipolar components between l_min
and l_max
to zero.
This is done for [ l_min
, l_max
] and all relevant m
in the spherical harmonics expansion basis. If ylm_output
is set to true
, then the spectral expansion basis of this
is left to be that of spherical harmonics.
Definition at line 107 of file scalar_manip.C.
References Valeur::base, Valeur::c, Valeur::c_cf, etat, Mtbl_cf::get_etat(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Base_val::give_quant_numbers(), Tensor::mp, Mtbl_cf::set(), set_etat_zero(), va, Valeur::ylm(), and Valeur::ylm_i().
Valeur ** Scalar::asymptot | ( | int | n, | |
const int | flag = 0 | |||
) | const |
Asymptotic expansion at r = infinity.
Determines the coefficients of the expansion
of *this
when .
n | order of the expansion | |
flag | : output |
n
+1 Valeur
s on mg->angu
describing the coefficients . This array is allocated by the routine. Definition at line 59 of file scalar_asymptot.C.
References Valeur::base, Valeur::c, dzpuis, Mg3d::get_angu(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, mult_r_ced(), Valeur::set(), Valeur::set_base(), Valeur::set_etat_c_qcq(), Tbl::set_etat_qcq(), Mtbl::set_etat_qcq(), Tbl::set_etat_zero(), Valeur::set_etat_zero(), set_grid_point(), Mtbl::t, va, and val_grid_point().
Performs one time-step integration (from ) of the scalar d'Alembert equation with *this
being the value of the function f at time J .
Works only with an affine mapping (class Map_af
) and, if the last domain is a compactified one, it simply copies the value of the field in this last domain at the time-step J to the last domain of the returned solution.
par | [input/output] possible parameters to control the resolution of the d'Alembert equation:
| |
fJm1 | [input] solution at time J-1 | |
source | [input] source of the d'Alembert equation . |
par.get_int(0)
. Definition at line 213 of file scalar_pde.C.
References Map::dalembert(), and Tensor::mp.
void Scalar::change_triad | ( | const Base_vect & | new_triad | ) | [virtual] |
bool Scalar::check_dzpuis | ( | int | dzi | ) | const |
Returns false
if the last domain is compactified and *this
is not zero in this domain and dzpuis
is not equal to dzi
, otherwise return true.
Definition at line 866 of file scalar.C.
References dz_nonzero(), and dzpuis.
Computes the Lie derivative of this
with respect to some vector field v
(protected method; the public interface is method derive_lie
).
Definition at line 335 of file tensor_calculus.C.
References Tensor::cmp, contract(), dec_dzpuis(), Tensor::derive_cov(), Map::flat_met_cart(), Map::flat_met_spher(), get_dzpuis(), Tensor::get_n_comp(), Tensor::get_triad(), Tensor::indices(), Tensor::mp, Tensor::n_comp, Tensor::operator()(), Tensor::set(), Itbl::set(), Tensor::triad, Tensor::type_indice, and Tensor::valence.
void Scalar::dec_dzpuis | ( | int | dec = 1 |
) | [virtual] |
Decreases by dec
units the value of dzpuis
and changes accordingly the values of the Scalar
in the compactified external domain (CED).
Reimplemented from Tensor.
Definition at line 414 of file scalar_r_manip.C.
References Map::dec2_dzpuis(), Map::dec_dzpuis(), del_deriv(), etat, and Tensor::mp.
void Scalar::del_deriv | ( | ) | const [protected, virtual] |
void Tensor::del_derive_met | ( | int | j | ) | const [protected, virtual, inherited] |
Logical destructor of the derivatives depending on the i-th element of met_depend
.
Reimplemented in Sym_tensor, and Vector.
Definition at line 410 of file tensor.C.
References Tensor::met_depend, Tensor::p_derive_con, Tensor::p_derive_cov, Tensor::p_divergence, Tensor::set_der_met_0x0(), and Metric::tensor_depend.
void Scalar::del_t | ( | ) | [protected] |
Logical destructor.
Definition at line 272 of file scalar.C.
References del_deriv(), Valeur::del_t(), Valeur::set_etat_nondef(), and va.
const Scalar & Scalar::deriv | ( | int | i | ) | const |
Returns the "contravariant" derivative of *this
with respect to some metric , by raising the index of the gradient (cf.
method derive_cov()
) with .
Reimplemented from Tensor.
Definition at line 401 of file scalar_deriv.C.
Returns the gradient (1-form = covariant vector) of *this
.
gam | metric components only used to get the triad with respect to which the components of the result are defined |
Reimplemented from Tensor.
Definition at line 389 of file scalar_deriv.C.
Computes the derivative of this
along a vector field v
.
Reimplemented from Tensor.
Definition at line 418 of file scalar_deriv.C.
References Tensor::compute_derive_lie(), and Tensor::mp.
void Scalar::div_cost | ( | ) |
Division by .
Definition at line 66 of file scalar_th_manip.C.
References del_deriv(), Map::div_cost(), and Tensor::mp.
void Scalar::div_r | ( | ) |
Division by r everywhere; dzpuis
is not changed.
Definition at line 129 of file scalar_r_manip.C.
References del_deriv(), Map::div_r(), and Tensor::mp.
void Scalar::div_r_ced | ( | ) |
Division by r in the compactified external domain (CED), the dzpuis
flag is not changed.
Definition at line 192 of file scalar_r_manip.C.
References del_deriv(), Map::div_r_zec(), and Tensor::mp.
void Scalar::div_r_dzpuis | ( | int | ced_mult_r | ) |
Division by r everywhere but with the output flag dzpuis
set to ced_mult_r
.
ced_mult_r | [input] value of dzpuis of the result. |
Definition at line 143 of file scalar_r_manip.C.
References allocate_all(), annule(), Valeur::base, dec_dzpuis(), del_deriv(), div_r(), domain(), dzpuis, etat, Valeur::get_base(), Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), inc_dzpuis(), Tensor::mp, set_domain(), set_etat_qcq(), set_spectral_base(), Base_val::sx(), and va.
void Scalar::div_rsint | ( | ) |
Division by everywhere; dzpuis
is not changed.
Definition at line 344 of file scalar_r_manip.C.
References del_deriv(), Map::div_rsint(), and Tensor::mp.
void Scalar::div_rsint_dzpuis | ( | int | ced_mult_r | ) |
Division by but with the output flag dzpuis
set to ced_mult_r
.
ced_mult_r | [input] value of dzpuis of the result. |
Definition at line 358 of file scalar_r_manip.C.
References allocate_all(), annule(), Valeur::base, dec_dzpuis(), del_deriv(), div_rsint(), div_sint(), domain(), dzpuis, etat, Valeur::get_base(), Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), inc_dzpuis(), Tensor::mp, set_domain(), set_dzpuis(), set_etat_qcq(), set_spectral_base(), Base_val::ssint(), Base_val::sx(), and va.
void Scalar::div_sint | ( | ) |
Division by .
Definition at line 94 of file scalar_th_manip.C.
References del_deriv(), Map::div_sint(), and Tensor::mp.
void Scalar::div_tant | ( | ) |
Division by .
Definition at line 107 of file scalar_th_manip.C.
References del_deriv(), Map::div_tant(), and Tensor::mp.
Computes the divergence of this
with respect to some metric .
The divergence is taken with respect of the last index of this
which thus must be contravariant. For instance if the tensor represented by this
is a twice contravariant tensor, whose components w.r.t. the triad are : , the divergence of w.r.t. is the vector
where denotes the connection associated with the metric .
gam | metric |
this
with respect to . Reimplemented in Sym_tensor, and Vector.
Definition at line 1051 of file tensor.C.
References Metric::connect(), Tensor::get_place_met(), Connection::p_divergence(), Tensor::p_divergence, and Tensor::set_dependance().
const Tbl& Scalar::domain | ( | int | l | ) | const [inline] |
Computes a new tensor by lowering an index of *this
.
ind | index to be lowered, with the following convention :
| |
gam | metric used to lower the index (contraction with the twice covariant form of the metric on the index ind ). |
Definition at line 261 of file tensor_calculus.C.
References contract(), Metric::cov(), Tensor::indices(), Tensor::mp, Tensor::n_comp, Tensor::set(), Itbl::set(), Tensor::triad, Tensor::type_indice, and Tensor::valence.
const Scalar & Scalar::dsdr | ( | ) | const |
Returns of *this
.
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 112 of file scalar_deriv.C.
References Map::dsdr(), dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_dsdr, Scalar(), set_dzpuis(), and set_etat_zero().
const Scalar & Scalar::dsdradial | ( | ) | const |
Returns of *this
if the mapping is affine (class Map_af
) and of *this
if the mapping is logarithmic (class Map_log
).
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 490 of file scalar_deriv.C.
References Map::dsdradial(), dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_dsdradial, Scalar(), set_dzpuis(), and set_etat_zero().
const Scalar & Scalar::dsdrho | ( | ) | const |
Returns of *this
.
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 521 of file scalar_deriv.C.
References dsdr(), dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), get_spectral_va(), Mg3d::get_type_r(), Tensor::mp, Valeur::mult_ct(), Valeur::mult_st(), p_dsdrho, Scalar(), set_dzpuis(), set_etat_zero(), and srdsdt().
const Scalar & Scalar::dsdt | ( | ) | const |
Returns of *this
.
Definition at line 207 of file scalar_deriv.C.
References Map::dsdt(), dzpuis, etat, Tensor::mp, p_dsdt, Scalar(), set_dzpuis(), and set_etat_zero().
const Scalar & Scalar::dsdx | ( | ) | const |
Returns of *this
, where .
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 265 of file scalar_deriv.C.
References Map::comp_x_from_spherical(), dsdr(), dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_dsdx, Scalar(), set_dzpuis(), set_etat_zero(), srdsdt(), and srstdsdp().
const Scalar & Scalar::dsdy | ( | ) | const |
Returns of *this
, where .
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 296 of file scalar_deriv.C.
References Map::comp_y_from_spherical(), dsdr(), dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_dsdy, Scalar(), set_dzpuis(), set_etat_zero(), srdsdt(), and srstdsdp().
const Scalar & Scalar::dsdz | ( | ) | const |
Returns of *this
, where .
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 327 of file scalar_deriv.C.
References Map::comp_z_from_spherical(), dsdr(), dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_dsdz, Scalar(), set_dzpuis(), set_etat_zero(), and srdsdt().
bool Scalar::dz_nonzero | ( | ) | const |
Returns true
if the last domain is compactified and *this
is not zero in this domain.
Definition at line 807 of file scalar.C.
References Valeur::c, Valeur::c_cf, Valeur::etat, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, and va.
void Scalar::exp_filter_r_all_domains | ( | Scalar & | ss, | |
int | p, | |||
double | alpha = -16. | |||
) |
Applies an exponential filter in radial direction in all domains.
(see Scalar:exponential_filter_r
). Note that this may cause regularity problems at the origin if applied in a nucleus.
void Scalar::exponential_filter_r | ( | int | lzmin, | |
int | lzmax, | |||
int | p, | |||
double | alpha = -16. | |||
) | [virtual] |
Applies an exponential filter to the spectral coefficients in the radial direction.
The filter is of the type: , with and N the number of radial coefficients.
lzmin,lzmax | [input] the indices of the domain where the filter is applied (in [lzmin , lzmax ]) | |
p | [input] the order of the filter | |
alpha | [input] appearing in the above formula. |
Reimplemented from Tensor.
Definition at line 56 of file scalar_exp_filter.C.
References Valeur::c, Valeur::c_cf, Valeur::coef(), del_deriv(), Valeur::del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, Mtbl_cf::set(), and va.
void Scalar::exponential_filter_ylm | ( | int | lzmin, | |
int | lzmax, | |||
int | p, | |||
double | alpha = -16. | |||
) | [virtual] |
Applies an exponential filter to the spectral coefficients in the angular directions.
The filter is of the type: , with defined for Scalar::exponential_filter_r
and the number of spherical harmonics used.
Reimplemented from Tensor.
Definition at line 138 of file scalar_exp_filter.C.
References Valeur::base, Valeur::c, Valeur::c_cf, del_deriv(), Valeur::del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Base_val::give_lmax(), Base_val::give_quant_numbers(), Tensor::mp, Mtbl_cf::set(), va, Valeur::ylm(), and Valeur::ylm_i().
void Scalar::filtre | ( | int | n | ) |
Sets the n
lasts coefficients in r to 0 in the external domain.
Definition at line 147 of file scalar_manip.C.
References Valeur::c_cf, Valeur::coef(), del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, Mtbl_cf::set(), Valeur::set_etat_cf_qcq(), and va.
void Scalar::filtre_phi | ( | int | n, | |
int | zone | |||
) |
Sets the n
lasts coefficients in to 0 in the domain zone
.
Definition at line 245 of file scalar_manip.C.
References Valeur::c_cf, Valeur::coef(), del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Tensor::mp, Mtbl_cf::set(), Valeur::set_etat_cf_qcq(), and va.
void Scalar::filtre_r | ( | int | n, | |
int | nzone | |||
) |
Sets the n
last coefficients in r to 0 in the domain nzone
.
Definition at line 214 of file scalar_manip.C.
References Valeur::c, Valeur::c_cf, Valeur::coef(), del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Tensor::mp, Mtbl_cf::set(), Valeur::set_etat_cf_qcq(), and va.
void Scalar::filtre_r | ( | int * | nn | ) |
Sets the n
lasts coefficients in r to 0 in all domains.
Definition at line 176 of file scalar_manip.C.
References Valeur::c, Valeur::c_cf, Valeur::coef(), del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, Mtbl_cf::set(), Valeur::set_etat_cf_qcq(), and va.
void Scalar::filtre_tp | ( | int | nn, | |
int | nz1, | |||
int | nz2 | |||
) |
Sets the n
lasts coefficients in to 0 in the domain nz1
to nz2
when expressed in spherical harmonics.
Definition at line 266 of file scalar_manip.C.
References Valeur::filtre_tp(), and va.
void Scalar::fixe_decroissance | ( | int | puis | ) |
Substracts all the components behaving like in the external domain, with n strictly lower than puis
, so that *this
decreases at least like at infinity.
Definition at line 342 of file scalar_manip.C.
References Valeur::base, Valeur::c_cf, Valeur::coef(), dzpuis, Map_af::get_alpha(), Base_val::get_base_r(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, mult_r_ced(), R_CHEBU, Mtbl_cf::set(), Valeur::set_etat_cf_qcq(), and va.
int Scalar::get_dzpuis | ( | ) | const [inline] |
int Scalar::get_etat | ( | ) | const [inline] |
Itbl Tensor::get_index_type | ( | ) | const [inline, inherited] |
Returns the types of all the indices.
Itbl
) of size valence
containing the type of each index, COV
for a covariant one and CON
for a contravariant one. Definition at line 892 of file tensor.h.
References Tensor::type_indice.
int Tensor::get_index_type | ( | int | i | ) | const [inline, inherited] |
Gives the type (covariant or contravariant) of the index number i
.
i
must be strictly lower than valence
and obey the following convention:
i
= 0 : first index i
= 1 : second index Definition at line 882 of file tensor.h.
References Tensor::type_indice.
const Map& Tensor::get_mp | ( | ) | const [inline, inherited] |
int Tensor::get_n_comp | ( | ) | const [inline, inherited] |
Returns the number of stored components.
Definition at line 868 of file tensor.h.
References Tensor::n_comp.
int Tensor::get_place_met | ( | const Metric & | metre | ) | const [protected, inherited] |
Returns the position of the pointer on metre
in the array met_depend
.
Definition at line 439 of file tensor.C.
References Tensor::met_depend.
const Base_val& Scalar::get_spectral_base | ( | ) | const [inline] |
Returns the spectral bases of the Valeur
va
.
Definition at line 1277 of file scalar.h.
References Valeur::base, and va.
const Valeur& Scalar::get_spectral_va | ( | ) | const [inline] |
const Base_vect* Tensor::get_triad | ( | ) | const [inline, inherited] |
Returns the vectorial basis (triad) on which the components are defined.
Definition at line 862 of file tensor.h.
References Tensor::triad.
int Tensor::get_valence | ( | ) | const [inline, inherited] |
void Scalar::import | ( | int | nzet, | |
const Scalar & | ci | |||
) |
Assignment to another Scalar
defined on a different mapping.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 76 of file scalar_import.C.
References Map::get_bvect_cart(), Tensor::get_mp(), import_align(), import_anti(), import_gal(), and Tensor::mp.
void Scalar::import | ( | const Scalar & | ci | ) |
Assignment to another Scalar
defined on a different mapping.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
ci | [input] Scalar to be imported. |
Definition at line 64 of file scalar_import.C.
References Map::get_mg(), Mg3d::get_nzone(), and Tensor::mp.
void Scalar::import_align | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings have aligned Cartesian axis.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 528 of file scalar_import.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::c_cf, Valeur::coef(), del_t(), Map::get_bvect_cart(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), get_spectral_va(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Mtbl_cf::val_point(), Map::x, Map::y, and Map::z.
void Scalar::import_align_asymy | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings have aligned Cartesian axis.
Case where the Scalar
is antisymmetric with respect to the plane y=0.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 374 of file scalar_import_asymy.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::c_cf, Valeur::coef(), del_t(), Map::get_bvect_cart(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), get_spectral_va(), Mg3d::get_type_p(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Mtbl_cf::val_point_asymy(), Map::x, Map::y, and Map::z.
void Scalar::import_align_symy | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings have aligned Cartesian axis.
Case where the Scalar
is symmetric with respect to the plane y=0.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 345 of file scalar_import_symy.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::c_cf, Valeur::coef(), del_t(), Map::get_bvect_cart(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), get_spectral_va(), Mg3d::get_type_p(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Mtbl_cf::val_point_symy(), Map::x, Map::y, and Map::z.
void Scalar::import_anti | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings have anti-aligned Cartesian axis (i.e.
, , ).
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 331 of file scalar_import.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::c_cf, Valeur::coef(), del_t(), Map::get_bvect_cart(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), get_spectral_va(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Mtbl_cf::val_point(), Map::x, Map::y, and Map::z.
void Scalar::import_anti_asymy | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings have anti-aligned Cartesian axis (i.e.
, , ). Case where the Scalar
is antisymmetric with respect to the plane y=0.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 127 of file scalar_import_asymy.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::c_cf, Valeur::coef(), del_t(), Map::get_bvect_cart(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), get_spectral_va(), Mg3d::get_type_p(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Mtbl_cf::val_point_asymy(), Map::x, Map::y, and Map::z.
void Scalar::import_anti_symy | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings have anti-aligned Cartesian axis (i.e.
, , ). Case where the Scalar
is symmetric with respect to the plane y=0.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 127 of file scalar_import_symy.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::c_cf, Valeur::coef(), del_t(), Map::get_bvect_cart(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), get_spectral_va(), Mg3d::get_type_p(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Mtbl_cf::val_point_symy(), Map::x, Map::y, and Map::z.
void Scalar::import_asymy | ( | int | nzet, | |
const Scalar & | ci | |||
) |
Assignment to another Scalar
defined on a different mapping.
Case where the Scalar
is antisymmetric with respect to the plane y=0. This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 80 of file scalar_import_asymy.C.
References Map::get_bvect_cart(), Tensor::get_mp(), import_align_asymy(), import_anti_asymy(), and Tensor::mp.
void Scalar::import_asymy | ( | const Scalar & | ci | ) |
Assignment to another Scalar
defined on a different mapping.
Case where the Scalar
is antisymmetric with respect to the plane y=0. This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
ci | [input] Scalar to be imported. |
Definition at line 68 of file scalar_import_asymy.C.
References Map::get_mg(), Mg3d::get_nzone(), and Tensor::mp.
void Scalar::import_gal | ( | int | nzet, | |
const Scalar & | ci | |||
) | [protected] |
Assignment to another Scalar
defined on a different mapping, when the two mappings do not have a particular relative orientation.
This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 127 of file scalar_import.C.
References Param::add_double(), Param::add_int(), Param::add_int_mod(), annule(), Valeur::c, Valeur::coef(), del_t(), get_dzpuis(), get_etat(), Map::get_mg(), Tensor::get_mp(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Map::get_ori_x(), Map::get_ori_y(), Map::get_ori_z(), Map::get_rot_phi(), get_spectral_va(), Tensor::mp, Map::phi, Map::r, set_dzpuis(), Valeur::set_etat_c_qcq(), set_etat_one(), Mtbl::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Tbl::t, Mtbl::t, Map::tet, va, Map::val_lx(), Map::xa, Map::ya, and Map::za.
void Scalar::import_symy | ( | int | nzet, | |
const Scalar & | ci | |||
) |
Assignment to another Scalar
defined on a different mapping.
Case where the Scalar
is symmetric with respect to the plane y=0. This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
nzet | [input] Number of domains of the destination mapping (i.e. this->mp ) where the importation is performed: the assignment is done for the domains whose indices are between 0 and nzet-1 . In the other domains, *this is set to zero. | |
ci | [input] Scalar to be imported. |
Definition at line 80 of file scalar_import_symy.C.
References Map::get_bvect_cart(), Tensor::get_mp(), import_align_symy(), import_anti_symy(), and Tensor::mp.
void Scalar::import_symy | ( | const Scalar & | ci | ) |
Assignment to another Scalar
defined on a different mapping.
Case where the Scalar
is symmetric with respect to the plane y=0. This assignment is performed point to point by means of the spectral expansion of the original Scalar
.
ci | [input] Scalar to be imported. |
Definition at line 68 of file scalar_import_symy.C.
References Map::get_mg(), Mg3d::get_nzone(), and Tensor::mp.
void Scalar::inc_dzpuis | ( | int | inc = 1 |
) | [virtual] |
Increases by inc
units the value of dzpuis
and changes accordingly the values of the Scalar
in the compactified external domain (CED).
Reimplemented from Tensor.
Definition at line 466 of file scalar_r_manip.C.
References del_deriv(), etat, Map::inc2_dzpuis(), Map::inc_dzpuis(), and Tensor::mp.
Itbl Tensor::indices | ( | int | pos | ) | const [virtual, inherited] |
Returns the indices of a component given by its position in the array cmp
.
pos | [input] position in the array cmp of the pointer to the Scalar representing a component |
Itbl
) of size valence
giving the value of each index for the component located at the position pos
in the array cmp
, with the following storage convention: Itbl(0)
: value of the first index (1, 2 or 3) Itbl(1)
: value of the second index (1, 2 or 3) Reimplemented in Tensor_sym, and Vector.
Definition at line 535 of file tensor.C.
References Tensor::n_comp, Itbl::set(), and Tensor::valence.
double Scalar::integrale | ( | ) | const |
Computes the integral over all space of *this
.
The computed quantity is (u being the field represented by *this
) . Note that in the compactified external domain (CED), dzpuis
must be 4 for the computation to take place.
Definition at line 57 of file scalar_integ.C.
References Map::get_mg(), Mg3d::get_nzone(), integrale_domains(), and Tensor::mp.
const Tbl & Scalar::integrale_domains | ( | ) | const |
Computes the integral in each domain of *this
.
The computed quantity is (u being the field represented by *this
) in each domain. The result is returned a Tbl
on the various domains. Note that in the compactified external domain (CED), dzpuis
must be 4 for the computation to take place.
Definition at line 75 of file scalar_integ.C.
References etat, Map::integrale(), Tensor::mp, and p_integ.
const Scalar & Scalar::lapang | ( | ) | const |
Returns the angular Laplacian of *this
, where .
Definition at line 460 of file scalar_deriv.C.
References dzpuis, etat, Map::lapang(), Tensor::mp, p_lapang, Scalar(), set_dzpuis(), and set_etat_zero().
const Scalar & Scalar::laplacian | ( | int | ced_mult_r = 4 |
) | const |
Returns the Laplacian of *this
.
ced_mult_r | [input] Determines the quantity computed in the compactified external domain (CED) (u in the field represented by *this ) :
|
Definition at line 435 of file scalar_deriv.C.
References etat, ind_lap, Map::laplacien(), Tensor::mp, p_lap, and Scalar().
Method for matching accross domains and imposing boundary condition.
Matching of the field represented by this
accross domains and imposition of the boundary condition using the collocation method.
Method for matching accross domains and imposing boundary condition.
Matching of the field represented by this
accross domains and imposition of the boundary condition using the tau method.
par_bc | [input] Param to control the boundary conditions par_bc must contain (at a minimum) a modifiable Tbl which specifies a physical boundary two integers, one specifying the domain that has the boundary the other specifying the number of conditions 1 -> Dirichlet 2 -> Robin (which may reduce to von Neumann, see below) two doubles, specifying the Robin BC parameters. If the first is zero, we see that Robin will reduce to von Neumann | |
par_mat | [input/output] optional Param in which the matching matrices are stored (together with their LU decomposition). |
Definition at line 62 of file scalar_match_tau.C.
References Param::add_matrice_mod(), Tbl::annule_hard(), Matrice::annule_hard(), annule_hard(), Valeur::c, Valeur::c_cf, Valeur::coef(), etat, Map_af::get_alpha(), Param::get_double_mod(), Param::get_int(), Param::get_matrice_mod(), Map::get_mg(), Param::get_n_double_mod(), Param::get_n_int(), Param::get_n_matrice_mod(), Param::get_n_tbl_mod(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), get_spectral_base(), Param::get_tbl_mod(), Mg3d::get_type_r(), Base_val::give_quant_numbers(), Matrice::inverse(), Tensor::mp, R_CHEBI, R_CHEBP, Mtbl_cf::set(), Tbl::set(), Matrice::set(), Matrice::set_lu(), va, and Valeur::ylm().
Method for matching accross domains and imposing boundary condition.
Matching of the field represented by this
accross domains and imposition of the boundary condition using the tau method.
void Scalar::mult_cost | ( | ) |
Multiplication by .
Definition at line 54 of file scalar_th_manip.C.
References del_deriv(), Tensor::mp, and Map::mult_cost().
void Scalar::mult_r | ( | ) |
Multiplication by r everywhere; dzpuis
is not changed.
Definition at line 204 of file scalar_r_manip.C.
References del_deriv(), Tensor::mp, and Map::mult_r().
void Scalar::mult_r_ced | ( | ) |
Multiplication by r in the compactified external domain (CED), the dzpuis
flag is not changed.
Definition at line 265 of file scalar_r_manip.C.
References del_deriv(), Tensor::mp, and Map::mult_r_zec().
void Scalar::mult_r_dzpuis | ( | int | ced_mult_r | ) |
Multiplication by r everywhere but with the output flag dzpuis
set to ced_mult_r
.
ced_mult_r | [input] value of dzpuis of the result. |
Definition at line 217 of file scalar_r_manip.C.
References allocate_all(), annule(), Valeur::base, dec_dzpuis(), del_deriv(), domain(), dzpuis, etat, Valeur::get_base(), Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), inc_dzpuis(), Tensor::mp, mult_r(), Base_val::mult_x(), set_domain(), set_etat_qcq(), set_spectral_base(), and va.
void Scalar::mult_rsint | ( | ) |
Multiplication by everywhere; dzpuis
is not changed.
Definition at line 277 of file scalar_r_manip.C.
References del_deriv(), Tensor::mp, and Map::mult_rsint().
void Scalar::mult_rsint_dzpuis | ( | int | ced_mult_r | ) |
Multiplication by but with the output flag dzpuis
set to ced_mult_r
.
ced_mult_r | [input] value of dzpuis of the result. |
Definition at line 290 of file scalar_r_manip.C.
References allocate_all(), annule(), Valeur::base, dec_dzpuis(), del_deriv(), domain(), dzpuis, etat, Valeur::get_base(), Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), inc_dzpuis(), Tensor::mp, mult_rsint(), Base_val::mult_sint(), mult_sint(), Base_val::mult_x(), set_domain(), set_dzpuis(), set_etat_qcq(), set_spectral_base(), and va.
void Scalar::mult_sint | ( | ) |
Multiplication by .
Definition at line 80 of file scalar_th_manip.C.
References del_deriv(), Tensor::mp, and Map::mult_sint().
Tbl Scalar::multipole_spectrum | ( | ) | const |
Gives the spectrum in terms of multipolar modes l .
Tbl
of size (nzone, lmax), where lmax is the maximal multipolar momentum over all domains. The l -th element contains the L1 norm of the l -th multipole (i.e. a sum over all m of the norms (coefficient space) of the component of a given . Definition at line 949 of file scalar.C.
References Tbl::annule_hard(), Mtbl_cf::base, Valeur::c_cf, Valeur::coef(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, Tbl::set(), Tbl::set_etat_zero(), va, and Valeur::ylm().
Scalar::operator Cmp | ( | ) | const |
Conversion to a Cmp
.
Definition at line 667 of file scalar.C.
References dzpuis, Tensor::mp, Cmp::set_dzpuis(), and va.
const Scalar & Tensor::operator() | ( | int | i1, | |
int | i2, | |||
int | i3, | |||
int | i4 | |||
) | const [inherited] |
Returns the value of a component for a tensor of valence 4 (read-only version).
i1 | value of the first index (1, 2 or 3) | |
i2 | value of the second index (1, 2 or 3) | |
i3 | value of the third index (1, 2 or 3) | |
i4 | value of the fourth index (1, 2 or 3) |
(i1,i2,i3,i4) Definition at line 779 of file tensor.C.
References Tensor::cmp, Tensor::position(), Itbl::set(), and Tensor::valence.
const Scalar & Tensor::operator() | ( | int | i1, | |
int | i2, | |||
int | i3 | |||
) | const [inherited] |
Returns the value of a component for a tensor of valence 3 (read-only version).
i1 | value of the first index (1, 2 or 3) | |
i2 | value of the second index (1, 2 or 3) | |
i3 | value of the third index (1, 2 or 3) |
(i1,i2,i3) Definition at line 767 of file tensor.C.
References Tensor::cmp, Tensor::position(), Itbl::set(), and Tensor::valence.
const Scalar & Tensor::operator() | ( | int | i1, | |
int | i2 | |||
) | const [inherited] |
Returns the value of a component for a tensor of valence 2 (read-only version).
i1 | value of the first index (1, 2 or 3) | |
i2 | value of the second index (1, 2 or 3) |
(i1,i2) Definition at line 756 of file tensor.C.
References Tensor::cmp, Tensor::position(), Itbl::set(), and Tensor::valence.
Returns the value of a component (read-only version).
ind | 1-D Itbl of size valence containing the values of each index specifing the component, with the following storage convention:
|
ind
Definition at line 794 of file tensor.C.
References Tensor::cmp, Itbl::get_dim(), Itbl::get_ndim(), Tensor::position(), and Tensor::valence.
void Scalar::operator*= | ( | const Scalar & | ci | ) |
*= Scalar
Definition at line 933 of file scalar_arithm.C.
References del_deriv(), dzpuis, etat, get_etat(), Tensor::get_mp(), Tensor::mp, operator=(), set_etat_nondef(), set_etat_zero(), and va.
void Scalar::operator+= | ( | const Scalar & | ci | ) |
+= Scalar
Reimplemented from Tensor.
Definition at line 829 of file scalar_arithm.C.
References del_deriv(), dz_nonzero(), dzpuis, etat, get_etat(), Tensor::get_mp(), Tensor::mp, set_dzpuis(), set_etat_nondef(), and va.
void Scalar::operator-= | ( | const Scalar & | ci | ) |
-= Scalar
Reimplemented from Tensor.
Definition at line 882 of file scalar_arithm.C.
References del_deriv(), dz_nonzero(), dzpuis, etat, get_etat(), Tensor::get_mp(), Tensor::mp, set_dzpuis(), set_etat_nondef(), and va.
void Scalar::operator= | ( | int | n | ) |
Assignment to an int
.
Definition at line 646 of file scalar.C.
References del_deriv(), dzpuis, set_etat_one(), set_etat_qcq(), set_etat_zero(), and va.
void Scalar::operator= | ( | double | x | ) |
Assignment to a double
.
Definition at line 625 of file scalar.C.
References del_deriv(), dzpuis, set_etat_one(), set_etat_qcq(), set_etat_zero(), and va.
void Scalar::operator= | ( | const Mtbl & | a | ) |
Assignment to a Mtbl
.
Definition at line 585 of file scalar.C.
References Valeur::c, del_deriv(), Valeur::del_t(), Mtbl::get_etat(), set_etat_qcq(), set_etat_zero(), and va.
void Scalar::operator= | ( | const Valeur & | a | ) |
Assignment to a Valeur
.
Definition at line 542 of file scalar.C.
References del_deriv(), Valeur::del_t(), Valeur::get_etat(), set_etat_qcq(), set_etat_zero(), and va.
void Scalar::operator= | ( | const Cmp & | a | ) |
Assignment to a Cmp
.
Definition at line 492 of file scalar.C.
References del_deriv(), Valeur::del_t(), dzpuis, Cmp::get_dzpuis(), Cmp::get_etat(), Cmp::get_mp(), Tensor::mp, set_etat_nondef(), set_etat_qcq(), set_etat_zero(), Cmp::va, and va.
void Scalar::operator= | ( | const Tensor & | a | ) | [virtual] |
Assignment to a Tensor
(of valence 0).
Definition at line 429 of file scalar.C.
References Tensor::cmp, operator=(), and Tensor::valence.
void Scalar::operator= | ( | const Scalar & | a | ) | [virtual] |
Assignment to another Scalar
defined on the same mapping.
Reimplemented from Tensor.
Definition at line 439 of file scalar.C.
References del_deriv(), Valeur::del_t(), dzpuis, etat, Valeur::get_base(), Tensor::mp, Valeur::set_base(), set_etat_nondef(), set_etat_one(), set_etat_qcq(), set_etat_zero(), and va.
Solves the scalar Poisson equation with *this
as a source (version with parameters to control the resolution).
The source of the equation is represented by the Scalar
*this
. Note that dzpuis
must be equal to 2 or 4, i.e. that the quantity stored in *this
is in fact or in the compactified external domain.
par | [input/output] possible parameters | |
uu | [input/output] solution u with the boundary condition u =0 at spatial infinity. |
Definition at line 148 of file scalar_pde.C.
References Tensor::mp, and Map::poisson().
Scalar Scalar::poisson | ( | ) | const |
Solves the scalar Poisson equation with *this
as a source.
The source of the equation is represented by the Scalar
*this
. Note that dzpuis
must be equal to 2 or 4, i.e. that the quantity stored in *this
is in fact or in the compactified external domain. The solution u with the boundary condition u =0 at spatial infinity is the returned Scalar
.
Definition at line 132 of file scalar_pde.C.
References Tensor::mp, and Map::poisson().
Scalar Scalar::poisson_angu | ( | double | lambda = 0 |
) | const |
Solves the (generalized) angular Poisson equation with *this
as source.
The generalized angular Poisson equation is , where .
lambda | [input] coefficient in the above equation (default value = 0) |
Definition at line 196 of file scalar_pde.C.
References Tensor::mp, and Map::poisson_angu().
Scalar Scalar::poisson_dir_neu | ( | const Valeur & | limite, | |
int | num, | |||
double | fact_dir, | |||
double | fact_neu | |||
) | const |
Is identicall to Scalar::poisson()
.
The regularity condition at the origin is replace by a mixed boundary condition (Dirichlet + Neumann).
limite | [input] : angular function. The boundary condition is given by limite [num] . | |
num | [input] : index of the boudary at which the condition is to be fullfilled. | |
fact_dir | [input] : double in front of (if is the variable solved). | |
fact_neu | [input] : double in front of the radial derivative of . |
More precisely we impose is equal to the source at the boundary between the domains num
and num+1
(the latter one being a shell).
Definition at line 76 of file scalar_pde_frontiere.C.
References Tensor::mp, Map::poisson_frontiere(), Cmp::set_dzpuis(), and Cmp::va.
Is identicall to Scalar::poisson()
.
The regularity condition at the origin is replace by a boundary condition of the Dirichlet type.
limite | [input] : angular function. The boundary condition is given by limite [num] . | |
num | [input] : index of the boudary at which the condition is to be fullfilled. |
More precisely we impose the solution is equal to limite
[num] at the boundary between the domains num
and num+1
(the latter one being a shell).
Definition at line 60 of file scalar_pde_frontiere.C.
References Tensor::mp, Map::poisson_frontiere(), Cmp::set_dzpuis(), and Cmp::va.
Scalar Scalar::poisson_frontiere_double | ( | const Valeur & | lim_func, | |
const Valeur & | lim_der, | |||
int | num_zone | |||
) | const |
Idem as Scalar::poisson_dirichlet
, the boundary condition being on both the function and its radial derivative.
The boundary condition at infinity is relaxed.
Definition at line 111 of file scalar_pde_frontiere.C.
References Tensor::mp.
Idem as Scalar::poisson_dirichlet
, the boundary condition being on the radial derivative of the solution.
Definition at line 94 of file scalar_pde_frontiere.C.
References Tensor::mp, Map::poisson_frontiere(), Cmp::set_dzpuis(), and Cmp::va.
void Scalar::poisson_regular | ( | int | k_div, | |
int | nzet, | |||
double | unsgam1, | |||
Param & | par, | |||
Scalar & | uu, | |||
Scalar & | uu_regu, | |||
Scalar & | uu_div, | |||
Tensor & | duu_div, | |||
Scalar & | source_regu, | |||
Scalar & | source_div | |||
) | const |
Solves the scalar Poisson equation with *this
as a source (version with parameters to control the resolution).
The source of the equation is represented by the Scalar
*this
. The regularized source is constructed and solved. Note that dzpuis
must be equal to 2 or 4, i.e. that the quantity stored in *this
is in fact or in the compactified external domain.
k_div | [input] regularization degree of the procedure | |
nzet | [input] number of domains covering the star | |
unsgam1 | [input] parameter where denotes the adiabatic index | |
par | [input/output] possible parameters | |
uu | [input/output] solution | |
uu_regu | [output] solution of the regular part of the source. | |
uu_div | [output] solution of the diverging part of the source. | |
duu_div | [output] derivative of the diverging potential. | |
source_regu | [output] regularized source | |
source_div | [output] diverging part of the source |
Definition at line 57 of file scalar_poisson_regu.C.
References Tensor::get_triad(), Tensor::mp, Map::poisson_regular(), Tensor::set(), Tenseur::set(), Itbl::set(), Itbl::set_etat_qcq(), and Tenseur::set_etat_qcq().
Solves the scalar Poisson equation with *this
as a source using a real Tau method (version with parameters to control the resolution) The source of the equation is represented by the Scalar
*this
.
Note that dzpuis
must be equal to 2, 3 or 4, i.e. that the quantity stored in *this
is in fact , or in the compactified external domain. The solution u with the boundary condition u =0 at spatial infinity is the returned Scalar
.
Definition at line 180 of file scalar_pde.C.
References Tensor::mp, and Map::poisson_tau().
Scalar Scalar::poisson_tau | ( | ) | const |
Solves the scalar Poisson equation with *this
as a source using a real Tau method The source of the equation is represented by the Scalar
*this
.
Note that dzpuis
must be equal to 2, 3 or 4, i.e. that the quantity stored in *this
is in fact , or in the compactified external domain. The solution u with the boundary condition u =0 at spatial infinity is the returned Scalar
.
Definition at line 165 of file scalar_pde.C.
References Tensor::mp, and Map::poisson_tau().
int Tensor::position | ( | const Itbl & | ind | ) | const [virtual, inherited] |
Returns the position in the array cmp
of a component given by its indices.
ind | [input] 1-D array of integers (class Itbl ) of size valence giving the values of each index specifing the component, with the following storage convention:
|
cmp
of the pointer to the Scalar
containing the component specified by ind
Reimplemented in Tensor_sym, and Vector.
Definition at line 521 of file tensor.C.
References Itbl::get_dim(), Itbl::get_ndim(), and Tensor::valence.
Scalar Scalar::primr | ( | bool | null_infty = true |
) | const |
Computes the radial primitive which vanishes for .
i.e. the function where f is the function represented by *this
(and must have a dzpuis
= 2).
null_infty | if true (default), the primitive is null at infinity (or on the grid boundary). F is null at the center otherwise |
Definition at line 97 of file scalar_integ.C.
References Tensor::mp, and Map::primr().
void Scalar::raccord | ( | int | n | ) |
Performs the matching of the nucleus with respect to the first shell.
Definition at line 58 of file scalar_raccord.C.
References Tbl::annule_hard(), Valeur::base, Valeur::c_cf, Valeur::coef(), etat, Map_af::get_alpha(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Mg3d::get_type_r(), Matrice::inverse(), Tensor::mp, R_CHEBI, R_CHEBP, Mtbl_cf::set(), Tbl::set(), Valeur::set_etat_cf_qcq(), Mtbl_cf::t, va, Valeur::ylm(), and Valeur::ylm_i().
void Scalar::raccord_c1_zec | ( | int | puis, | |
int | nbre, | |||
int | lmax | |||
) |
Performs the matching of the external domain with respect to the last shell using function like with for each spherical harmonics with .
Definition at line 73 of file scalar_raccord_zec.C.
References Tbl::annule_hard(), Valeur::base, Valeur::c_cf, Valeur::coef(), dec_dzpuis(), dsdr(), etat, Map_af::get_alpha(), Base_val::get_base_r(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Matrice::inverse(), Tensor::mp, R_CHEB, Mtbl_cf::set(), Matrice::set(), Tbl::set(), Matrice::set_band(), set_dzpuis(), Valeur::set_etat_cf_qcq(), Matrice::set_etat_qcq(), Tbl::set_etat_qcq(), Matrice::set_lu(), Mtbl_cf::t, va, Valeur::ylm(), and Valeur::ylm_i().
void Scalar::raccord_externe | ( | int | puis, | |
int | nbre, | |||
int | lmax | |||
) |
Matching of the external domain with the outermost shell.
Definition at line 62 of file scalar_raccord_externe.C.
References Tbl::annule_hard(), Valeur::base, Valeur::c_cf, Valeur::coef(), Map_af::get_alpha(), Map_af::get_beta(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, Mtbl_cf::set(), Tbl::set(), Matrice::set(), set_dzpuis(), Valeur::set_etat_cf_qcq(), Mtbl_cf::set_etat_qcq(), Tbl::set_etat_qcq(), Matrice::set_etat_qcq(), Mtbl_cf::t, va, Valeur::ylm(), and Valeur::ylm_i().
void Scalar::sarra_filter_r | ( | int | lzmin, | |
int | lzmax, | |||
double | p, | |||
double | alpha = -1E-16 | |||
) |
Applies an exponential filter to the spectral coefficients in the radial direction.
The filter is of the type: , with and N the number of radial coefficients.
lzmin,lzmax | [input] the indices of the domain where the filter is applied (in [lzmin , lzmax ]) | |
p | [input] the order of the filter | |
alpha | [input] appearing in the above formula. |
Definition at line 92 of file scalar_exp_filter.C.
References Valeur::c, Valeur::c_cf, Valeur::coef(), del_deriv(), Valeur::del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Tensor::mp, Mtbl_cf::set(), and va.
void Scalar::sarra_filter_r_all_domains | ( | double | p, | |
double | alpha = 1E-16 | |||
) |
Applies an exponential filter in radial direction in all domains for the case where p is a double (see Scalar:sarra_filter_r
).
Note that this may cause regularity problems at the origin if applied in a nucleus.
Definition at line 132 of file scalar_exp_filter.C.
References Map::get_mg(), Tensor::get_mp(), Mg3d::get_nzone(), and sarra_filter_r().
void Scalar::sauve | ( | FILE * | fd | ) | const [virtual] |
Save in a file.
Reimplemented from Tensor.
Definition at line 679 of file scalar.C.
References dzpuis, etat, fwrite_be(), Valeur::sauve(), and va.
Scalar Scalar::scalar_out_bound | ( | int | n, | |
bool | leave_ylm = false | |||
) |
Returns the Scalar
containing the values of angular coefficients at the outer boundary.
l_dom | [input] domain index | |
leave_ylm | [input] flag to decide whether the coefficients are expressed in spherical harmonics or Fourier base |
Definition at line 463 of file scalar_manip.C.
References Valeur::base, Valeur::c_cf, Valeur::coef(), etat, Base_val::get_base_t(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Tensor::mp, Map::mp_angu(), Base_val::set_base_t(), va, Mtbl_cf::val_out_bound_jk(), and Valeur::ylm().
Scalar & Tensor::set | ( | int | i1, | |
int | i2, | |||
int | i3, | |||
int | i4 | |||
) | [inherited] |
Returns the value of a component for a tensor of valence 4 (read/write version).
i1 | value of the first index (1, 2 or 3) | |
i2 | value of the second index (1, 2 or 3) | |
i3 | value of the third index (1, 2 or 3) | |
i4 | value of the fourth index (1, 2 or 3) |
(i1,i2,i3,i4) Definition at line 633 of file tensor.C.
References Tensor::cmp, Tensor::del_deriv(), Tensor::position(), Itbl::set(), and Tensor::valence.
Scalar & Tensor::set | ( | int | i1, | |
int | i2, | |||
int | i3 | |||
) | [inherited] |
Returns the value of a component for a tensor of valence 3 (read/write version).
i1 | value of the first index (1, 2 or 3) | |
i2 | value of the second index (1, 2 or 3) | |
i3 | value of the third index (1, 2 or 3) |
(i1,i2,i3) Definition at line 617 of file tensor.C.
References Tensor::cmp, Tensor::del_deriv(), Tensor::position(), Itbl::set(), and Tensor::valence.
Scalar & Tensor::set | ( | int | i1, | |
int | i2 | |||
) | [inherited] |
Returns the value of a component for a tensor of valence 2 (read/write version).
i1 | value of the first index (1, 2 or 3) | |
i2 | value of the second index (1, 2 or 3) |
(i1,i2) Definition at line 602 of file tensor.C.
References Tensor::cmp, Tensor::del_deriv(), Tensor::position(), Itbl::set(), and Tensor::valence.
Returns the value of a component (read/write version).
ind | 1-D Itbl of size valence containing the values of each index specifing the component, with the following storage convention:
|
ind
Definition at line 650 of file tensor.C.
References Tensor::cmp, Tensor::del_deriv(), Itbl::get_dim(), Itbl::get_ndim(), Tensor::position(), and Tensor::valence.
void Tensor::set_dependance | ( | const Metric & | met | ) | const [protected, inherited] |
To be used to describe the fact that the derivatives members have been calculated with met
.
First it sets a null element of met_depend
to &met
and puts this
in the list of the dependancies of met
.
Definition at line 449 of file tensor.C.
References Tensor::met_depend, and Metric::tensor_depend.
void Scalar::set_der_0x0 | ( | ) | const [protected] |
void Tensor::set_der_met_0x0 | ( | int | i | ) | const [protected, inherited] |
Sets all the i-th components of met_depend
, p_derive_cov
, etc.
.. to 0x0.
Reimplemented in Sym_tensor, and Vector.
Definition at line 429 of file tensor.C.
References Tensor::met_depend, Tensor::p_derive_con, Tensor::p_derive_cov, and Tensor::p_divergence.
Tbl& Scalar::set_domain | ( | int | l | ) | [inline] |
Read/write of the value in a given domain.
This method should be used only to set the value in a given domain (it performs a call to del_deriv
); for reading the value in a domain without changing it, the method domain(int )
is preferable.
l | [input] domain index |
Tbl
containing the value of the field in domain l
. Definition at line 597 of file scalar.h.
References del_deriv(), etat, Valeur::set(), and va.
void Scalar::set_dzpuis | ( | int | dzi | ) |
Modifies the dzpuis
flag.
NB: this method does not change the field values stored in the compactified external domain (use methods dec_dzpuis()
, etc... for this purpose).
Definition at line 801 of file scalar.C.
References dzpuis.
void Scalar::set_etat_nondef | ( | ) | [virtual] |
void Scalar::set_etat_one | ( | ) |
void Scalar::set_etat_qcq | ( | ) | [virtual] |
Sets the logical state to ETATQCQ
(ordinary state).
If the state is already ETATQCQ
, this function does nothing. Otherwise, it calls the logical destructor of the Valeur
va
and deallocates the memory occupied by all the derivatives.
Reimplemented from Tensor.
void Scalar::set_etat_zero | ( | ) | [virtual] |
Sets the logical state to ETATZERO
(zero).
Calls the logical destructor of the Valeur
va
and deallocates the memory occupied by all the derivatives.
Reimplemented from Tensor.
Definition at line 317 of file scalar.C.
References del_deriv(), etat, Valeur::set_etat_zero(), and va.
double& Scalar::set_grid_point | ( | int | l, | |
int | k, | |||
int | j, | |||
int | i | |||
) | [inline] |
Setting the value of the field at a given grid point.
CAUTION: to gain in efficiency (especially when this method is invoqued inside a loop), the method del_deriv()
(to delete the derived members) is not called by set_grid_point
. It must thus be invoqued by the user, after all the calls to set_grid_point
have been performed.
l | [input] domain index | |
k | [input] index | |
j | [input] index | |
i | [input] r () index |
Definition at line 666 of file scalar.h.
References etat, Valeur::set(), and va.
Itbl& Tensor::set_index_type | ( | ) | [inline, inherited] |
Sets the types of all the indices.
Itbl
) of size valence
that can be modified (COV
for a covariant one and CON
for a contravariant one) Definition at line 914 of file tensor.h.
References Tensor::type_indice.
int& Tensor::set_index_type | ( | int | i | ) | [inline, inherited] |
Sets the type of the index number i
.
i
must be strictly lower than valence
and obey the following convention:
i
= 0 : first index i
= 1 : second index COV
for a covariant index, CON
for a contravariant one) Definition at line 905 of file tensor.h.
References Itbl::set(), and Tensor::type_indice.
void Scalar::set_inner_boundary | ( | int | l, | |
double | x | |||
) |
Sets the value of the Scalar
at the inner boundary of a given domain.
l | [input] domain index | |
x | [input] (constant) value at the inner boundary of domain no. l |
Definition at line 279 of file scalar_manip.C.
References annule_hard(), Valeur::coef_i(), del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Tensor::mp, Valeur::set(), Valeur::set_etat_c_qcq(), and va.
void Scalar::set_outer_boundary | ( | int | l, | |
double | x | |||
) |
Sets the value of the Scalar
at the outer boundary of a given domain.
l | [input] domain index | |
x | [input] (constant) value at the outer boundary of domain no. l |
Definition at line 311 of file scalar_manip.C.
References annule_hard(), Valeur::coef_i(), del_deriv(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Tensor::mp, Valeur::set(), Valeur::set_etat_c_qcq(), and va.
void Scalar::set_spectral_base | ( | const Base_val & | bi | ) |
Sets the spectral bases of the Valeur
va
.
Definition at line 790 of file scalar.C.
References Valeur::set_base(), and va.
Valeur& Scalar::set_spectral_va | ( | ) | [inline] |
void Tensor::set_triad | ( | const Base_vect & | new_triad | ) | [inherited] |
Assigns a new vectorial basis (triad) of decomposition.
NB: this function modifies only the member triad
and leave unchanged the components (member cmp
). In order to change them coherently with the new basis, the function change_triad
(const Base_vect& ) must be called instead.
Definition at line 515 of file tensor.C.
References Tensor::triad.
void Scalar::smooth_decay | ( | int | k, | |
int | n | |||
) |
Performs a matching of the last non-compactified shell with a decaying function where is the spherical harmonic index and n is some specifiable parameter.
Definition at line 209 of file scalar_raccord_zec.C.
References Valeur::c, Valeur::c_cf, Valeur::coef(), dsdr(), etat, Valeur::get_base(), get_dzpuis(), get_etat(), Valeur::get_etat(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), Mg3d::get_type_r(), Matrice::inverse(), Tensor::mp, Map::r, Matrice::set(), Tbl::set(), Matrice::set_band(), Valeur::set_base(), Valeur::set_base_r(), Valeur::set_etat_cf_qcq(), Matrice::set_etat_qcq(), Mtbl_cf::set_etat_qcq(), set_etat_qcq(), Tbl::set_etat_qcq(), Tbl::set_etat_zero(), Matrice::set_lu(), std_spectral_base(), Mtbl_cf::t, va, val_grid_point(), Valeur::ylm(), and Valeur::ylm_i().
Scalar Scalar::sol_divergence | ( | int | n | ) | const |
Resolution of a divergence-like equation.
The equation solved reads: , with the unknown and the source represented by this
.
n | [input] the coefficient in front of the 1/r term. |
Definition at line 64 of file scalar_sol_div.C.
References Tbl::annule_hard(), Mtbl_cf::annule_hard(), Valeur::c_cf, check_dzpuis(), etat, Map_af::get_alpha(), Map_af::get_beta(), Diff_sx::get_matrice(), Diff_dsdx::get_matrice(), Diff_xdsdx::get_matrice(), Mtbl_cf::get_mg(), Map::get_mg(), Mg3d::get_np(), Mg3d::get_nr(), Mg3d::get_nt(), Mg3d::get_nzone(), get_spectral_base(), Mg3d::get_type_r(), Base_val::give_quant_numbers(), Matrice::inverse(), Tensor::mp, Base_val::mult_x(), R_CHEB, R_CHEBI, R_CHEBP, R_CHEBU, Mtbl_cf::set(), Matrice::set(), Tbl::set(), Valeur::set_etat_cf_qcq(), Tbl::set_etat_qcq(), Matrice::set_etat_qcq(), set_etat_qcq(), set_etat_zero(), Matrice::set_lu(), set_spectral_base(), set_spectral_va(), va, and Valeur::ylm_i().
Scalar Scalar::sol_elliptic | ( | Param_elliptic & | params | ) | const |
Resolution of a general elliptic equation, putting zero at infinity.
params | [input] the operators and variables to be used. |
Definition at line 230 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic().
Scalar Scalar::sol_elliptic_2d | ( | Param_elliptic & | ope_var | ) | const |
Solves the scalar 2-dimensional elliptic equation with *this
as a source.
Note that dzpuis
must be equal to 2, 3 or 4, i.e. The solution u with the boundary condition u =0 at spatial infinity is the returned Scalar
.
Definition at line 405 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_2d().
Scalar Scalar::sol_elliptic_boundary | ( | Param_elliptic & | params, | |
const Scalar & | bound, | |||
double | fact_dir, | |||
double | fact_neu | |||
) | const |
Resolution of general elliptic equation, with inner boundary conditions as Scalars on mono-domain angulare grids.
Definition at line 278 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_boundary().
Scalar Scalar::sol_elliptic_boundary | ( | Param_elliptic & | params, | |
const Mtbl_cf & | bound, | |||
double | fact_dir, | |||
double | fact_neu | |||
) | const |
Resolution of a general elliptic equation, putting zero at infinity and with inner boundary conditions.
params | [input] the operators and variables to be used. | |
bound | [input] : the boundary condition | |
fact_dir | : 1 Dirchlet condition, 0 Neumann condition | |
fact_neu | : 0 Dirchlet condition, 1 Neumann condition |
Definition at line 252 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_boundary().
Scalar Scalar::sol_elliptic_fixe_der_zero | ( | double | val, | |
Param_elliptic & | params | |||
) | const |
Resolution of a general elliptic equation fixing the dericative at the origin and relaxing one continuity condition.
val | [input] value of the derivative. | |
params | [input] the operators and variables to be used. |
Definition at line 382 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_fixe_der_zero().
Scalar Scalar::sol_elliptic_no_zec | ( | Param_elliptic & | params, | |
double | val = 0 | |||
) | const |
Resolution of a general elliptic equation, putting a given value at the outermost shell and not solving in the compactified domain.
params | [input] the operators and variables to be used. | |
val | [input] value at the last shell. |
Definition at line 310 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_no_zec().
Scalar Scalar::sol_elliptic_only_zec | ( | Param_elliptic & | params, | |
double | val | |||
) | const |
Resolution of a general elliptic equation solving in the compactified domain and putting a given value at the inner boundary.
params | [input] the operators and variables to be used. | |
val | [input] value at the inner boundary of the external domain. |
Definition at line 337 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_only_zec().
Scalar Scalar::sol_elliptic_pseudo_1d | ( | Param_elliptic & | ope_var | ) | const |
Solves a pseudo-1d elliptic equation with *this
as a source.
Definition at line 426 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_pseudo_1d().
Scalar Scalar::sol_elliptic_sin_zec | ( | Param_elliptic & | params, | |
double * | coefs, | |||
double * | phases | |||
) | const |
General elliptic solver.
The equation is not solved in the compactified domain and the matching is done with an homogeneous solution.
params | [input] the operators and variables to be used. | |
coef | [output] : coefficient of the oscillatory solution in the external domain. | |
phases | [output] : phases (i.e. choice of the homogeneous solution to match with). |
Definition at line 359 of file scalar_pde.C.
References Tensor::mp, set_etat_qcq(), and Map_af::sol_elliptic_sin_zec().
void Scalar::spectral_display | ( | const char * | comment = 0x0 , |
|
double | threshold = 1.e-7 , |
|||
int | precision = 4 , |
|||
ostream & | ostr = cout | |||
) | const [virtual] |
Displays the spectral coefficients and the associated basis functions.
This function shows only the values greater than a given threshold.
comment | comment to be printed at top of the display (default: 0x0 = nothing printed) | |
threshold | [input] Value above which a coefficient is printed (default: 1.e-7) | |
precision | [input] Number of printed digits (default: 4) | |
ostr | [input] Output stream used for the printing (default: cout) |
Reimplemented from Tensor.
Definition at line 734 of file scalar.C.
References Valeur::display_coef(), dzpuis, etat, and va.
const Scalar & Scalar::srdsdt | ( | ) | const |
Returns of *this
.
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 144 of file scalar_deriv.C.
References dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_srdsdt, Scalar(), set_dzpuis(), set_etat_zero(), and Map::srdsdt().
const Scalar & Scalar::srstdsdp | ( | ) | const |
Returns of *this
.
If dzpuis
is zero, then the returned Scalar
has dzpuis
= 2. It is increased by 1 otherwise.
Definition at line 176 of file scalar_deriv.C.
References dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Mg3d::get_type_r(), Tensor::mp, p_srstdsdp, Scalar(), set_dzpuis(), set_etat_zero(), and Map::srstdsdp().
void Scalar::std_spectral_base | ( | ) | [virtual] |
void Scalar::std_spectral_base_odd | ( | ) | [virtual] |
const Scalar & Scalar::stdsdp | ( | ) | const |
Returns of *this
.
Definition at line 237 of file scalar_deriv.C.
References dzpuis, etat, Tensor::mp, p_stdsdp, Scalar(), set_dzpuis(), set_etat_zero(), and Map::stdsdp().
Tbl Scalar::tbl_in_bound | ( | int | n, | |
bool | leave_ylm = false | |||
) |
Returns the Tbl
containing the values of angular coefficients at the inner boundary.
l_dom | [input] domain index | |
leave_ylm | [input] flag to decide whether the coefficients are expressed in spherical harmonics or Fourier base |
Definition at line 444 of file scalar_manip.C.
References Valeur::c_cf, Valeur::coef(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Mg3d::get_type_r(), Tensor::mp, Tbl::set(), Tbl::set_etat_qcq(), Tbl::set_etat_zero(), va, Mtbl_cf::val_in_bound_jk(), and Valeur::ylm().
Tbl Scalar::tbl_out_bound | ( | int | l_dom, | |
bool | leave_ylm = false | |||
) |
Returns the Tbl
containing the values of angular coefficients at the outer boundary.
l_dom | [input] domain index | |
leave_ylm | [input] flag to decide whether the coefficients are expressed in spherical harmonics or Fourier base |
Definition at line 426 of file scalar_manip.C.
References Valeur::c_cf, Valeur::coef(), etat, Map::get_mg(), Mg3d::get_np(), Mg3d::get_nt(), Tensor::mp, Tbl::set(), Tbl::set_etat_qcq(), Tbl::set_etat_zero(), va, Mtbl_cf::val_out_bound_jk(), and Valeur::ylm().
Checks if a Poisson equation with *this
as a source has been correctly solved.
uu | [input] Solution u of the Poisson equation , being represented by the Scalar *this . | |
ostr | [input/output] Output stream used for displaying err . | |
detail | [input]
|
Tbl
err
decribing the errors in each domain: err(0,l)
: Relative error in domain no. l
, defined as the maximum value of in that domain divided by M , where M is the maximum value of over all domains if dzpuis
= 0 or is zero in the compactified external domain (CED). If dzpuis
!= 0 and does not vanish in the CED, the value of M used in the non-compactified domains is the maximum value over these domains, whereas the value of M used in the compactified external domain is the maximum value on that particular domain. err(1,l)
: Maximum value of the absolute error in domain no. l
err(2,l)
: Maximum value of in domain no. l
Definition at line 56 of file scalar_test_poisson.C.
References check_dzpuis(), dzpuis, Map::get_mg(), Tensor::get_mp(), Mg3d::get_nzone(), laplacian(), Tensor::mp, Tbl::set(), and Tbl::set_etat_qcq().
Trace with respect to a given metric for a valence 2 tensor.
gam | metric used to raise or lower one of the indices, in order to take the trace |
Definition at line 193 of file tensor_calculus.C.
References Metric::con(), contract(), Metric::cov(), Tensor::trace(), Tensor::type_indice, and Tensor::valence.
Scalar Tensor::trace | ( | ) | const [inherited] |
Trace on two different type indices for a valence 2 tensor.
Definition at line 176 of file tensor_calculus.C.
References Tensor::mp, Tensor::operator()(), set_etat_zero(), Tensor::type_indice, and Tensor::valence.
Trace with respect to a given metric.
ind1 | first index for the contraction, with the following convention :
| |
ind2 | second index for the contraction | |
gam | metric used to raise or lower ind1 in order that it has a opposite type than ind2 |
Definition at line 149 of file tensor_calculus.C.
References Metric::con(), contract(), Metric::cov(), Tensor::trace(), Tensor::type_indice, and Tensor::valence.
Tensor Tensor::trace | ( | int | ind1, | |
int | ind2 | |||
) | const [inherited] |
Trace on two different type indices.
ind1 | first index for the contraction, with the following convention :
| |
ind2 | second index for the contraction |
Definition at line 90 of file tensor_calculus.C.
References Tensor::cmp, Tensor::get_n_comp(), Tensor::indices(), Tensor::mp, Tensor::position(), Tensor::set(), Itbl::set(), set_etat_zero(), Tensor::triad, Tensor::type_indice, and Tensor::valence.
Computes a new tensor by raising an index of *this
.
ind | index to be raised, with the following convention :
| |
gam | metric used to raise the index (contraction with the twice contravariant form of the metric on the index ind ). |
Definition at line 221 of file tensor_calculus.C.
References Metric::con(), contract(), Tensor::indices(), Tensor::mp, Tensor::n_comp, Tensor::set(), Itbl::set(), Tensor::triad, Tensor::type_indice, and Tensor::valence.
Computes a new tensor by raising or lowering all the indices of *this
.
gam | metric used to lower the contravariant indices and raising the covariant ones. |
Definition at line 301 of file tensor_calculus.C.
References Tensor::down(), Tensor::Tensor(), Tensor::type_indice, Tensor::up(), and Tensor::valence.
double Scalar::val_grid_point | ( | int | l, | |
int | k, | |||
int | j, | |||
int | i | |||
) | const [inline] |
double Scalar::val_point | ( | double | r, | |
double | theta, | |||
double | phi | |||
) | const |
Computes the value of the field at an arbitrary point , by means of the spectral expansion.
NB: if is a point of the spectral grid, the method val_grid_point
is to be preferred, being much more efficient.
r | [input] value of the coordinate r | |
theta | [input] value of the coordinate | |
phi | [input] value of the coordinate |
*this
. NB: in the compactified external domain, the returned value is the actual value of the field, i.e. the stored value divided by . Definition at line 883 of file scalar.C.
References dzpuis, etat, Map::get_mg(), Mg3d::get_nzone(), Tensor::mp, va, Map::val_lx(), and Valeur::val_point().
void Scalar::visu_box | ( | double | xmin, | |
double | xmax, | |||
double | ymin, | |||
double | ymax, | |||
double | zmin, | |||
double | zmax, | |||
const char * | title0 = 0x0 , |
|||
const char * | filename0 = 0x0 , |
|||
bool | start_dx = true , |
|||
int | nx = 40 , |
|||
int | ny = 40 , |
|||
int | nz = 40 | |||
) | const |
3D visualization (volume rendering) via OpenDX.
Prepares files for visualization by OpenDX of the values of the field in some rectangular box.
xmin | [input] defines with xmax the x range of the visualization box | |
xmax | [input] defines with xmin the x range of the visualization box | |
ymin | [input] defines with ymax the y range of the visualization box | |
ymax | [input] defines with ymin the y range of the visualization box | |
zmin | [input] defines with zmax the z range of the visualization box | |
zmax | [input] defines with zmin the z range of the visualization box | |
title | [input] title for the graph (for OpenDX legend) | |
filename | [input] name for the file which will be the input for OpenDX; the default 0x0 is transformed into "scalar_box" | |
start_dx | [input] determines whether OpenDX must be launched (as a subprocess) to view the field; if set to false , only input files for future usage of OpenDX are created | |
nx | [input] number of points in the x direction (uniform sampling) | |
ny | [input] number of points in the y direction (uniform sampling) | |
nz | [input] number of points in the z direction (uniform sampling) |
Definition at line 338 of file scalar_visu.C.
References Valeur::c_cf, check_dzpuis(), Valeur::coef(), Map::convert_absolute(), dec_dzpuis(), dzpuis, Tensor::mp, Scalar(), va, Map::val_lx(), and Mtbl_cf::val_point().
void Scalar::visu_section | ( | const Tbl & | plane, | |
double | umin, | |||
double | umax, | |||
double | vmin, | |||
double | vmax, | |||
const char * | title = 0x0 , |
|||
const char * | filename = 0x0 , |
|||
bool | start_dx = true , |
|||
int | nu = 200 , |
|||
int | nv = 200 | |||
) | const |
3D visualization via a plane section.
Prepares files for visualization by OpenDX of the values of the field in any given plane.
plane | [input] : 2D Tbl defining the section plane: plane must of dimension 3x3 with the following content:
| |
umin | [input] defines with umax the range of the plane coordinate u | |
umax | [input] defines with umin the range of the plane coordinate u | |
vmin | [input] defines with vmax the range of the plane coordinate v | |
vmax | [input] defines with vmin the range of the plane coordinate v | |
title | [input] title for the graph (for OpenDX legend) | |
filename | [input] name for the file which will be the input for OpenDX; the default 0x0 is transformed into "scalar_section" | |
start_dx | [input] determines whether OpenDX must be launched (as a subprocess) to view the field; if set to false , only input files for future usage of OpenDX are created | |
nu | [input] number of points in the u direction (uniform sampling) | |
nv | [input] number of points in the v direction (uniform sampling) |
Definition at line 152 of file scalar_visu.C.
References Valeur::c_cf, Valeur::coef(), Map::convert_absolute(), Tensor::mp, va, Map::val_lx(), and Mtbl_cf::val_point().
void Scalar::visu_section | ( | const char | section_type, | |
double | aa, | |||
double | umin, | |||
double | umax, | |||
double | vmin, | |||
double | vmax, | |||
const char * | title = 0x0 , |
|||
const char * | filename = 0x0 , |
|||
bool | start_dx = true , |
|||
int | nu = 200 , |
|||
int | nv = 200 | |||
) | const |
3D visualization via a plane section.
Prepares files for visualization by OpenDX of the values of the field in a plane x=const, y=const or z=const
section_type | [input] defines the type of section :
| |
aa | [input] constant a defining the section plane | |
umin | [input] defines with umax the range of the plane coordinate u | |
umax | [input] defines with umin the range of the plane coordinate u | |
vmin | [input] defines with vmax the range of the plane coordinate v | |
vmax | [input] defines with vmin the range of the plane coordinate v | |
title | [input] title for the graph (for OpenDX legend) | |
filename | [input] name for the file which will be the input for OpenDX; the default 0x0 is transformed into "scalar_section" | |
start_dx | [input] determines whether OpenDX must be launched (as a subprocess) to view the field; if set to false , only input files for future usage of OpenDX are created | |
nu | [input] number of points in the u direction (uniform sampling) | |
nv | [input] number of points in the v direction (uniform sampling) |
Definition at line 74 of file scalar_visu.C.
References Tbl::set(), and Tbl::set_etat_qcq().
void Scalar::visu_section_anim | ( | const char | section_type, | |
double | aa, | |||
double | umin, | |||
double | umax, | |||
double | vmin, | |||
double | vmax, | |||
int | jtime, | |||
double | ttime, | |||
int | jgraph = 1 , |
|||
const char * | title = 0x0 , |
|||
const char * | filename_root = 0x0 , |
|||
bool | start_dx = false , |
|||
int | nu = 200 , |
|||
int | nv = 200 | |||
) | const |
3D visualization via time evolving plane section (animation).
Prepares files for visualization by OpenDX of the values of the field in a plane x=const, y=const or z=const at successive time steps
section_type | [input] defines the type of section :
| |
aa | [input] constant a defining the section plane | |
umin | [input] defines with umax the range of the plane coordinate u | |
umax | [input] defines with umin the range of the plane coordinate u | |
vmin | [input] defines with vmax the range of the plane coordinate v | |
vmax | [input] defines with vmin the range of the plane coordinate v | |
jtime | [input] time step label | |
ttime | [input] time t corresponding to jtime | |
jgraph | [input] number of time steps between two graphs: the graph will be generated only if jtime is a multiple of jgraph | |
title | [input] title for the graph (for OpenDX legend) | |
filename_root | [input] beginning of the names for the files which will be the input for OpenDX (the end of names will be automatically generated from the time steps); the default 0x0 is transformed into "anim" | |
start_dx | [input] determines whether OpenDX must be launched (as a subprocess) to view the field; if set to false , only input files for future usage of OpenDX are created | |
nu | [input] number of points in the u direction (uniform sampling) | |
nv | [input] number of points in the v direction (uniform sampling) |
Definition at line 531 of file scalar_visu.C.
References visu_section().
Tensor_sym operator* | ( | const Tensor_sym & | a, | |
const Tensor_sym & | b | |||
) | [friend, inherited] |
Tensorial product of two symmetric tensors.
NB: the output is an object of class Tensor_sym
, with the two symmetric indices corresponding to the symmetric indices of tensor a
. This means that the symmetries of tensor b
indices are not used in the storage, since there is currently no class in Lorene to manage tensors with more than two symmetric indices.
Definition at line 147 of file tensor_sym_calculus.C.
Scalar** Tensor::cmp [protected, inherited] |
int Scalar::dzpuis [protected] |
int Scalar::etat [protected] |
int Scalar::ind_lap [mutable, protected] |
const Metric* Tensor::met_depend[N_MET_MAX] [mutable, protected, inherited] |
const Map* const Tensor::mp [protected, inherited] |
int Tensor::n_comp [protected, inherited] |
Tensor* Tensor::p_derive_con[N_MET_MAX] [mutable, protected, inherited] |
Array of pointers on the contravariant derivatives of this
with respect to various metrics.
See the comments of met_depend
. See also the comments of method derive_con()
for a precise definition of a "contravariant" derivative.
Tensor* Tensor::p_derive_cov[N_MET_MAX] [mutable, protected, inherited] |
Array of pointers on the covariant derivatives of this
with respect to various metrics.
See the comments of met_depend
. See also the comments of method derive_cov()
for the index convention of the covariant derivation.
Tensor* Tensor::p_divergence[N_MET_MAX] [mutable, protected, inherited] |
Array of pointers on the divergence of this
with respect to various metrics.
See the comments of met_depend
. See also the comments of method divergence()
for a precise definition of a the divergence with respect to a given metric.
Scalar* Scalar::p_dsdr [mutable, protected] |
Scalar* Scalar::p_dsdradial [mutable, protected] |
Scalar* Scalar::p_dsdrho [mutable, protected] |
Scalar* Scalar::p_dsdt [mutable, protected] |
Scalar* Scalar::p_dsdx [mutable, protected] |
Scalar* Scalar::p_dsdy [mutable, protected] |
Scalar* Scalar::p_dsdz [mutable, protected] |
Tbl* Scalar::p_integ [mutable, protected] |
Scalar* Scalar::p_lap [mutable, protected] |
Scalar* Scalar::p_lapang [mutable, protected] |
Scalar* Scalar::p_srdsdt [mutable, protected] |
Scalar* Scalar::p_srstdsdp [mutable, protected] |
Scalar* Scalar::p_stdsdp [mutable, protected] |
const Base_vect* Tensor::triad [protected, inherited] |
Itbl Tensor::type_indice [protected, inherited] |
Valeur Scalar::va [protected] |
int Tensor::valence [protected, inherited] |