LORENE
Lorene::Tenseur_sym Class Reference

Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class Tensor_sym instead ***. More...

#include <tenseur.h>

Inheritance diagram for Lorene::Tenseur_sym:
Lorene::Tenseur

Public Member Functions

 Tenseur_sym (const Map &map, int val, const Itbl &tipe, const Base_vect &triad_i, const Metrique *met=0x0, double weight=0)
 Standard constructor. More...
 
 Tenseur_sym (const Map &map, int val, int tipe, const Base_vect &triad_i, const Metrique *met=0x0, double weight=0)
 Standard constructor when all the indices are of the same type. More...
 
 Tenseur_sym (const Tenseur_sym &)
 Copy constructor. More...
 
 Tenseur_sym (const Tenseur &)
 Constructor from a Tenseur . More...
 
 Tenseur_sym (const Map &map, const Base_vect &triad_i, FILE *fich, const Metrique *met=0x0)
 Constructor from a file (see sauve(FILE*) ). More...
 
virtual ~Tenseur_sym ()
 Destructor. More...
 
virtual void operator= (const Tenseur &)
 Assignment from a Tenseur . More...
 
void operator= (const Tenseur_sym &)
 
virtual int donne_place (const Itbl &idx) const
 Returns the position in the Cmp 1-D array c of a component given by its indices. More...
 
virtual Itbl donne_indices (int place) const
 Returns the indices of a component given by its position in the Cmp 1-D array c . More...
 
void set_etat_nondef ()
 Sets the logical state to ETATNONDEF (undefined state). More...
 
void set_etat_zero ()
 Sets the logical state to ETATZERO (zero state). More...
 
void set_etat_qcq ()
 Sets the logical state to ETATQCQ (ordinary state). More...
 
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. More...
 
void change_triad (const Base_vect &new_triad)
 Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly. More...
 
void set_triad (const Base_vect &new_triad)
 Assigns a new vectorial basis (triad) of decomposition. More...
 
void set_poids (double weight)
 Sets the weight for a tensor density. More...
 
void set_metric (const Metrique &met)
 Sets the pointer on the metric for a tensor density. More...
 
Cmpset ()
 Read/write for a scalar (see also operator=(const Cmp&) ). More...
 
Cmpset (int)
 Read/write for a vector. More...
 
Cmpset (int, int)
 Read/write for a tensor of valence 2. More...
 
Cmpset (int, int, int)
 Read/write for a tensor of valence 3. More...
 
Cmpset (const Itbl &)
 Read/write in the general case. More...
 
void annule (int l)
 Sets the Tenseur to zero in a given domain. More...
 
void annule (int l_min, int l_max)
 Sets the Tenseur to zero in several domains. More...
 
void set_std_base ()
 Set the standard spectal basis of decomposition for each component. More...
 
void dec_dzpuis ()
 dzpuis -= 1 ; More...
 
void inc_dzpuis ()
 dzpuis += 1 ; More...
 
void dec2_dzpuis ()
 dzpuis -= 2 ; More...
 
void inc2_dzpuis ()
 dzpuis += 2 ; More...
 
void mult_r_zec ()
 Multiplication by r in the external zone. More...
 
Tenseur inverse_poisson_vect (double lambda) const
 Compute $\Delta + \lambda \nabla\nabla$ of *this , *this being of valence 1. More...
 
const Mapget_mp () const
 Returns pointer on the mapping. More...
 
const Base_vectget_triad () const
 Returns the vectorial basis (triad) on which the components are defined. More...
 
int get_etat () const
 Returns the logical state. More...
 
int get_valence () const
 Returns the valence. More...
 
int get_n_comp () const
 Returns the number of components. More...
 
int get_type_indice (int i) const
 Returns the type of the index number i . More...
 
Itbl get_type_indice () const
 Returns the types of all the indices. More...
 
double get_poids () const
 Returns the weight. More...
 
const Metrique * get_metric () const
 Returns a pointer on the metric defining the conformal factor for tensor densities. More...
 
const Cmpoperator() () const
 Read only for a scalar. More...
 
const Cmpoperator() (int) const
 Read only for a vector. More...
 
const Cmpoperator() (int, int) const
 Read only for a tensor of valence 2. More...
 
const Cmpoperator() (int, int, int) const
 Read only for a tensor of valence 3. More...
 
const Cmpoperator() (const Itbl &) const
 Read only in the general case. More...
 
void sauve (FILE *) const
 Save in a file. More...
 
const Tenseurgradient () const
 Returns the gradient of *this (Cartesian coordinates) More...
 
const Tenseurgradient_spher () const
 Returns the gradient of *this (Spherical coordinates) (scalar field only). More...
 
const Tenseurderive_cov (const Metrique &met) const
 Returns the covariant derivative of *this , with respect to met . More...
 
const Tenseurderive_con (const Metrique &) const
 Returns the contravariant derivative of *this , with respect to met . More...
 
const Tenseurcarre_scal (const Metrique &) const
 Returns the scalar square of *this , with respect to met . More...
 
void poisson_vect (double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal) const
 Solves the vectorial Poisson equation : $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$. More...
 
