32char tensor_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Tensor/tensor.C,v 1.43 2014/10/13 08:53:44 j_novak Exp $" ;
193#include "utilitaires.h"
204 : mp(&map), valence(val), triad(&
triad_i), type_indice(
tipe),
205 n_comp(
int(
pow(3., val))){
227 : mp(&map), valence(val), triad(
triad_i), type_indice(
tipe),
228 n_comp(
int(
pow(3., val))){
252 : mp(&map), valence(val), triad(&
triad_i), type_indice(val),
253 n_comp(
int(
pow(3., val))){
273 type_indice(
source.type_indice) {
323 type_indice(0), n_comp(1) {
361 mp(&map), valence(val), triad(&
triad_i), type_indice(val), n_comp(
compo)
400 for (
int i=0;
i<N_MET_MAX;
i++)
409 for (
int i=0;
i<N_MET_MAX;
i++)
419 for (
int i=0 ;
i<N_TENSOR_DEPEND ;
i++)
445 for (
int i=0;
i<N_MET_MAX;
i++)
457 for (
int i=0;
i<N_MET_MAX;
i++) {
461 if (
nmet == N_MET_MAX) {
462 cout <<
"Too many metrics in Tensor::set_dependances" <<
endl ;
467 while ((
conte < N_TENSOR_DEPEND) && (
met.tensor_depend[
conte] != 0x0))
470 if (
conte == N_TENSOR_DEPEND) {
471 cout <<
"Too many dependancies in Tensor::set_dependances " <<
endl ;
674 if ( (
l_min == 0) && (
l_max ==
mp->get_mg()->get_nzone()-1) ) {
695 int nz =
mp->get_mg()->get_nzone() ;
697 if ((2*
deg+1) >=
mp->get_mg()->get_nr(
lrac))
698 cout <<
"Tensor::annule_extern_cn : \n"
700 <<
"The number of coefficients in r is too low \n"
701 <<
"to do a clean matching..." <<
endl ;
709 binom.annule_hard() ;
711 for (
int n=1;
n<=
deg;
n++) {
713 for (
int k=1;
k<=
n;
k++)
722 for (
int i=
deg-1;
i>=0;
i--) {
731 double aa = coef(
deg) ;
732 for (
int i =
deg-1;
i>=0;
i--)
742 for (
int i=
deg-1;
i>=0;
i--)
751 rac.std_spectral_base() ;
829 flux <<
"Lorene class : " <<
typeid(
source).name()
830 <<
" Valence : " <<
source.valence <<
'\n' ;
832 if (
source.get_triad() != 0x0) {
833 flux <<
"Vectorial basis (triad) on which the components are defined :"
835 flux << *(
source.get_triad()) <<
'\n' ;
839 flux <<
"Type of the indices : " ;
840 for (
int i=0 ; i<source.
valence ; i++) {
841 flux <<
"index " << i <<
" : " ;
843 flux <<
" contravariant." <<
'\n' ;
845 flux <<
" covariant." <<
'\n' ;
846 if ( i < source.
valence-1 ) flux <<
" " ;
851 for (
int i=0 ; i<source.
n_comp ; i++) {
855 "===================== Scalar field ========================= \n" ;
858 flux <<
"================ Component " ;
859 Itbl num_indices = source.
indices(i) ;
860 for (
int j=0 ; j<source.
valence ; j++) {
861 flux <<
" " << num_indices(j) ;
863 flux <<
" ================ \n" ;
867 flux << *source.
cmp[i] <<
'\n' ;
881 ost <<
"Lorene class : " <<
typeid(*this).name()
882 <<
" Valence : " <<
valence <<
'\n' ;
888 "===================== Scalar field ========================= \n" ;
891 ost <<
"================ Component " ;
896 ost <<
" ================ \n" ;
937 "Tensor::std_spectral_base: should not be called on a Tensor"
938 <<
" of valence 1 but on a Vector !" <<
endl ;
946 if(
triad->identify() == (
mp->get_bvect_cart()).identify() ) {
947 bases =
mp->get_mg()->std_base_vect_cart() ;
950 assert(
triad->identify() == (
mp->get_bvect_spher()).identify()) ;
951 bases =
mp->get_mg()->std_base_vect_spher() ;
961 for (
int i=0 ;
i<3 ;
i++) {
972 cout <<
"Tensor::std_spectral_base: the case valence = " <<
valence
973 <<
" is not treated yet !" <<
endl ;
993 cout <<
"Tensor::std_spectral_base_odd: the case valence = " <<
valence
994 <<
" is not treated yet !" <<
endl ;
1033 tsym->sym_index1(),
tsym->sym_index2()) ;
1068 if(
triad->identify() == (
mp->get_bvect_cart()).identify() )
1072 cout <<
"Tensor::exponential_filter_r : " <<
endl ;
1073 cout <<
"Only Cartesian triad is implemented!" <<
endl ;
1081 if(
triad->identify() == (
mp->get_bvect_cart()).identify() )
1085 cout <<
"Tensor::exponential_filter_ylm : " <<
endl ;
1086 cout <<
"Only Cartesian triad is implemented!" <<
endl ;
Bases of the spectral expansions.
Vectorial bases (triads) with respect to which the tensorial components are defined.
static Base_vect * bvect_from_file(FILE *)
Construction of a vectorial basis from a file (see sauve(FILE* ) ).
Time evolution with partial storage (*** under development ***).
int position(int j) const
Gives the position in the arrays step, the_time and val corresponding to the time step j.
Basic integer array class.
void sauve(FILE *) const
Save in a file.
Base class for coordinate mappings.
Metric for tensor calculation.
Tensor field of valence 0 (or component of a tensorial field).
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
virtual void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
virtual void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
virtual void std_spectral_base_odd()
Sets the spectral bases of the Valeur va to the standard odd ones for a scalar field.
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 set_spectral_base(const Base_val &)
Sets the spectral bases of the Valeur va
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
double & set(int i)
Read/write of a particular element (index i) (1D case)
Symmetric tensors (with respect to two of their arguments).
Tensor field of valence 1.
Cmp pow(const Cmp &, int)
Power .
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
virtual void exponential_filter_ylm(int lzmin, int lzmax, int p, double alpha=-16.)
Applies exponential filters to all components (see Scalar::exponential_filter_ylm ).
virtual void operator=(const Tensor &)
Assignment to a Tensor.
virtual void sauve(FILE *) const
Save in a binary file.
virtual void del_derive_met(int) const
Logical destructor of the derivatives depending on the i-th element of met_depend .
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values in the compactified ext...
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
virtual void allocate_all()
Performs the memory allocation of all the elements, down to the double arrays of the Tbl s.
virtual void std_spectral_base_odd()
Sets the standard odd spectal bases of decomposition for each component.
Tensor * p_divergence[N_MET_MAX]
Array of pointers on the divergence of this with respect to various metrics.
void annule_extern_cn(int l_0, int deg)
Performs a smooth (C^n) transition in a given domain to zero.
int get_place_met(const Metric &) const
Returns the position of the pointer on metre in the array met_depend .
void set_der_0x0() const
Sets the pointers on derived quantities to 0x0.
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 of each component.
const Tensor & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
Tensor * p_derive_cov[N_MET_MAX]
Array of pointers on the covariant derivatives of this with respect to various metrics.
const Metric * met_depend[N_MET_MAX]
Array on the Metric 's which were used to compute derived quantities, like p_derive_cov ,...
virtual Itbl indices(int pos) const
Returns the indices of a component given by its position in the array cmp .
virtual int position(const Itbl &ind) const
Returns the position in the array cmp of a component given by its indices.
Tensor(const Map &map, int val, const Itbl &tipe, const Base_vect &triad_i)
Standard constructor.
int valence
Valence of the tensor (0 = scalar, 1 = vector, etc...)
virtual void set_etat_nondef()
Sets the logical state of all components to ETATNONDEF (undefined state).
virtual void annule(int l_min, int l_max)
Sets the Tensor to zero in several domains.
virtual ~Tensor()
Destructor.
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
const Scalar & operator()(const Itbl &ind) const
Returns the value of a component (read-only version).
void set_dependance(const Metric &) const
To be used to describe the fact that the derivatives members have been calculated with met .
Scalar ** cmp
Array of size n_comp of pointers onto the components.
int n_comp
Number of stored components, depending on the symmetry.
void operator-=(const Tensor &)
-= Tensor
void set_der_met_0x0(int) const
Sets all the i-th components of met_depend , p_derive_cov , etc... to 0x0.
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
virtual void del_deriv() const
Deletes the derived quantities.
void operator+=(const Tensor &)
+= Tensor
const Tensor & divergence(const Metric &gam) const
Computes the divergence of this with respect to some metric .
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
const Tensor & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
Itbl type_indice
1D array of integers (class Itbl ) of size valence containing the type of each index: COV for a cov...
virtual void exponential_filter_r(int lzmin, int lzmax, int p, double alpha=-16.)
Applies exponential filters to all components (see Scalar::exponential_filter_r ).
Tensor * p_derive_con[N_MET_MAX]
Array of pointers on the contravariant derivatives of this with respect to various metrics.
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .