LORENE
sym_tensor_decomp.C
1/*
2 * Methods transverse( ) and longit_pot( ) of class Sym_tensor
3 *
4 * (see file sym_tensor.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2003-2004 Eric Gourgoulhon & Jerome Novak
10 *
11 * This file is part of LORENE.
12 *
13 * LORENE is free software; you can redistribute it and/or modify
14 * it under the terms of the GNU General Public License version 2
15 * as published by the Free Software Foundation.
16 *
17 * LORENE is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
21 *
22 * You should have received a copy of the GNU General Public License
23 * along with LORENE; if not, write to the Free Software
24 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 *
26 */
27
28char sym_tensor_decomp_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/sym_tensor_decomp.C,v 1.14 2014/10/13 08:53:43 j_novak Exp $" ;
29
30/*
31 * $Id: sym_tensor_decomp.C,v 1.14 2014/10/13 08:53:43 j_novak Exp $
32 * $Log: sym_tensor_decomp.C,v $
33 * Revision 1.14 2014/10/13 08:53:43 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.13 2008/12/05 08:45:12 j_novak
37 * Modified dzpuis treatment.
38 *
39 * Revision 1.12 2008/12/03 10:20:00 j_novak
40 * Modified output.
41 *
42 * Revision 1.11 2004/06/17 06:56:42 e_gourgoulhon
43 * Simplified code for method transverse (use of Vector::ope_killing).
44 * Slight modif. of output in method longit_pot.
45 *
46 * Revision 1.10 2004/06/14 20:45:41 e_gourgoulhon
47 * Added argument method_poisson to methods longit_pot and transverse.
48 *
49 * Revision 1.9 2004/05/25 15:05:57 f_limousin
50 * Add parameters in argument of the functions transverse, longit_pot
51 * for the case of Map_et.
52 *
53 * Revision 1.8 2004/03/30 14:01:19 j_novak
54 * Copy constructors and operator= now copy the "derived" members.
55 *
56 * Revision 1.7 2004/03/30 08:01:16 j_novak
57 * Better treatment of dzpuis in mutators.
58 *
59 * Revision 1.6 2004/03/29 16:13:07 j_novak
60 * New methods set_longit_trans and set_tt_trace .
61 *
62 * Revision 1.5 2004/02/09 12:56:27 e_gourgoulhon
63 * Method longit_pot: added test of the vector Poisson equation.
64 *
65 * Revision 1.4 2004/02/02 09:18:11 e_gourgoulhon
66 * Method longit_pot: treatment of case divergence dzpuis = 5.
67 *
68 * Revision 1.3 2003/12/10 10:17:54 e_gourgoulhon
69 * First operational version.
70 *
71 * Revision 1.2 2003/11/27 16:01:47 e_gourgoulhon
72 * First implmentation.
73 *
74 * Revision 1.1 2003/11/26 21:57:03 e_gourgoulhon
75 * First version; not ready yet.
76 *
77 *
78 * $Header: /cvsroot/Lorene/C++/Source/Tensor/sym_tensor_decomp.C,v 1.14 2014/10/13 08:53:43 j_novak Exp $
79 *
80 */
81
82
83// Lorene headers
84#include "metric.h"
85#include "param.h"
86
87namespace Lorene {
89 const Sym_tensor_trans& ht ) {
90
91 assert ( v_pot.get_index_type(0) == CON ) ;
92
93 const Metric& metre = ht.get_met_div() ;
94
95 *this = ht + v_pot.ope_killing(metre) ; // this has dzpuis = 2, if v_pot not 0
96 if ((*this)(1,1).get_dzpuis() == 2)
97 dec_dzpuis(2) ; // which is decreased so to add *this to a flat metric
98
99 del_deriv() ;
100
102 int jp = get_place_met(metre) ;
103 assert ((jp>=0) && (jp<N_MET_MAX)) ;
104
106 p_longit_pot[jp] = new Vector( v_pot ) ;
107
108}
109
111 int method_poisson) const {
112
114 int jp = get_place_met(metre) ;
115 assert ((jp>=0) && (jp<N_MET_MAX)) ;
116
117 if (p_transverse[jp] == 0x0) { // a new computation is necessary
118
119 assert( (type_indice(0) == CON) && (type_indice(1) == CON) ) ;
120
121 for (int ic=0; ic<n_comp; ic++) {
122 assert(cmp[ic]->check_dzpuis(4)) ; // dzpuis=4 is assumed
123 }
124
125 const Vector& ww = longit_pot(metre, par, method_poisson) ;
126
127 Sym_tensor lww = ww.ope_killing(metre) ; // D^i W^j + D^j W^i
128
129 lww.inc_dzpuis(2) ;
130
132
133 *(p_transverse[jp]) = *this - lww ;
134
135 }
136
137 return *p_transverse[jp] ;
138
139
140}
141
142
144 int method_poisson) const {
145
147 int jp = get_place_met(metre) ;
148 assert ((jp>=0) && (jp<N_MET_MAX)) ;
149
150 if (p_longit_pot[jp] == 0x0) { // a new computation is necessary
151
152 const Metric_flat* metf = dynamic_cast<const Metric_flat*>(&metre) ;
153 if (metf == 0x0) {
154 cout << "Sym_tensor::longit_pot : the case of a non flat metric"
155 << endl <<" is not treated yet !" << endl ;
156 abort() ;
157 }
158
159 Vector hhh = divergence(metre) ;
160
161
162 // If dpzuis is 5, it should be decreased to 4 for the Poisson equation:
163 bool dzp5 = false ;
164 for (int i=1; i<=3; i++) {
165 dzp5 = dzp5 || hhh(i).check_dzpuis(5) ;
166 }
167 if (dzp5) hhh.dec_dzpuis() ;
168
169 if (dynamic_cast<const Map_af*>(mp) != 0x0)
170 p_longit_pot[jp] = new Vector( hhh.poisson(double(1),
171 method_poisson) ) ;
172 else
173 p_longit_pot[jp] = new Vector( hhh.poisson(double(1), *par,
174 method_poisson) ) ;
175
176
177 // Test of resolution of the vector Poisson equation:
178 const Vector& ww = *(p_longit_pot[jp]) ;
179
180 hhh.dec_dzpuis() ;
181
182 Vector lapw = ww.derive_con(metre).divergence(metre)
184#ifndef NDEBUG
185 cout << "## Sym_tensor::longit_pot : test of Poisson : \n" ;
186 cout <<
187 " Max absolute error in each domain on the vector Poisson equation:\n" ;
188 maxabs(lapw - hhh) ;
189
190 int nz = mp->get_mg()->get_nzone() ; // total number of domains
191
192 cout << " Relative error in each domain on the vector Poisson equation:\n" ;
193 for (int i=1; i<=3; i++){
194 cout << " Comp. " << i << " : " ;
195 for (int l=0; l<nz; l++){
196 cout << diffrel(lapw(i),hhh(i) )(l) << " " ;
197 }
198 cout << endl ;
199 }
200 cout << endl ;
201#endif
202 }
203
204 return *p_longit_pot[jp] ;
205
206
207}
208
209
210}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Affine radial mapping.
Definition map.h:2027
Flat metric for tensor calculation.
Definition metric.h:261
Metric for tensor calculation.
Definition metric.h:90
Parameter storage.
Definition param.h:125
Transverse symmetric tensors of rank 2.
Definition sym_tensor.h:608
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
void set_longit_trans(const Vector &v, const Sym_tensor_trans &a)
Assigns the derived members p_longit_pot and p_transverse and updates the components accordingly.
const Sym_tensor_trans & transverse(const Metric &gam, Param *par=0x0, int method_poisson=6) const
Computes the transverse part of the tensor with respect to a given metric, transverse meaning diverg...
const Vector & longit_pot(const Metric &gam, Param *par=0x0, int method_poisson=6) const
Computes the vector potential of longitudinal part of the tensor (see documentation of method transv...
Sym_tensor_trans * p_transverse[N_MET_MAX]
Array of the transverse part of the tensor with respect to various metrics, transverse meaning diver...
Definition sym_tensor.h:239
virtual void del_deriv() const
Deletes the derived quantities.
Definition sym_tensor.C:286
Vector * p_longit_pot[N_MET_MAX]
Array of the vector potential of the longitudinal part of the tensor with respect to various metrics ...
Definition sym_tensor.h:246
const Vector & divergence(const Metric &) const
Returns the divergence of this with respect to a Metric .
Definition sym_tensor.C:349
Tensor field of valence 1.
Definition vector.h:188
Sym_tensor ope_killing(const Metric &gam) const
Computes the Killing operator associated with a given metric.
Definition vector.C:438
Vector poisson(double lambda, int method=6) const
Solves the vector Poisson equation with *this as a source.
const Scalar & divergence(const Metric &) const
The divergence of this with respect to a Metric .
Definition vector.C:381
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition cmp_math.C:504
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
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition tensor.h:295
int get_place_met(const Metric &) const
Returns the position of the pointer on metre in the array met_depend .
Definition tensor.C:443
const Tensor & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
Definition tensor.C:1014
void set_dependance(const Metric &) const
To be used to describe the fact that the derivatives members have been calculated with met .
Definition tensor.C:453
Scalar ** cmp
Array of size n_comp of pointers onto the components.
Definition tensor.h:315
int n_comp
Number of stored components, depending on the symmetry.
Definition tensor.h:312
const Tensor_sym & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
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
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Definition tensor.h:303
Tbl maxabs(const Tensor &aa, const char *comment=0x0, ostream &ost=cout, bool verb=true)
Maxima in each domain of the absolute values of the tensor components.
Lorene prototypes.
Definition app_hor.h:64