Tenseur poisson_vect (double lambda, Tenseur &vect, Tenseur &scal) const
 Solves the vectorial Poisson equation $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$. More...
 
void poisson_vect_tau (double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal) const
 
Tenseur poisson_vect_tau (double lambda, Tenseur &vect, Tenseur &scal) const
 
void poisson_vect_falloff (double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal, int *k_falloff) const
 
Tenseur poisson_vect_falloff (double lambda, Tenseur &vect, Tenseur &scal, int *k_falloff) const
 
void poisson_vect_ylm (double lambda, Param &para, Tenseur &shift, Tenseur &vecteur, Tenseur &scalaire, int nylm, double *intvec) const
 
Tenseur poisson_vect_ylm (double lambda, Tenseur &vecteur, Tenseur &scalaire, int nylm, double *intvec) const
 
void poisson_vect_oohara (double lambda, Param &par, Tenseur &shift, Tenseur &scal) const
 Solves the vectorial Poisson equation $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$. More...
 
Tenseur poisson_vect_oohara (double lambda, Tenseur &scal) const
 Solves the vectorial Poisson equation $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$. More...
 
void poisson_vect_oohara_tau (double lambda, Param &par, Tenseur &shift, Tenseur &scal) const
 
Tenseur poisson_vect_oohara_tau (double lambda, Tenseur &scal) const
 
void poisson_vect_regu (int k_div, int nzet, double unsgam1, double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal) const
 Solves the vectorial Poisson equation : $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$. More...
 

Protected Member Functions

virtual void fait_gradient () const
 Calculates, if needed, the gradient of *this . More...
 
virtual void fait_derive_cov (const Metrique &met, int i) const
 Calculates, if needed, the covariant derivative of *this , with respect to met . More...
 
virtual void fait_derive_con (const Metrique &, int i) const
 Calculates, if needed, the contravariant derivative of *this , with respect to met . More...
 
bool verif () const
 Returns false for a tensor density without a defined metric. More...
 
void new_der_met ()
 Builds the arrays met_depend , p_derive_cov , p_derive_con and p_carre_scal and fills them with null pointers. More...
 
void del_t ()
 Logical destructor. More...
 
void del_derive_met (int i) const
 Logical destructor of the derivatives depending on the i-th element of *met_depend . More...
 
void del_derive () const
 Logical destructor of all the derivatives. More...
 
void set_der_met_0x0 (int i) const
 Sets the pointers of the derivatives depending on the i-th element of *met_depend to zero (as well as that i-th element). More...
 
void set_der_0x0 () const
 Sets the pointers of all the derivatives to zero. More...
 
void fait_gradient_spher () const
 Calculates, if needed, the gradient of *this in a spherical orthonormal basis (scalar field only). More...
 
void fait_carre_scal (const Metrique &, int i) const
 Calculates, if needed, the scalar square of *this , with respect to met . More...
 
void set_dependance (const Metrique &met) const
 To be used to describe the fact that the derivatives members have been calculated with met . More...
 
int get_place_met (const Metrique &metre) const
 Returns the position of the pointer on metre in the array met_depend . More...
 

Protected Attributes

const Map *const mp
 Reference mapping. More...
 
int valence
 Valence. More...
 
const Base_vecttriad
 Vectorial basis (triad) with respect to which the tensor components are defined. More...
 
Itbl type_indice
 Array of size valence contening the type of each index, COV for a covariant one and CON for a contravariant one. More...
 
int n_comp
 Number of components, depending on the symmetry. More...
 
int etat
 Logical state ETATZERO , ETATQCQ or ETATNONDEF. More...
 
Cmp ** c
 The components. More...
 
double poids
 For tensor densities: the weight. More...
 
const Metrique * metric
 For tensor densities: the metric defining the conformal factor. More...
 
const Metrique ** met_depend
 Array of pointers on the Metrique 's used to calculate derivatives members. More...
 
Tenseurp_gradient
 Pointer on the gradient of *this . More...
 
Tenseurp_gradient_spher
 Pointer on the gradient of *this in a spherical orthonormal basis (scalar field only). More...
 
Tenseur ** p_derive_cov
 Array of pointers on the covariant derivatives of *this
with respect to the corresponding metric in *met_depend . More...
 
Tenseur ** p_derive_con
 Array of pointers on the contravariant derivatives of *this
with respect to the corresponding metric in *met_depend . More...
 
Tenseur ** p_carre_scal
 Array of pointers on the scalar squares of *this
with respect to the corresponding metric in *met_depend . More...
 

Friends

Tenseur_sym operator* (const Tenseur &, const Tenseur_sym &)
 Tensorial product. More...
 
Tenseur_sym manipule (const Tenseur_sym &, const Metrique &)
 Raise or lower all the indices, depending on their type, using the given Metrique . More...
 
Tenseur lie_derive (const Tenseur &, const Tenseur &, const Metrique *)
 Lie Derivative of t with respect to x . More...
 

Detailed Description

Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class Tensor_sym instead ***.

The storage and the calculations are different and quicker than with an usual Tenseur .()

The valence must be >1.

Definition at line 1256 of file tenseur.h.

Constructor & Destructor Documentation

