Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class Tensor_sym
instead ***.
More...
#include <tenseur.h>
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. | |
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. | |
Tenseur_sym (const Tenseur_sym &) | |
Copy constructor. | |
Tenseur_sym (const Tenseur &) | |
Constructor from a Tenseur . | |
Tenseur_sym (const Map &map, const Base_vect &triad_i, FILE *fich, const Metrique *met=0x0) | |
Constructor from a file (see sauve(FILE*) ). | |
virtual | ~Tenseur_sym () |
Destructor. | |
virtual void | operator= (const Tenseur &) |
Assignment from a Tenseur . | |
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. | |
virtual Itbl | donne_indices (int place) const |
Returns the indices of a component given by its position in the Cmp 1-D array c . | |
void | set_etat_nondef () |
Sets the logical state to ETATNONDEF (undefined state). | |
void | set_etat_zero () |
Sets the logical state to ETATZERO (zero state). | |
void | set_etat_qcq () |
Sets the logical state to ETATQCQ (ordinary state). | |
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 | change_triad (const Base_vect &new_triad) |
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly. | |
void | set_triad (const Base_vect &new_triad) |
Assigns a new vectorial basis (triad) of decomposition. | |
void | set_poids (double weight) |
Sets the weight for a tensor density. | |
void | set_metric (const Metrique &met) |
Sets the pointer on the metric for a tensor density. | |
Cmp & | set () |
Read/write for a scalar (see also operator= (const Cmp& ) ). | |
Cmp & | set (int) |
Read/write for a vector. | |
Cmp & | set (int, int) |
Read/write for a tensor of valence 2. | |
Cmp & | set (int, int, int) |
Read/write for a tensor of valence 3. | |
Cmp & | set (const Itbl &) |
Read/write in the general case. | |
void | annule (int l) |
Sets the Tenseur to zero in a given domain. | |
void | annule (int l_min, int l_max) |
Sets the Tenseur to zero in several domains. | |
void | set_std_base () |
Set the standard spectal basis of decomposition for each component. | |
void | dec_dzpuis () |
dzpuis -= 1 ; | |
void | inc_dzpuis () |
dzpuis += 1 ; | |
void | dec2_dzpuis () |
dzpuis -= 2 ; | |
void | inc2_dzpuis () |
dzpuis += 2 ; | |
void | mult_r_zec () |
Multiplication by r in the external zone. | |
Tenseur | inverse_poisson_vect (double lambda) const |
Compute of *this , *this being of valence 1. | |
const Map * | get_mp () const |
Returns pointer on the mapping. | |
const Base_vect * | get_triad () const |
Returns the vectorial basis (triad) on which the components are defined. | |
int | get_etat () const |
Returns the logical state. | |
int | get_valence () const |
Returns the valence. | |
int | get_n_comp () const |
Returns the number of components. | |
int | get_type_indice (int i) const |
Returns the type of the index number i . | |
Itbl | get_type_indice () const |
Returns the types of all the indices. | |
double | get_poids () const |
Returns the weight. | |
const Metrique * | get_metric () const |
Returns a pointer on the metric defining the conformal factor for tensor densities. | |
const Cmp & | operator() () const |
Read only for a scalar. | |
const Cmp & | operator() (int) const |
Read only for a vector. | |
const Cmp & | operator() (int, int) const |
Read only for a tensor of valence 2. | |
const Cmp & | operator() (int, int, int) const |
Read only for a tensor of valence 3. | |
const Cmp & | operator() (const Itbl &) const |
Read only in the general case. | |
void | sauve (FILE *) const |
Save in a file. | |
const Tenseur & | gradient () const |
Returns the gradient of *this (Cartesian coordinates). | |
const Tenseur & | gradient_spher () const |
Returns the gradient of *this (Spherical coordinates) (scalar field only). | |
const Tenseur & | derive_cov (const Metrique &met) const |
Returns the covariant derivative of *this , with respect to met . | |
const Tenseur & | derive_con (const Metrique &) const |
Returns the contravariant derivative of *this , with respect to met . | |
const Tenseur & | carre_scal (const Metrique &) const |
Returns the scalar square of *this , with respect to met . | |
void | poisson_vect (double lambda, Param &par, Tenseur &shift, Tenseur &vect, Tenseur &scal) const |
Solves the vectorial Poisson equation : . | |
Tenseur | poisson_vect (double lambda, Tenseur &vect, Tenseur &scal) const |
Solves the vectorial Poisson equation . | |
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 ¶, 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 . | |
Tenseur | poisson_vect_oohara (double lambda, Tenseur &scal) const |
Solves the vectorial Poisson equation . | |
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 : . | |
void | compare (const Tenseur &tens, const char *name) |
Functions to compare the values of tensors. | |
void | compare (FILE *fich, const char *name_i) |
Protected Member Functions | |
virtual void | fait_gradient () const |
Calculates, if needed, the gradient of *this . | |
virtual void | fait_derive_cov (const Metrique &met, int i) const |
Calculates, if needed, the covariant derivative of *this , with respect to met . | |
virtual void | fait_derive_con (const Metrique &, int i) const |
Calculates, if needed, the contravariant derivative of *this , with respect to met . | |
bool | verif () const |
Returns false for a tensor density without a defined metric. | |
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. | |
void | del_t () |
Logical destructor. | |
void | del_derive_met (int i) const |
Logical destructor of the derivatives depending on the i-th element of *met_depend . | |
void | del_derive () const |
Logical destructor of all the derivatives. | |
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). | |
void | set_der_0x0 () const |
Sets the pointers of all the derivatives to zero. | |
void | fait_gradient_spher () const |
Calculates, if needed, the gradient of *this in a spherical orthonormal basis (scalar field only). | |
void | fait_carre_scal (const Metrique &, int i) const |
Calculates, if needed, the scalar square of *this , with respect to met . | |
void | set_dependance (const Metrique &met) const |
To be used to describe the fact that the derivatives members have been calculated with met . | |
int | get_place_met (const Metrique &metre) const |
Returns the position of the pointer on metre in the array met_depend . | |
Protected Attributes | |
const Map *const | mp |
Reference mapping. | |
int | valence |
Valence. | |
const Base_vect * | triad |
Vectorial basis (triad) with respect to which the tensor components are defined. | |
Itbl | type_indice |
Array of size valence contening the type of each index, COV for a covariant one and CON for a contravariant one. | |
int | n_comp |
Number of components, depending on the symmetry. | |
int | etat |
Logical state ETATZERO , ETATQCQ or ETATNONDEF . | |
Cmp ** | c |
The components. | |
double | poids |
For tensor densities: the weight. | |
const Metrique * | metric |
For tensor densities: the metric defining the conformal factor. | |
const Metrique ** | met_depend |
Array of pointers on the Metrique 's used to calculate derivatives members. | |
Tenseur * | p_gradient |
Pointer on the gradient of *this . | |
Tenseur * | p_gradient_spher |
Pointer on the gradient of *this in a spherical orthonormal basis (scalar field only). | |
Tenseur ** | p_derive_cov |
Array of pointers on the covariant derivatives of *this with respect to the corresponding metric in *met_depend . | |
Tenseur ** | p_derive_con |
Array of pointers on the contravariant derivatives of *this with respect to the corresponding metric in *met_depend . | |
Tenseur ** | p_carre_scal |
Array of pointers on the scalar squares of *this with respect to the corresponding metric in *met_depend . | |
Friends | |
Tenseur_sym | operator* (const Tenseur &, const Tenseur_sym &) |
Tenseur_sym | manipule (const Tenseur_sym &, const Metrique &) |
Tenseur | lie_derive (const Tenseur &, const Tenseur &, const Metrique *) |
Lie Derivative of t with respect to x . | |
class | Tenseur_sym |
class | Metrique |
ostream & | operator<< (ostream &, const Tenseur &) |
Tenseur | operator* (const Tenseur &, const Tenseur &) |
Tenseur | operator% (const Tenseur &, const Tenseur &) |
Tenseur | contract (const Tenseur &, int id1, int id2) |
Tenseur | contract (const Tenseur &, int id1, const Tenseur &, int id2) |
Tenseur | contract_desal (const Tenseur &, int id1, const Tenseur &, int id2) |
Tenseur | flat_scalar_prod (const Tenseur &t1, const Tenseur &t2) |
Tenseur | flat_scalar_prod_desal (const Tenseur &t1, const Tenseur &t2) |
Tenseur | manipule (const Tenseur &, const Metrique &, int idx) |
Tenseur | skxk (const Tenseur &) |
Tenseur | lie_derive (const Tenseur &, const Tenseur &, const Metrique *) |
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 1253 of file tenseur.h.
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.
map | the mapping | |
val | valence of the tensor; must be greater or equal to 2. | |
tipe | 1-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:
| |
triad_i | vectorial basis (triad) with respect to which the tensor components are defined |
Definition at line 92 of file tenseur_sym.C.
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.
map | the mapping | |
val | valence of the tensor; must be greater or equal to 2. | |
tipe | the type of the indices. | |
triad_i | vectorial basis (triad) with respect to which the tensor components are defined |
Definition at line 103 of file tenseur_sym.C.
Tenseur_sym::Tenseur_sym | ( | const Tenseur_sym & | source | ) |
Copy constructor.
Definition at line 114 of file tenseur_sym.C.
References Tenseur::c, donne_indices(), donne_place(), Tenseur::etat, Tenseur::met_depend, Tenseur::n_comp, Tenseur::p_carre_scal, Tenseur::p_derive_con, Tenseur::p_derive_cov, Tenseur::p_gradient, Tenseur::set_dependance(), Tenseur::Tenseur(), Tenseur_sym(), and Tenseur::valence.
Tenseur_sym::Tenseur_sym | ( | const Tenseur & | source | ) | [explicit] |
Constructor from a Tenseur
.
The symmetry is assumed to be true but not checked.
Definition at line 155 of file tenseur_sym.C.
References Tenseur::c, donne_indices(), Tenseur::donne_place(), Tenseur::etat, Tenseur::met_depend, Tenseur::n_comp, Tenseur::p_carre_scal, Tenseur::p_derive_con, Tenseur::p_derive_cov, Tenseur::p_gradient, Tenseur::Tenseur(), and Tenseur::valence.
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*)
).
map | the mapping | |
triad_i | vectorial 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. | |
fich | file which has been created by the function sauve(FILE*) . |
Definition at line 185 of file tenseur_sym.C.
References Tenseur::n_comp, pow(), and Tenseur::valence.
Tenseur_sym::~Tenseur_sym | ( | ) | [virtual] |
Destructor.
Definition at line 197 of file tenseur_sym.C.
void 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 652 of file tenseur.C.
References Cmp::allocate_all(), Tenseur::c, Tenseur::n_comp, and Tenseur::set_etat_qcq().
void Tenseur::annule | ( | int | l_min, | |
int | l_max | |||
) | [inherited] |
Sets the Tenseur
to zero in several domains.
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 annule
(0,nz-1), where nz
is the total number of domains, is equivalent to set_etat_zero()
.
Definition at line 900 of file tenseur.C.
References Tenseur::annule(), Cmp::annule(), Tenseur::c, Tenseur::etat, Map::get_mg(), Mg3d::get_nzone(), Tenseur::mp, Tenseur::n_comp, Tenseur::p_carre_scal, Tenseur::p_derive_con, Tenseur::p_derive_cov, Tenseur::p_gradient, Tenseur::p_gradient_spher, and Tenseur::set_etat_zero().
void Tenseur::annule | ( | int | l | ) | [inherited] |
const Tenseur & Tenseur::carre_scal | ( | const Metrique & | metre | ) | const [inherited] |
Returns the scalar square of *this
, with respect to met
.
Definition at line 1572 of file tenseur.C.
References Tenseur::fait_carre_scal(), Tenseur::get_place_met(), Tenseur::p_carre_scal, and Tenseur::set_dependance().
void 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 663 of file tenseur.C.
References Base_vect::change_basis().
void Tenseur::compare | ( | const Tenseur & | tens, | |
const char * | name | |||
) | [inherited] |
Functions to compare the values of tensors.
Definition at line 53 of file tenseur_compare.C.
References Cmp::compare(), Tenseur::get_valence(), and Tenseur::valence.
void Tenseur::dec2_dzpuis | ( | ) | [inherited] |
dzpuis -= 2 ;
Definition at line 1125 of file tenseur.C.
References Tenseur::c, Cmp::dec2_dzpuis(), Tenseur::etat, and Tenseur::n_comp.
void Tenseur::dec_dzpuis | ( | ) | [inherited] |
dzpuis -= 1 ;
Definition at line 1099 of file tenseur.C.
References Tenseur::c, Cmp::dec_dzpuis(), Tenseur::etat, and Tenseur::n_comp.
void Tenseur::del_derive | ( | ) | const [protected, inherited] |
Logical destructor of all the derivatives.
Definition at line 568 of file tenseur.C.
References Tenseur::del_derive_met(), Tenseur::p_gradient, Tenseur::p_gradient_spher, and Tenseur::set_der_0x0().
void Tenseur::del_derive_met | ( | int | i | ) | const [protected, inherited] |
Logical destructor of the derivatives depending on the i-th element of *met_depend
.
Definition at line 549 of file tenseur.C.
References Tenseur::met_depend, Tenseur::p_carre_scal, Tenseur::p_derive_con, Tenseur::p_derive_cov, and Tenseur::set_der_met_0x0().
void Tenseur::del_t | ( | ) | [protected, inherited] |
Logical destructor.
Definition at line 540 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), and Tenseur::n_comp.
const Tenseur & Tenseur::derive_con | ( | const Metrique & | metre | ) | const [inherited] |
Returns the contravariant derivative of *this
, with respect to met
.
Definition at line 1563 of file tenseur.C.
References Tenseur::fait_derive_con(), Tenseur::get_place_met(), Tenseur::p_derive_con, and Tenseur::set_dependance().
const Tenseur & Tenseur::derive_cov | ( | const Metrique & | met | ) | const [inherited] |
Returns the covariant derivative of *this
, with respect to met
.
Definition at line 1549 of file tenseur.C.
References Tenseur::fait_derive_cov(), Tenseur::get_place_met(), Tenseur::gradient(), Tenseur::p_derive_cov, Tenseur::set_dependance(), and Tenseur::valence.
Itbl 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
.
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 Tenseur.
Definition at line 248 of file tenseur_sym.C.
References Tenseur::n_comp, Itbl::set(), Itbl::set_etat_qcq(), and Tenseur::valence.
int 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.
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 Tenseur.
Definition at line 203 of file tenseur_sym.C.
References Itbl::get_dim(), Itbl::get_ndim(), and Tenseur::valence.
void Tenseur::fait_carre_scal | ( | const Metrique & | met, | |
int | i | |||
) | const [protected, inherited] |
Calculates, if needed, the scalar square of *this
, with respect to met
.
The result is in *p_carre_scal
[i]
Definition at line 1512 of file tenseur.C.
References Tenseur::p_carre_scal, Tenseur::Tenseur(), and Tenseur::valence.
void Tenseur_sym::fait_derive_con | ( | const Metrique & | metre, | |
int | i | |||
) | const [protected, virtual] |
Calculates, if needed, the contravariant derivative of *this
, with respect to met
.
The result is in *p_derive_con
[i]
Reimplemented from Tenseur.
Definition at line 437 of file tenseur_sym.C.
References Tenseur::derive_cov(), Tenseur::gradient(), Tenseur::p_derive_con, Tenseur_sym(), and Tenseur::valence.
void Tenseur_sym::fait_derive_cov | ( | const Metrique & | met, | |
int | i | |||
) | const [protected, virtual] |
Calculates, if needed, the covariant derivative of *this
, with respect to met
.
The result is in *p_derive_cov
[i]
Reimplemented from Tenseur.
Definition at line 368 of file tenseur_sym.C.
References donne_indices(), Tenseur::etat, Tenseur::gradient(), Tenseur::n_comp, Tenseur::p_derive_cov, Tenseur::poids, Tenseur::set(), Itbl::set(), Itbl::set_etat_qcq(), Tenseur::Tenseur(), Tenseur_sym(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
void Tenseur_sym::fait_gradient | ( | ) | const [protected, virtual] |
Calculates, if needed, the gradient of *this
.
The result is in *p_gradient
Reimplemented from Tenseur.
Definition at line 321 of file tenseur_sym.C.
References Tenseur::etat, Map::get_bvect_cart(), Tenseur::metric, Tenseur::mp, Tenseur::p_gradient, Tenseur::poids, Itbl::set(), Itbl::set_etat_qcq(), Tenseur_sym(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
void Tenseur::fait_gradient_spher | ( | ) | const [protected, inherited] |
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 1387 of file tenseur.C.
References Tenseur::c, Cmp::dsdr(), Tenseur::etat, Map::get_bvect_spher(), Tenseur::metric, Tenseur::mp, Tenseur::p_gradient_spher, Tenseur::poids, Tenseur::set(), Tenseur::set_etat_qcq(), Tenseur::set_etat_zero(), Cmp::srdsdt(), Cmp::srstdsdp(), Tenseur::Tenseur(), and Tenseur::valence.
int Tenseur::get_etat | ( | ) | const [inline, inherited] |
const Metrique* Tenseur::get_metric | ( | ) | const [inline, inherited] |
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 738 of file tenseur.h.
References Tenseur::metric.
const Map* Tenseur::get_mp | ( | ) | const [inline, inherited] |
int Tenseur::get_n_comp | ( | ) | const [inline, inherited] |
Returns the number of components.
Definition at line 706 of file tenseur.h.
References Tenseur::n_comp.
int Tenseur::get_place_met | ( | const Metrique & | metre | ) | const [protected, inherited] |
Returns the position of the pointer on metre
in the array met_depend
.
Definition at line 593 of file tenseur.C.
References Tenseur::met_depend.
double Tenseur::get_poids | ( | ) | const [inline, inherited] |
const Base_vect* Tenseur::get_triad | ( | ) | const [inline, inherited] |
Returns the vectorial basis (triad) on which the components are defined.
Definition at line 697 of file tenseur.h.
References Tenseur::triad.
Itbl Tenseur::get_type_indice | ( | ) | 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 728 of file tenseur.h.
References Tenseur::type_indice.
int Tenseur::get_type_indice | ( | int | i | ) | const [inline, inherited] |
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 Definition at line 719 of file tenseur.h.
References Tenseur::type_indice.
int Tenseur::get_valence | ( | ) | const [inline, inherited] |
const Tenseur & Tenseur::gradient | ( | ) | const [inherited] |
Returns the gradient of *this
(Cartesian coordinates).
Definition at line 1537 of file tenseur.C.
References Tenseur::fait_gradient(), and Tenseur::p_gradient.
const Tenseur & Tenseur::gradient_spher | ( | ) | const [inherited] |
Returns the gradient of *this
(Spherical coordinates) (scalar field only).
Definition at line 1543 of file tenseur.C.
References Tenseur::fait_gradient_spher(), and Tenseur::p_gradient_spher.
void Tenseur::inc2_dzpuis | ( | ) | [inherited] |
dzpuis += 2 ;
Definition at line 1138 of file tenseur.C.
References Tenseur::c, Tenseur::etat, Cmp::inc2_dzpuis(), and Tenseur::n_comp.
void Tenseur::inc_dzpuis | ( | ) | [inherited] |
dzpuis += 1 ;
Definition at line 1112 of file tenseur.C.
References Tenseur::c, Tenseur::etat, Cmp::inc_dzpuis(), and Tenseur::n_comp.
Tenseur Tenseur::inverse_poisson_vect | ( | double | lambda | ) | const [inherited] |
Compute of *this
, *this
being of valence 1.
Definition at line 55 of file tenseur_inv_pois_vect.C.
References Tenseur::dec2_dzpuis(), Tenseur::etat, Tenseur::get_triad(), Tenseur::gradient(), Tenseur::metric, Tenseur::mp, Tenseur::poids, Tenseur::set(), Tenseur::set_etat_qcq(), and Tenseur::valence.
void Tenseur::mult_r_zec | ( | ) | [inherited] |
Multiplication by r in the external zone.
Definition at line 1151 of file tenseur.C.
References Tenseur::c, Tenseur::etat, Cmp::mult_r_zec(), and Tenseur::n_comp.
void Tenseur::new_der_met | ( | ) | [protected, inherited] |
Builds the arrays met_depend
, p_derive_cov
, p_derive_con
and p_carre_scal
and fills them with null pointers.
Definition at line 191 of file tenseur.C.
References Tenseur::met_depend, Tenseur::p_carre_scal, Tenseur::p_derive_con, Tenseur::p_derive_cov, and Tenseur::set_der_0x0().
Read only in the general case.
Definition at line 1066 of file tenseur.C.
References Tenseur::c, Map::cmp_zero(), Tenseur::donne_place(), Tenseur::etat, Itbl::get_dim(), Itbl::get_ndim(), Tenseur::mp, and Tenseur::valence.
const Cmp & Tenseur::operator() | ( | int | indice1, | |
int | indice2, | |||
int | indice3 | |||
) | const [inherited] |
Read only for a tensor of valence 3.
Definition at line 1029 of file tenseur.C.
References Tenseur::c, Map::cmp_zero(), Tenseur::donne_place(), Tenseur::etat, Tenseur::mp, Itbl::set(), Itbl::set_etat_qcq(), and Tenseur::valence.
const Cmp & Tenseur::operator() | ( | int | indice1, | |
int | indice2 | |||
) | const [inherited] |
Read only for a tensor of valence 2.
Definition at line 995 of file tenseur.C.
References Tenseur::c, Map::cmp_zero(), Tenseur::donne_place(), Tenseur::etat, Tenseur::mp, Itbl::set(), Itbl::set_etat_qcq(), and Tenseur::valence.
const Cmp & Tenseur::operator() | ( | int | indice | ) | const [inherited] |
Read only for a vector.
Definition at line 967 of file tenseur.C.
References Tenseur::c, Map::cmp_zero(), Tenseur::etat, Tenseur::mp, and Tenseur::valence.
const Cmp & Tenseur::operator() | ( | ) | const [inherited] |
Read only for a scalar.
Definition at line 938 of file tenseur.C.
References Tenseur::c, Map::cmp_zero(), Tenseur::etat, Tenseur::mp, and Tenseur::valence.
void Tenseur_sym::operator= | ( | const Tenseur & | t | ) | [virtual] |
Assignment from a Tenseur
.
The symmetry is assumed but not checked.
Definition at line 279 of file tenseur_sym.C.
References Tenseur::c, donne_indices(), Tenseur::donne_place(), Tenseur::get_etat(), Tenseur::get_valence(), Tenseur::metric, Tenseur::n_comp, Tenseur::poids, Tenseur::set_etat_nondef(), Tenseur::set_etat_qcq(), Tenseur::set_etat_zero(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
Solves the vectorial Poisson equation .
with .
*this
must be given with dzpuis
= 4.
It uses the Shibata scheme, where is given by :
with and .
This version is to be used only with an affine mapping.
lambda | [input] . | |
vect | [input] at the previous step. Zero if nothing is known. | |
vect | [output] at this step. | |
scal | [input/output] the same thing than for shift but for . |
Definition at line 196 of file tenseur_pde.C.
References Tenseur::metric, Tenseur::mp, Tenseur::poids, Tenseur::poisson_vect(), Tenseur::set_etat_qcq(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
void Tenseur::poisson_vect | ( | double | lambda, | |
Param & | par, | |||
Tenseur & | shift, | |||
Tenseur & | vect, | |||
Tenseur & | scal | |||
) | const [inherited] |
Solves the vectorial Poisson equation : .
with .
*this
must be given with dzpuis
= 4.
It uses the Shibata scheme, where is given by :
with and .
lambda | [input] . | |
par | [input/output] see Map::donne_para_poisson_vect. | |
shift | [input] solution 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 . | |
scal | [input/output] the same thing than for shift but for . |
Definition at line 114 of file tenseur_pde.C.
References Tenseur::dec_dzpuis(), Map::donne_para_poisson_vect(), Tenseur::etat, Tenseur::get_type_indice(), Tenseur::get_valence(), Tenseur::gradient(), Tenseur::mp, Tenseur::set(), Tenseur::set_etat_qcq(), Tenseur::set_etat_zero(), Tenseur::set_triad(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
Solves the vectorial Poisson equation .
with .
*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 is given by :
with solution of :
This version is to be used only with an affine mapping.
lambda | [input] . | |
scal | [input] at the previous step. Zero if nothing is known. | |
scal | [output] at this step. |
Definition at line 286 of file tenseur_pde.C.
References Tenseur::metric, Tenseur::mp, Tenseur::poids, Tenseur::poisson_vect_oohara(), Tenseur::set_etat_qcq(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
void Tenseur::poisson_vect_oohara | ( | double | lambda, | |
Param & | par, | |||
Tenseur & | shift, | |||
Tenseur & | scal | |||
) | const [inherited] |
Solves the vectorial Poisson equation .
with .
*this
must be given with dzpuis
= 3 or 4 and be continuous.
It uses the Oohara scheme, where is given by
with solution of :
lambda | [input] . | |
par | [input/output] see Map::donne_para_poisson_vect. | |
shift | [input] solution 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 . |
Definition at line 214 of file tenseur_pde.C.
References Tenseur::dec2_dzpuis(), Tenseur::dec_dzpuis(), Map::donne_para_poisson_vect(), Tenseur::etat, Tenseur::get_etat(), Tenseur::get_triad(), Tenseur::get_type_indice(), Tenseur::get_valence(), Tenseur::gradient(), Tenseur::mp, Tenseur::set(), Tenseur::set_etat_qcq(), Tenseur::set_etat_zero(), Tenseur::set_triad(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
void 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 : .
with by regularizing the source term.
*this
must be given with dzpuis
= 4.
It uses the Shibata scheme, where is given by :
with and .
k_div | [input] regularization degree. | |
nzet | [input] number of domains covering a star. | |
unsgam1 | [input] . | |
lambda | [input] . | |
par | [input/output] see Map::donne_para_poisson_vect. | |
shift | [input] solution 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 . | |
scal | [input/output] the same thing than for shift but for . |
Definition at line 67 of file tenseur_pde_regu.C.
References Tenseur::dec_dzpuis(), Map::donne_para_poisson_vect(), Tenseur::etat, Map::get_bvect_cart(), Tenseur::get_type_indice(), Tenseur::get_valence(), Tenseur::gradient(), Tenseur::metric, Tenseur::mp, Tenseur::poids, Tbl::set(), Cmp::set(), Tenseur::set(), Cmp::set_etat_qcq(), Tenseur::set_etat_qcq(), Tenseur::set_etat_zero(), Tenseur::set_triad(), Cmp::std_base_scal(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
void Tenseur::sauve | ( | FILE * | fd | ) | const [inherited] |
Save in a file.
Definition at line 1320 of file tenseur.C.
References Tenseur::c, Tenseur::etat, fwrite_be(), Tenseur::n_comp, Tenseur::poids, Base_vect::sauve(), Itbl::sauve(), Tenseur::triad, Tenseur::type_indice, and Tenseur::valence.
Read/write in the general case.
Definition at line 879 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), Tenseur::donne_place(), Tenseur::etat, Itbl::get_dim(), Itbl::get_ndim(), and Tenseur::valence.
Cmp & Tenseur::set | ( | int | ind1, | |
int | ind2, | |||
int | ind3 | |||
) | [inherited] |
Read/write for a tensor of valence 3.
Definition at line 859 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), Tenseur::donne_place(), Tenseur::etat, Itbl::set(), Itbl::set_etat_qcq(), and Tenseur::valence.
Cmp & Tenseur::set | ( | int | ind1, | |
int | ind2 | |||
) | [inherited] |
Read/write for a tensor of valence 2.
Definition at line 840 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), Tenseur::donne_place(), Tenseur::etat, Itbl::set(), Itbl::set_etat_qcq(), and Tenseur::valence.
Cmp & Tenseur::set | ( | int | ind | ) | [inherited] |
Read/write for a vector.
Definition at line 829 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), Tenseur::etat, and Tenseur::valence.
Cmp & Tenseur::set | ( | ) | [inherited] |
Read/write for a scalar (see also operator=
(const Cmp&
) ).
Definition at line 819 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), Tenseur::etat, and Tenseur::valence.
void Tenseur::set_dependance | ( | const Metrique & | 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 603 of file tenseur.C.
References Tenseur::met_depend.
void Tenseur::set_der_0x0 | ( | ) | const [protected, inherited] |
Sets the pointers of all the derivatives to zero.
Definition at line 586 of file tenseur.C.
References Tenseur::p_gradient, Tenseur::p_gradient_spher, and Tenseur::set_der_met_0x0().
void Tenseur::set_der_met_0x0 | ( | int | i | ) | const [protected, inherited] |
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 578 of file tenseur.C.
References Tenseur::met_depend, Tenseur::p_carre_scal, Tenseur::p_derive_con, and Tenseur::p_derive_cov.
void Tenseur::set_etat_nondef | ( | ) | [inherited] |
Sets the logical state to ETATNONDEF
(undefined state).
The components are not allocated.
Definition at line 645 of file tenseur.C.
References Tenseur::del_t(), and Tenseur::etat.
void 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 631 of file tenseur.C.
References Tenseur::c, Tenseur::del_derive(), Tenseur::etat, Tenseur::mp, and Tenseur::n_comp.
void Tenseur::set_etat_zero | ( | ) | [inherited] |
Sets the logical state to ETATZERO
(zero state).
The components are not allocated.
Definition at line 640 of file tenseur.C.
References Tenseur::del_t(), and Tenseur::etat.
void Tenseur::set_metric | ( | const Metrique & | met | ) | [inherited] |
Sets the pointer on the metric for a tensor density.
Definition at line 680 of file tenseur.C.
References Tenseur::metric.
void Tenseur::set_poids | ( | double | weight | ) | [inherited] |
Sets the weight for a tensor density.
Definition at line 675 of file tenseur.C.
References Tenseur::poids.
void 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 1165 of file tenseur.C.
References Tenseur::c, Tenseur::donne_indices(), Tenseur::etat, Map::get_bvect_cart(), Map::get_bvect_spher(), Map::get_mg(), Base_vect::identify(), Tenseur::mp, Tenseur::n_comp, Itbl::set_etat_qcq(), Cmp::std_base_scal(), Mg3d::std_base_vect_cart(), Mg3d::std_base_vect_spher(), Tenseur::triad, Cmp::va, and Tenseur::valence.
void 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 669 of file tenseur.C.
References Tenseur::triad.
bool Tenseur::verif | ( | ) | const [protected, inherited] |
Returns false for a tensor density without a defined metric.
Definition at line 187 of file tenseur.C.
References Tenseur::metric, and Tenseur::poids.
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 812 of file tenseur_operateur.C.
Cmp** Tenseur::c [protected, inherited] |
int Tenseur::etat [protected, inherited] |
const Metrique** Tenseur::met_depend [protected, inherited] |
const Metrique* Tenseur::metric [protected, inherited] |
const Map* const Tenseur::mp [protected, inherited] |
int Tenseur::n_comp [protected, inherited] |
Tenseur** Tenseur::p_carre_scal [protected, inherited] |
Tenseur** Tenseur::p_derive_con [protected, inherited] |
Tenseur** Tenseur::p_derive_cov [protected, inherited] |
Tenseur* Tenseur::p_gradient [mutable, protected, inherited] |
Tenseur* Tenseur::p_gradient_spher [mutable, protected, inherited] |
double Tenseur::poids [protected, inherited] |
const Base_vect* Tenseur::triad [protected, inherited] |
Itbl Tenseur::type_indice [protected, inherited] |
int Tenseur::valence [protected, inherited] |