sym_tensor_trans.C
/*
* Methods of class Sym_tensor_trans
*
* (see file sym_tensor.h for documentation)
*
*/
/*
* Copyright (c) 2003-2004 Eric Gourgoulhon & Jerome Novak
*
* This file is part of LORENE.
*
* LORENE is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* LORENE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with LORENE; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
*/
char sym_tensor_trans_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/sym_tensor_trans.C,v 1.18 2014/10/13 08:53:43 j_novak Exp $" ;
/*
* $Id: sym_tensor_trans.C,v 1.18 2014/10/13 08:53:43 j_novak Exp $
* $Log: sym_tensor_trans.C,v $
* Revision 1.18 2014/10/13 08:53:43 j_novak
* Lorene classes and functions now belong to the namespace Lorene.
*
* Revision 1.17 2014/10/06 15:13:19 j_novak
* Modified #include directives to use c++ syntax.
*
* Revision 1.16 2008/12/03 10:20:00 j_novak
* Modified output.
*
* Revision 1.15 2006/09/15 08:48:01 j_novak
* Suppression of some messages in the optimized version.
*
* Revision 1.14 2005/01/03 12:37:08 j_novak
* Sym_tensor_trans::trace_from_det_one : modified the test on hijtt to
* be compatible with older compilers.
*
* Revision 1.13 2005/01/03 08:36:36 f_limousin
* Come back to the previous version.
*
* Revision 1.12 2005/01/03 08:13:26 f_limousin
* The first argument of the function trace_from_det_one(...) is now
* a Sym_tensor_trans instead of a Sym_tensor_tt (because of a
* compilation error with some compilators).
*
* Revision 1.11 2004/12/28 14:21:48 j_novak
* Added the method Sym_tensor_trans::trace_from_det_one
*
* Revision 1.10 2004/05/25 15:07:12 f_limousin
* Add parameters in argument of the function tt_part for the case
* of a Map_et.
*
* Revision 1.9 2004/03/30 14:01:19 j_novak
* Copy constructors and operator= now copy the "derived" members.
*
* Revision 1.8 2004/03/30 08:01:16 j_novak
* Better treatment of dzpuis in mutators.
*
* Revision 1.7 2004/03/29 16:13:07 j_novak
* New methods set_longit_trans and set_tt_trace .
*
* Revision 1.6 2004/03/03 13:22:14 j_novak
* The case where dzpuis = 0 is treated in tt_part().
*
* Revision 1.5 2004/02/18 18:48:39 e_gourgoulhon
* Method trace() renamed the_trace() to avoid any confusion with
* the new method Tensor::trace().
*
* Revision 1.4 2004/02/09 12:57:13 e_gourgoulhon
* First implementation of method tt_part().
*
* Revision 1.3 2004/01/04 20:52:45 e_gourgoulhon
* Added assignement (operator=) to a Tensor_sym.
*
* Revision 1.2 2003/10/28 21:24:52 e_gourgoulhon
* Added new methods trace() and tt_part().
*
* Revision 1.1 2003/10/27 10:50:54 e_gourgoulhon
* First version.
*
*
*
* $Header: /cvsroot/Lorene/C++/Source/Tensor/sym_tensor_trans.C,v 1.18 2014/10/13 08:53:43 j_novak Exp $
*
*/
// Headers C
#include <cstdlib>
// Headers Lorene
#include "metric.h"
#include "param.h"
//--------------//
// Constructors //
//--------------//
// Standard constructor
// --------------------
namespace Lorene {
Sym_tensor_trans::Sym_tensor_trans(const Map& map, const Base_vect& triad_i,
const Metric& met)
: Sym_tensor(map, CON, triad_i),
met_div(&met) {
set_der_0x0() ;
}
// Copy constructor
// ----------------
Sym_tensor_trans::Sym_tensor_trans (const Sym_tensor_trans& source)
: Sym_tensor(source),
met_div(source.met_div) {
set_der_0x0() ;
if (source.p_trace != 0x0) p_trace = new Scalar( *(source.p_trace) ) ;
if (source.p_tt != 0x0) p_tt = new Sym_tensor_tt( *(source.p_tt) ) ;
}
// Constructor from a file
// -----------------------
Sym_tensor_trans::Sym_tensor_trans(const Map& mapping, const Base_vect& triad_i,
const Metric& met, FILE* fd)
: Sym_tensor(mapping, triad_i, fd),
met_div(&met) {
set_der_0x0() ;
}
//--------------//
// Destructor //
//--------------//
Sym_tensor_trans::~Sym_tensor_trans() {
Sym_tensor_trans::del_deriv() ; // in order not to follow the virtual aspect
// of del_deriv()
}
//-------------------//
// Memory managment //
//-------------------//
void Sym_tensor_trans::del_deriv() const {
if (p_trace != 0x0) delete p_trace ;
if (p_tt != 0x0) delete p_tt ;
set_der_0x0() ;
Sym_tensor::del_deriv() ;
}
void Sym_tensor_trans::set_der_0x0() const {
p_trace = 0x0 ;
p_tt = 0x0 ;
}
//--------------//
// Assignment //
//--------------//
void Sym_tensor_trans::operator=(const Sym_tensor_trans& source) {
// Assignment of quantities common to all the derived classes of Sym_tensor
Sym_tensor::operator=(source) ;
// Assignment of proper quantities of class Sym_tensor_trans
assert(met_div == source.met_div) ;
del_deriv() ;
if (source.p_trace != 0x0) p_trace = new Scalar( *(source.p_trace) ) ;
if (source.p_tt != 0x0) p_tt = new Sym_tensor_tt( *(source.p_tt) ) ;
}
void Sym_tensor_trans::operator=(const Sym_tensor& source) {
// Assignment of quantities common to all the derived classes of Sym_tensor
Sym_tensor::operator=(source) ;
// The metric which was set by the constructor is kept
del_deriv() ;
}
void Sym_tensor_trans::operator=(const Tensor_sym& source) {
// Assignment of quantities common to all the derived classes of Sym_tensor
Sym_tensor::operator=(source) ;
// The metric which was set by the constructor is kept
del_deriv() ;
}
void Sym_tensor_trans::operator=(const Tensor& source) {
// Assignment of quantities common to all the derived classes of Sym_tensor
Sym_tensor::operator=(source) ;
// The metric which was set by the constructor is kept
del_deriv() ;
}
void Sym_tensor_trans::set_tt_trace(const Sym_tensor_tt& htt,
const Scalar& htrace, Param* par ) {
assert (met_div == &htt.get_met_div() ) ;
assert (htrace.check_dzpuis(4)) ;
Scalar pot (*mp) ;
if (dynamic_cast<const Map_af*>(mp) != 0x0) {
pot = htrace.poisson() ;
}
else {
pot = 0. ;
htrace.poisson(*par, pot) ;
}
Sym_tensor tmp = (pot.derive_con(*met_div)).derive_con(*met_div) ;
tmp.dec_dzpuis() ;
*this = htrace * met_div->con() ;
dec_dzpuis(2) ; // this has now dzpuis = 2
*this = htt + 0.5 * ( *this - tmp ) ;
del_deriv() ;
p_trace = new Scalar( htrace ) ;
p_tt = new Sym_tensor_tt( htt ) ;
}
//-----------------------------//
// Computational methods //
//-----------------------------//
const Scalar& Sym_tensor_trans::the_trace() const {
if (p_trace == 0x0) { // a new computation is necessary
assert( (type_indice(0)==CON) && (type_indice(1)==CON) ) ;
p_trace = new Scalar( trace(*met_div) ) ;
}
return *p_trace ;
}
const Sym_tensor_tt& Sym_tensor_trans::tt_part(Param* par) const {
if (p_tt == 0x0) { // a new computation is necessary
int dzp = the_trace().get_dzpuis() ;
assert((dzp == 0) || (dzp == 4)) ;
p_tt = new Sym_tensor_tt(*mp, *triad, *met_div) ;
Scalar pot (*mp) ;
if (dynamic_cast<const Map_af*>(mp) != 0x0) {
pot = the_trace().poisson() ;
}
else {
pot = 0. ;
the_trace().poisson(*par, pot) ;
}
Sym_tensor tmp = (pot.derive_con(*met_div)).derive_con(*met_div) ;
(dzp == 4) ? tmp.inc_dzpuis() : tmp.dec_dzpuis(3) ; //## to be improved ?
*p_tt = *this - 0.5 * ( the_trace() * met_div->con() - tmp ) ;
}
return *p_tt ;
}
void Sym_tensor_trans::trace_from_det_one(const Sym_tensor_tt& hijtt,
double precis, int it_max) {
#ifndef NDEBUG
const Sym_tensor_trans* ptmp =
dynamic_cast<const Sym_tensor_trans*>(&hijtt) ;
assert (ptmp != 0x0) ;
assert (ptmp != this) ;
for (int i=0; i<hijtt.get_n_comp(); i++)
assert(hijtt.cmp[i]->check_dzpuis(2)) ;
#endif
assert( (precis > 0.) && (it_max > 0) ) ;
assert (met_div == &hijtt.get_met_div() ) ;
Sym_tensor_trans& hij = *this ;
hij = hijtt ; //initialization
// The trace h = f_{ij} h^{ij} :
Scalar htrace(*mp) ;
// Value of h at previous step of the iterative procedure below :
Scalar htrace_prev(*mp) ;
htrace_prev.set_etat_zero() ; // initialisation to zero
for (int it=0; it<=it_max; it++) {
// Trace h from the condition det(f^{ij} + h^{ij}) = det f^{ij} :
htrace = hij(1,1) * hij(2,3) * hij(2,3)
+ hij(2,2) * hij(1,3) * hij(1,3) + hij(3,3) * hij(1,2) * hij(1,2)
- 2.* hij(1,2) * hij(1,3) * hij(2,3)
- hij(1,1) * hij(2,2) * hij(3,3) ;
htrace.dec_dzpuis(2) ; // dzpuis: 6 --> 4
htrace += hij(1,2) * hij(1,2) + hij(1,3) * hij(1,3)
+ hij(2,3) * hij(2,3) - hij(1,1) * hij(2,2)
- hij(1,1) * hij(3,3) - hij(2,2) * hij(3,3) ;
// New value of hij from htrace and hijtt (obtained by solving
// the Poisson equation for Phi) :
set_tt_trace(hijtt, htrace) ;
double diff = max(max(abs(htrace - htrace_prev))) ;
#ifndef NDEBUG
cout << "Sym_tensor_trans::trace_from_det_one : "
<< "iteration : " << it << " convergence on trace(h): " << diff << endl ;
#endif
if (diff < precis) break ;
else htrace_prev = htrace ;
if (it == it_max) {
cout << "Sym_tensor_trans::trace_from_det_one : convergence not reached \n" ;
cout << " for the required accuracy (" << precis << ") ! " << endl ;
abort() ;
}
}
}
}
Lorene::Base_vect
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition base_vect.h:105
Lorene::Map_af
Affine radial mapping.
Definition map.h:2027
Lorene::Param
Parameter storage.
Definition param.h:125
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Lorene::Scalar::poisson
Scalar poisson() const
Solves the scalar Poisson equation with *this as a source.
Definition scalar_pde.C:136
Lorene::Scalar::get_dzpuis
int get_dzpuis() const
Returns dzpuis.
Definition scalar.h:557
Lorene::Scalar::check_dzpuis
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 n...
Definition scalar.C:873
Lorene::Scalar::set_etat_zero
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
Lorene::Scalar::derive_con
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Definition scalar_deriv.C:402
Lorene::Scalar::dec_dzpuis
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 co...
Definition scalar_r_manip.C:418
Lorene::Sym_tensor_trans::p_tt
Sym_tensor_tt * p_tt
Traceless part with respect to the metric *met_div.
Definition sym_tensor.h:620
Lorene::Sym_tensor_trans::Sym_tensor_trans
Sym_tensor_trans(const Map &map, const Base_vect &triad_i, const Metric &met)
Standard constructor.
Definition sym_tensor_trans.C:114
Lorene::Sym_tensor_trans::trace_from_det_one
void trace_from_det_one(const Sym_tensor_tt &htt, double precis=1.e-14, int it_max=100)
Assigns the derived member p_tt and computes the trace so that *this + the flat metric has a determin...
Definition sym_tensor_trans.C:315
Lorene::Sym_tensor_trans::met_div
const Metric *const met_div
Metric with respect to which the divergence and the trace are defined.
Definition sym_tensor.h:614
Lorene::Sym_tensor_trans::p_trace
Scalar * p_trace
Trace with respect to the metric *met_div.
Definition sym_tensor.h:617
Lorene::Sym_tensor_trans::the_trace
const Scalar & the_trace() const
Returns the trace of the tensor with respect to metric *met_div.
Definition sym_tensor_trans.C:270
Lorene::Sym_tensor_trans::set_der_0x0
void set_der_0x0() const
Sets the pointers on derived quantities to 0x0.
Definition sym_tensor_trans.C:174
Lorene::Sym_tensor_trans::get_met_div
const Metric & get_met_div() const
Returns the metric with respect to which the divergence and the trace are defined.
Definition sym_tensor.h:669
Lorene::Sym_tensor_trans::tt_part
const Sym_tensor_tt & tt_part(Param *par=0x0) const
Returns the transverse traceless part of the tensor, the trace being defined with respect to metric *...
Definition sym_tensor_trans.C:284
Lorene::Sym_tensor_trans::~Sym_tensor_trans
virtual ~Sym_tensor_trans()
Destructor.
Definition sym_tensor_trans.C:151
Lorene::Sym_tensor_trans::set_tt_trace
void set_tt_trace(const Sym_tensor_tt &a, const Scalar &h, Param *par=0x0)
Assigns the derived members p_tt and p_trace and updates the components accordingly.
Definition sym_tensor_trans.C:235
Lorene::Sym_tensor_trans::operator=
virtual void operator=(const Sym_tensor_trans &a)
Assignment to another Sym_tensor_trans.
Definition sym_tensor_trans.C:186
Lorene::Sym_tensor_trans::del_deriv
virtual void del_deriv() const
Deletes the derived quantities.
Definition sym_tensor_trans.C:164
Lorene::Sym_tensor_tt
Transverse and traceless symmetric tensors of rank 2.
Definition sym_tensor.h:938
Lorene::Sym_tensor::operator=
virtual void operator=(const Sym_tensor &a)
Assignment to another Sym_tensor.
Definition sym_tensor.C:234
Lorene::Sym_tensor::del_deriv
virtual void del_deriv() const
Deletes the derived quantities.
Definition sym_tensor.C:286
Lorene::Sym_tensor::Sym_tensor
Sym_tensor(const Map &map, const Itbl &tipe, const Base_vect &triad_i)
Standard constructor.
Definition sym_tensor.C:142
Lorene::Tensor_sym
Symmetric tensors (with respect to two of their arguments).
Definition tensor.h:1037
Lorene::Tensor
Tensor handling.
Definition tensor.h:288
Lorene::max
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition cmp_math.C:435
Lorene::abs
Cmp abs(const Cmp &)
Absolute value.
Definition cmp_math.C:410
Lorene::Tensor::dec_dzpuis
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition tensor.C:808
Lorene::Tensor::mp
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition tensor.h:295
Lorene::Tensor::trace
Scalar trace() const
Trace on two different type indices for a valence 2 tensor.
Definition tensor_calculus.C:180
Lorene::Tensor::get_n_comp
int get_n_comp() const
Returns the number of stored components.
Definition tensor.h:872
Lorene::Tensor::inc_dzpuis
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition tensor.C:816
Lorene::Tensor::cmp
Scalar ** cmp
Array of size n_comp of pointers onto the components.
Definition tensor.h:315
Lorene::Tensor::type_indice
Itbl type_indice
1D array of integers (class Itbl ) of size valence containing the type of each index: COV for a cov...
Definition tensor.h:310
Lorene::Tensor::triad
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Definition tensor.h:303
Lorene
Lorene prototypes.
Definition app_hor.h:64
Lorene::Map
Map(const Mg3d &)
Constructor from a multi-domain 3D grid.
Definition map.C:139
Lorene::dec_dzpuis
virtual void dec_dzpuis(Scalar &) const =0
Decreases by 1 the value of dzpuis of a Scalar and changes accordingly its values in the compacti...