◆ Tenseur_sym() [1/5]

Lorene::Tenseur_sym::Tenseur_sym ( const Map map,
int  val,
const Itbl tipe,
const Base_vect triad_i,
const Metrique *  met = 0x0,
double  weight = 0 
)

Standard constructor.

Parameters
mapthe mapping
valvalence of the tensor; must be greater or equal to 2.
tipe1-D Itbl of size valence containing the type of each index, COV for a covariant one and CON for a contravariant one, with the following storage convention:
  • tipe(0) : type of the first index
  • tipe(1) : type of the second index
  • and so on...
triad_ivectorial basis (triad) with respect to which the tensor components are defined

Definition at line 105 of file tenseur_sym.C.

◆ Tenseur_sym() [2/5]

Lorene::Tenseur_sym::Tenseur_sym ( const Map map,
int  val,
int  tipe,
const Base_vect triad_i,
const Metrique *  met = 0x0,
double  weight = 0 
)

Standard constructor when all the indices are of the same type.

Parameters
mapthe mapping
valvalence of the tensor; must be greater or equal to 2.
tipethe type of the indices.
triad_ivectorial basis (triad) with respect to which the tensor components are defined

Definition at line 116 of file tenseur_sym.C.

◆ Tenseur_sym() [3/5]

◆ Tenseur_sym() [4/5]

◆ Tenseur_sym() [5/5]

Lorene::Tenseur_sym::Tenseur_sym ( const Map map,
const Base_vect triad_i,
FILE *  fich,
const Metrique *  met = 0x0 
)

Constructor from a file (see sauve(FILE*) ).

Parameters
mapthe mapping
triad_ivectorial basis (triad) with respect to which the tensor components are defined. It will be checked that it coincides with the basis saved in the file.
fichfile which has been created by the function sauve(FILE*) .

Definition at line 198 of file tenseur_sym.C.

References Lorene::Tenseur::n_comp, Lorene::pow(), and Lorene::Tenseur::valence.

◆ ~Tenseur_sym()

Lorene::Tenseur_sym::~Tenseur_sym ( )
virtual

Destructor.

Definition at line 210 of file tenseur_sym.C.

Member Function Documentation

◆ allocate_all()

void Lorene::Tenseur::allocate_all ( )
inherited

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 Tenseur -> Cmp -> Valeur -> Mtbl -> Tbl .

Definition at line 673 of file tenseur.C.

References Lorene::Cmp::allocate_all(), Lorene::Tenseur::c, Lorene::Tenseur::n_comp, and Lorene::Tenseur::set_etat_qcq().

◆ annule() [1/2]

void Lorene::Tenseur::annule ( int  l)
inherited

Sets the Tenseur to zero in a given domain.

Parameters
l[input] Index of the domain in which the Tenseur will be set (logically) to zero.

Definition at line 916 of file tenseur.C.

◆ annule() [2/2]

void Lorene::Tenseur::annule ( int  l_min,
int  l_max 
)
inherited

Sets the Tenseur to zero in several domains.

Parameters
l_min[input] The Tenseur 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 \c annule(0,nz-1), where \c nz  is the total number
of domains, is equivalent to set_etat_zero() .

Definition at line 921 of file tenseur.C.

References Lorene::Tenseur::etat, Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), Lorene::Tenseur::mp, and Lorene::Tenseur::set_etat_zero().

◆ carre_scal()

const Tenseur & Lorene::Tenseur::carre_scal ( const Metrique &  metre) const
inherited

Returns the scalar square of *this , with respect to met .

Definition at line 1593 of file tenseur.C.

References Lorene::Tenseur::fait_carre_scal(), Lorene::Tenseur::get_place_met(), Lorene::Tenseur::p_carre_scal, and Lorene::Tenseur::set_dependance().

◆ change_triad()

void Lorene::Tenseur::change_triad ( const Base_vect new_triad)
inherited

Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.

Definition at line 684 of file tenseur.C.

References Lorene::Base_vect::change_basis().

◆ dec2_dzpuis()

void Lorene::Tenseur::dec2_dzpuis ( )
inherited

dzpuis -= 2 ;

Definition at line 1146 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ dec_dzpuis()

void Lorene::Tenseur::dec_dzpuis ( )
inherited

dzpuis -= 1 ;

Definition at line 1120 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ del_derive()

void Lorene::Tenseur::del_derive ( ) const
protectedinherited

Logical destructor of all the derivatives.

Definition at line 589 of file tenseur.C.

References Lorene::Tenseur::del_derive_met(), Lorene::Tenseur::p_gradient, Lorene::Tenseur::p_gradient_spher, and Lorene::Tenseur::set_der_0x0().

◆ del_derive_met()

void Lorene::Tenseur::del_derive_met ( int  i) const
protectedinherited

Logical destructor of the derivatives depending on the i-th element of *met_depend .

Definition at line 570 of file tenseur.C.

References Lorene::Tenseur::met_depend, Lorene::Tenseur::p_carre_scal, Lorene::Tenseur::p_derive_con, Lorene::Tenseur::p_derive_cov, and Lorene::Tenseur::set_der_met_0x0().

◆ del_t()

void Lorene::Tenseur::del_t ( )
protectedinherited

Logical destructor.

Definition at line 561 of file tenseur.C.

References Lorene::Tenseur::c, Lorene::Tenseur::del_derive(), and Lorene::Tenseur::n_comp.

◆ derive_con()

const Tenseur & Lorene::Tenseur::derive_con ( const Metrique &  metre) const
inherited

Returns the contravariant derivative of *this , with respect to met .

Definition at line 1584 of file tenseur.C.

References Lorene::Tenseur::fait_derive_con(), Lorene::Tenseur::get_place_met(), Lorene::Tenseur::p_derive_con, and Lorene::Tenseur::set_dependance().

◆ derive_cov()

const Tenseur & Lorene::Tenseur::derive_cov ( const Metrique &  met) const
inherited

◆ donne_indices()

Itbl Lorene::Tenseur_sym::donne_indices ( int  place) const
virtual

Returns the indices of a component given by its position in the Cmp 1-D array c .

Returns
1-D array of integers (Itbl ) of size valence giving the value of each index for the component located at the position place in the Cmp 1-D array c . Each element of this Itbl is 0, 1 or 2, which corresponds to spatial indices 1, 2 or 3 respectively.

Reimplemented from Lorene::Tenseur.

Definition at line 261 of file tenseur_sym.C.

References Lorene::Tenseur::n_comp, Lorene::Itbl::set(), Lorene::Itbl::set_etat_qcq(), and Lorene::Tenseur::valence.

◆ donne_place()

int Lorene::Tenseur_sym::donne_place ( const Itbl idx) const
virtual

Returns the position in the Cmp 1-D array c of a component given by its indices.

Returns
position in the Cmp 1-D array c
corresponding to the indices given in idx . idx must be a 1-D Itbl of size valence , each element of which must be 0, 1 or 2, corresponding to spatial indices 1, 2 or 3 respectively.

Reimplemented from Lorene::Tenseur.

Definition at line 216 of file tenseur_sym.C.

References Lorene::Itbl::get_dim(), Lorene::Itbl::get_ndim(), and Lorene::Tenseur::valence.

◆ fait_carre_scal()

void Lorene::Tenseur::fait_carre_scal ( const Metrique &  met,
int  i 
) const
protectedinherited

Calculates, if needed, the scalar square of *this , with respect to met .

The result is in *p_carre_scal[i]

Definition at line 1533 of file tenseur.C.

References Lorene::Tenseur::contract, Lorene::Tenseur::manipule, Lorene::Tenseur::p_carre_scal, Lorene::Tenseur::Tenseur(), and Lorene::Tenseur::valence.

◆ fait_derive_con()

void Lorene::Tenseur_sym::fait_derive_con ( const Metrique &  metre,
int  i 
) const
protectedvirtual

Calculates, if needed, the contravariant derivative of *this , with respect to met .

The result is in *p_derive_con[i]

Reimplemented from Lorene::Tenseur.

Definition at line 454 of file tenseur_sym.C.

References Lorene::Tenseur::contract, Lorene::Tenseur::derive_cov(), Lorene::Tenseur::gradient(), Lorene::Tenseur::p_derive_con, Tenseur_sym(), and Lorene::Tenseur::valence.

◆ fait_derive_cov()

void Lorene::Tenseur_sym::fait_derive_cov ( const Metrique &  met,
int  i 
) const
protectedvirtual

Calculates, if needed, the covariant derivative of *this , with respect to met .

The result is in *p_derive_cov[i]

Reimplemented from Lorene::Tenseur.

Definition at line 385 of file tenseur_sym.C.

References Lorene::Tenseur::etat.

◆ fait_gradient()

void Lorene::Tenseur_sym::fait_gradient ( ) const
protectedvirtual

Calculates, if needed, the gradient of *this .

The result is in *p_gradient

Reimplemented from Lorene::Tenseur.

Definition at line 338 of file tenseur_sym.C.

References Lorene::Tenseur::etat.

◆ fait_gradient_spher()

void Lorene::Tenseur::fait_gradient_spher ( ) const
protectedinherited

Calculates, if needed, the gradient of *this in a spherical orthonormal basis (scalar field only).

The result is in *p_gradient_spher

Definition at line 1408 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ get_etat()

int Lorene::Tenseur::get_etat ( ) const
inlineinherited

Returns the logical state.

Definition at line 710 of file tenseur.h.

References Lorene::Tenseur::etat.

◆ get_metric()

const Metrique* Lorene::Tenseur::get_metric ( ) const
inlineinherited

Returns a pointer on the metric defining the conformal factor for tensor densities.

Otherwise (case of a pure tensor), it returns 0x0.

Definition at line 748 of file tenseur.h.

References Lorene::Tenseur::metric.

◆ get_mp()

const Map* Lorene::Tenseur::get_mp ( ) const
inlineinherited

Returns pointer on the mapping.

Definition at line 702 of file tenseur.h.

References Lorene::Tenseur::mp.

◆ get_n_comp()

int Lorene::Tenseur::get_n_comp ( ) const
inlineinherited

Returns the number of components.

Definition at line 716 of file tenseur.h.

References Lorene::Tenseur::n_comp.

◆ get_place_met()

int Lorene::Tenseur::get_place_met ( const Metrique &  metre) const
protectedinherited

Returns the position of the pointer on metre in the array met_depend .

Definition at line 614 of file tenseur.C.

References Lorene::Tenseur::met_depend.

◆ get_poids()

double Lorene::Tenseur::get_poids ( ) const
inlineinherited

Returns the weight.

Definition at line 741 of file tenseur.h.

References Lorene::Tenseur::poids.

◆ get_triad()

const Base_vect* Lorene::Tenseur::get_triad ( ) const
inlineinherited

Returns the vectorial basis (triad) on which the components are defined.

Definition at line 707 of file tenseur.h.

References Lorene::Tenseur::triad.

◆ get_type_indice() [1/2]

int Lorene::Tenseur::get_type_indice ( int  i) const
inlineinherited

Returns 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
  • and so on...
Returns
COV for a covariant index, CON for a contravariant one.

Definition at line 729 of file tenseur.h.

References Lorene::Tenseur::type_indice.

◆ get_type_indice() [2/2]

Itbl Lorene::Tenseur::get_type_indice ( ) const
inlineinherited

Returns the types of all the indices.

Returns
1-D Itbl of size valence containing the type of each index, COV for a covariant one and CON
for a contravariant one.

Definition at line 738 of file tenseur.h.

References Lorene::Tenseur::type_indice.

◆ get_valence()

int Lorene::Tenseur::get_valence ( ) const
inlineinherited

Returns the valence.

Definition at line 713 of file tenseur.h.

References Lorene::Tenseur::valence.

◆ gradient()

const Tenseur & Lorene::Tenseur::gradient ( ) const
inherited

Returns the gradient of *this (Cartesian coordinates)

Definition at line 1558 of file tenseur.C.

References Lorene::Tenseur::fait_gradient(), and Lorene::Tenseur::p_gradient.

◆ gradient_spher()

const Tenseur & Lorene::Tenseur::gradient_spher ( ) const
inherited

Returns the gradient of *this (Spherical coordinates) (scalar field only).

Definition at line 1564 of file tenseur.C.

References Lorene::Tenseur::fait_gradient_spher(), and Lorene::Tenseur::p_gradient_spher.

◆ inc2_dzpuis()

void Lorene::Tenseur::inc2_dzpuis ( )
inherited

dzpuis += 2 ;

Definition at line 1159 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ inc_dzpuis()

void Lorene::Tenseur::inc_dzpuis ( )
inherited

dzpuis += 1 ;

Definition at line 1133 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ inverse_poisson_vect()

Tenseur Lorene::Tenseur::inverse_poisson_vect ( double  lambda) const
inherited

Compute $\Delta + \lambda \nabla\nabla$ of *this , *this being of valence 1.

Definition at line 62 of file tenseur_inv_pois_vect.C.

References Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ mult_r_zec()

void Lorene::Tenseur::mult_r_zec ( )
inherited

Multiplication by r in the external zone.

Definition at line 1172 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ new_der_met()

void Lorene::Tenseur::new_der_met ( )
protectedinherited

Builds the arrays met_depend , p_derive_cov , p_derive_con and p_carre_scal and fills them with null pointers.

Definition at line 212 of file tenseur.C.

References Lorene::Tenseur::met_depend, Lorene::Tenseur::p_carre_scal, Lorene::Tenseur::p_derive_con, Lorene::Tenseur::p_derive_cov, and Lorene::Tenseur::set_der_0x0().

◆ operator()() [1/5]

const Cmp & Lorene::Tenseur::operator() ( ) const
inherited

Read only for a scalar.

Definition at line 959 of file tenseur.C.

References Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ operator()() [2/5]

const Cmp & Lorene::Tenseur::operator() ( int  indice) const
inherited

Read only for a vector.

Definition at line 988 of file tenseur.C.

References Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ operator()() [3/5]

const Cmp & Lorene::Tenseur::operator() ( int  indice1,
int  indice2 
) const
inherited

Read only for a tensor of valence 2.

Definition at line 1016 of file tenseur.C.

References Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ operator()() [4/5]

const Cmp & Lorene::Tenseur::operator() ( int  indice1,
int  indice2,
int  indice3 
) const
inherited

Read only for a tensor of valence 3.

Definition at line 1050 of file tenseur.C.

References Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ operator()() [5/5]

const Cmp & Lorene::Tenseur::operator() ( const Itbl indices) const
inherited

Read only in the general case.

Definition at line 1087 of file tenseur.C.

References Lorene::Tenseur::etat, Lorene::Itbl::get_dim(), Lorene::Itbl::get_ndim(), and Lorene::Tenseur::valence.

◆ operator=()

void Lorene::Tenseur_sym::operator= ( const Tenseur t)
virtual

Assignment from a Tenseur .

The symmetry is assumed but not checked.

Reimplemented from Lorene::Tenseur.

Definition at line 292 of file tenseur_sym.C.

References Lorene::Tenseur::get_etat(), Lorene::Tenseur::get_valence(), Lorene::Tenseur::metric, Lorene::Tenseur::poids, Lorene::Tenseur::triad, Lorene::Tenseur::type_indice, and Lorene::Tenseur::valence.

◆ poisson_vect() [1/2]

void Lorene::Tenseur::poisson_vect ( double  lambda,
Param par,
Tenseur shift,
Tenseur vect,
Tenseur scal 
) const
inherited

Solves the vectorial Poisson equation : $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$.

with $\lambda \not= 1$.

*this must be given with dzpuis = 4.

It uses the Shibata scheme, where $N^i$ is given by :

\[ N^i = \frac{1}{2}\frac{\lambda+2}{\lambda+1}W^i-\frac{1}{2} \frac{\lambda}{\lambda+1}\left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

with $\Delta W^i = S^i$ and $\Delta \chi = -x_kS^k$.

Parameters
lambda[input] $\lambda$.
par[input/output] see Map::donne_para_poisson_vect.
shift[input] solution $N^i$ at the previous step. Zero if nothing is known.
shift[output] solution at this step.
vect[input/output] the same thing than for shift but for $W^i$.
scal[input/output] the same thing than for shift but for $\chi$.

Definition at line 121 of file tenseur_pde.C.

References Lorene::Tenseur::etat, Lorene::Tenseur::get_type_indice(), Lorene::Tenseur::get_valence(), Lorene::Tenseur::type_indice, and Lorene::Tenseur::valence.

◆ poisson_vect() [2/2]

Tenseur Lorene::Tenseur::poisson_vect ( double  lambda,
Tenseur vect,
Tenseur scal 
) const
inherited

Solves the vectorial Poisson equation $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$.

with $\lambda \not= 1$.

*this must be given with dzpuis = 4.

It uses the Shibata scheme, where $N^i$ is given by :

\[ N^i = \frac{1}{2}\frac{\lambda+2}{\lambda+1}W^i-\frac{1}{2} \frac{\lambda}{\lambda+1}\left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

with $\Delta W^i = S^i$ and $\Delta \chi = -x_kS^k$.

This version is to be used only with an affine mapping.

Parameters
lambda[input] $\lambda$.
vect[input] $W^i$ at the previous step. Zero if nothing is known.
vect[output] $W^i$ at this step.
scal[input/output] the same thing than for shift but for $\chi$.
Returns
the solution $N^i$.

Definition at line 203 of file tenseur_pde.C.

References Lorene::Tenseur::metric, Lorene::Tenseur::mp, Lorene::Tenseur::poids, Lorene::Tenseur::poisson_vect(), Lorene::Tenseur::set_etat_qcq(), Lorene::Tenseur::triad, Lorene::Tenseur::type_indice, and Lorene::Tenseur::valence.

◆ poisson_vect_oohara() [1/2]

void Lorene::Tenseur::poisson_vect_oohara ( double  lambda,
Param par,
Tenseur shift,
Tenseur scal 
) const
inherited

Solves the vectorial Poisson equation $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$.

with $\lambda \not= 1$.

*this must be given with dzpuis = 3 or 4 and be continuous.

It uses the Oohara scheme, where $N^i$ is given by

\[ \Delta N^i = S^i-\lambda \nabla^i \chi \]

with $\chi$ solution of :

\[ \Delta \chi = \frac{1}{\lambda+1}\nabla_k S^k \]

Parameters
lambda[input] $\lambda$.
par[input/output] see Map::donne_para_poisson_vect.
shift[input] solution $N^i$ at the previous step. Zero if nothing is known.
shift[output] solution at this step.
scal[input/output] the same thing than for shift but for $\chi$.

Definition at line 221 of file tenseur_pde.C.

References Lorene::Tenseur::etat, Lorene::Tenseur::get_type_indice(), Lorene::Tenseur::get_valence(), Lorene::Tenseur::type_indice, and Lorene::Tenseur::valence.

◆ poisson_vect_oohara() [2/2]

Tenseur Lorene::Tenseur::poisson_vect_oohara ( double  lambda,
Tenseur scal 
) const
inherited

Solves the vectorial Poisson equation $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$.

with $\lambda \not= 1$.

*this must be given with dzpuis = 3 or 4 and be continuous.

This version is to be used only with an affine mapping.

It uses the Oohara scheme, where $N^i$ is given by :

\[ \Delta N^i = S^i-\lambda \nabla^i \chi \]

with $\chi$ solution of :

\[ \Delta \chi = \frac{1}{\lambda+1}\nabla_k S^k \]

This version is to be used only with an affine mapping.

Parameters
lambda[input] $\lambda$.
scal[input] $\chi$ at the previous step. Zero if nothing is known.
scal[output] $\chi$ at this step.
Returns
the solution $N^i$.

Definition at line 293 of file tenseur_pde.C.

References Lorene::Tenseur::metric, Lorene::Tenseur::mp, Lorene::Tenseur::poids, Lorene::Tenseur::poisson_vect_oohara(), Lorene::Tenseur::set_etat_qcq(), Lorene::Tenseur::triad, Lorene::Tenseur::type_indice, and Lorene::Tenseur::valence.

◆ poisson_vect_regu()

void Lorene::Tenseur::poisson_vect_regu ( int  k_div,
int  nzet,
double  unsgam1,
double  lambda,
Param par,
Tenseur shift,
Tenseur vect,
Tenseur scal 
) const
inherited

Solves the vectorial Poisson equation : $\Delta N^i +\lambda \nabla^i \nabla_k N^k = S^i$.

with $\lambda \not= 1$ by regularizing the source term.

*this must be given with dzpuis = 4.

It uses the Shibata scheme, where $N^i$ is given by :

\[ N^i = \frac{1}{2}\frac{\lambda+2}{\lambda+1}W^i-\frac{1}{2} \frac{\lambda}{\lambda+1}\left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

with $\Delta W^i = S^i$ and $\Delta \chi = -x_kS^k$.

Parameters
k_div[input] regularization degree.
nzet[input] number of domains covering a star.
unsgam1[input] $1/(\gamma - 1)$.
lambda[input] $\lambda$.
par[input/output] see Map::donne_para_poisson_vect.
shift[input] solution $N^i$ at the previous step. Zero if nothing is known.
shift[output] solution at this step.
vect[input/output] the same thing than for shift but for $W^i$.
scal[input/output] the same thing than for shift but for $\chi$.

Definition at line 74 of file tenseur_pde_regu.C.

References Lorene::Tenseur::etat, Lorene::Tenseur::get_type_indice(), Lorene::Tenseur::get_valence(), Lorene::Tenseur::type_indice, and Lorene::Tenseur::valence.

◆ sauve()

void Lorene::Tenseur::sauve ( FILE *  fd) const
inherited

◆ set() [1/5]

Cmp & Lorene::Tenseur::set ( )
inherited

Read/write for a scalar (see also operator=(const Cmp&) ).

Definition at line 840 of file tenseur.C.

References Lorene::Tenseur::del_derive(), and Lorene::Tenseur::etat.

◆ set() [2/5]

Cmp & Lorene::Tenseur::set ( int  ind)
inherited

Read/write for a vector.

Definition at line 850 of file tenseur.C.

References Lorene::Tenseur::del_derive(), Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ set() [3/5]

Cmp & Lorene::Tenseur::set ( int  ind1,
int  ind2 
)
inherited

Read/write for a tensor of valence 2.

Definition at line 861 of file tenseur.C.

References Lorene::Tenseur::del_derive(), Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ set() [4/5]

Cmp & Lorene::Tenseur::set ( int  ind1,
int  ind2,
int  ind3 
)
inherited

Read/write for a tensor of valence 3.

Definition at line 880 of file tenseur.C.

References Lorene::Tenseur::del_derive(), Lorene::Tenseur::etat, and Lorene::Tenseur::valence.

◆ set() [5/5]

Cmp & Lorene::Tenseur::set ( const Itbl indices)
inherited

◆ set_dependance()

void Lorene::Tenseur::set_dependance ( const Metrique &  met) const
protectedinherited

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 624 of file tenseur.C.

References Lorene::Tenseur::met_depend.

◆ set_der_0x0()

void Lorene::Tenseur::set_der_0x0 ( ) const
protectedinherited

Sets the pointers of all the derivatives to zero.

Definition at line 607 of file tenseur.C.

References Lorene::Tenseur::p_gradient, Lorene::Tenseur::p_gradient_spher, and Lorene::Tenseur::set_der_met_0x0().

◆ set_der_met_0x0()

void Lorene::Tenseur::set_der_met_0x0 ( int  i) const
protectedinherited

Sets the pointers of the derivatives depending on the i-th element of *met_depend to zero (as well as that i-th element).

Definition at line 599 of file tenseur.C.

References Lorene::Tenseur::met_depend, Lorene::Tenseur::p_carre_scal, Lorene::Tenseur::p_derive_con, and Lorene::Tenseur::p_derive_cov.

◆ set_etat_nondef()

void Lorene::Tenseur::set_etat_nondef ( )
inherited

Sets the logical state to ETATNONDEF (undefined state).

The components are not allocated.

Definition at line 666 of file tenseur.C.

References Lorene::Tenseur::del_t(), and Lorene::Tenseur::etat.

◆ set_etat_qcq()

void Lorene::Tenseur::set_etat_qcq ( )
inherited

Sets the logical state to ETATQCQ (ordinary state).

The components are now allocated and set to ETATNONDEF .

Definition at line 652 of file tenseur.C.

References Lorene::Tenseur::c, Lorene::Tenseur::del_derive(), Lorene::Tenseur::etat, Lorene::Tenseur::mp, and Lorene::Tenseur::n_comp.

◆ set_etat_zero()

void Lorene::Tenseur::set_etat_zero ( )
inherited

Sets the logical state to ETATZERO (zero state).

The components are not allocated.

Definition at line 661 of file tenseur.C.

References Lorene::Tenseur::del_t(), and Lorene::Tenseur::etat.

◆ set_metric()

void Lorene::Tenseur::set_metric ( const Metrique &  met)
inherited

Sets the pointer on the metric for a tensor density.

Definition at line 701 of file tenseur.C.

References Lorene::Tenseur::metric.

◆ set_poids()

void Lorene::Tenseur::set_poids ( double  weight)
inherited

Sets the weight for a tensor density.

Definition at line 696 of file tenseur.C.

References Lorene::Tenseur::poids.

◆ set_std_base()

void Lorene::Tenseur::set_std_base ( )
inherited

Set the standard spectal basis of decomposition for each component.

To be used only with valence strictly lower than 3.

Definition at line 1186 of file tenseur.C.

References Lorene::Tenseur::etat.

◆ set_triad()

void Lorene::Tenseur::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 c ). In order to change them coherently with the new basis, the function change_triad(const Base_vect&) must be called instead.

Definition at line 690 of file tenseur.C.

References Lorene::Tenseur::triad.

◆ verif()

bool Lorene::Tenseur::verif ( ) const
protectedinherited

Returns false for a tensor density without a defined metric.

Definition at line 208 of file tenseur.C.

References Lorene::Tenseur::metric, and Lorene::Tenseur::poids.

Friends And Related Function Documentation

◆ lie_derive

Tenseur lie_derive ( const Tenseur t,
const Tenseur x,
const Metrique *  met = 0x0 
)
friend

Lie Derivative of t with respect to x .

If no other argument is given, it uses partial derivatives with respect to cartesian coordinates to calculate the result (this is the default). Otherwise, it uses the covariant derivative associated to the metric given as last argument.

Definition at line 819 of file tenseur_operateur.C.

◆ manipule

Tenseur_sym manipule ( const Tenseur_sym t1,
const Metrique &  met 
)
friend

Raise or lower all the indices, depending on their type, using the given Metrique .

Definition at line 140 of file tenseur_sym_operateur.C.

◆ operator*

Tenseur_sym operator* ( const Tenseur t1,
const Tenseur_sym t2 
)
friend

Tensorial product.

Definition at line 87 of file tenseur_sym_operateur.C.

Member Data Documentation

◆ c

Cmp** Lorene::Tenseur::c
protectedinherited

The components.

Definition at line 325 of file tenseur.h.

◆ etat

int Lorene::Tenseur::etat
protectedinherited

Logical state ETATZERO , ETATQCQ or ETATNONDEF.

Definition at line 324 of file tenseur.h.

◆ met_depend

const Metrique** Lorene::Tenseur::met_depend
protectedinherited

Array of pointers on the Metrique 's used to calculate derivatives members.

If no such member has been calculated the pointer is set to zero. The array has the size N_MET_MAX abd the i-th corresponds to the i-th ones in p_derive_cov and p_derive_con.

Definition at line 340 of file tenseur.h.

◆ metric

const Metrique* Lorene::Tenseur::metric
protectedinherited

For tensor densities: the metric defining the conformal factor.

Definition at line 328 of file tenseur.h.

◆ mp

const Map* const Lorene::Tenseur::mp
protectedinherited

Reference mapping.

Definition at line 309 of file tenseur.h.

◆ n_comp

int Lorene::Tenseur::n_comp
protectedinherited

Number of components, depending on the symmetry.

Definition at line 323 of file tenseur.h.

◆ p_carre_scal

Tenseur** Lorene::Tenseur::p_carre_scal
protectedinherited

Array of pointers on the scalar squares of *this
with respect to the corresponding metric in *met_depend .

It is set to zero if it has not been calculated yet.

Definition at line 375 of file tenseur.h.

◆ p_derive_con

Tenseur** Lorene::Tenseur::p_derive_con
protectedinherited

Array of pointers on the contravariant derivatives of *this
with respect to the corresponding metric in *met_depend .

It is set to zero if it has not been calculated yet.

Definition at line 368 of file tenseur.h.

◆ p_derive_cov

Tenseur** Lorene::Tenseur::p_derive_cov
protectedinherited

Array of pointers on the covariant derivatives of *this
with respect to the corresponding metric in *met_depend .

It is set to zero if it has not been calculated yet, or for a scalar field.

Definition at line 361 of file tenseur.h.

◆ p_gradient

Tenseur* Lorene::Tenseur::p_gradient
mutableprotectedinherited

Pointer on the gradient of *this .

It is set to zero if it has not been calculated yet.

Definition at line 346 of file tenseur.h.

◆ p_gradient_spher

Tenseur* Lorene::Tenseur::p_gradient_spher
mutableprotectedinherited

Pointer on the gradient of *this in a spherical orthonormal basis (scalar field only).

It is set to zero if it has not been calculated yet.

Definition at line 353 of file tenseur.h.

◆ poids

double Lorene::Tenseur::poids
protectedinherited

For tensor densities: the weight.

Definition at line 326 of file tenseur.h.

◆ triad

const Base_vect* Lorene::Tenseur::triad
protectedinherited

Vectorial basis (triad) with respect to which the tensor components are defined.

Definition at line 315 of file tenseur.h.

◆ type_indice

Itbl Lorene::Tenseur::type_indice
protectedinherited

Array of size valence contening the type of each index, COV for a covariant one and CON for a contravariant one.

Definition at line 321 of file tenseur.h.

◆ valence

int Lorene::Tenseur::valence
protectedinherited

Valence.

Definition at line 310 of file tenseur.h.


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