LORENE
tslice_adm_mass.C
1 /*
2  * Virtual methods of class Time_slice and derived classes to
3  * compute the ADM mass
4  *
5  * (see file time_slice.h for documentation).
6  *
7  */
8 
9 /*
10  * Copyright (c) 2004 Eric Gourgoulhon & Jerome Novak
11  *
12  * This file is part of LORENE.
13  *
14  * LORENE is free software; you can redistribute it and/or modify
15  * it under the terms of the GNU General Public License version 2
16  * as published by the Free Software Foundation.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 char tslice_adm_mass_C[] = "$Header: /cvsroot/Lorene/C++/Source/Time_slice/tslice_adm_mass.C,v 1.7 2014/10/13 08:53:47 j_novak Exp $" ;
30 
31 /*
32  * $Id: tslice_adm_mass.C,v 1.7 2014/10/13 08:53:47 j_novak Exp $
33  * $Log: tslice_adm_mass.C,v $
34  * Revision 1.7 2014/10/13 08:53:47 j_novak
35  * Lorene classes and functions now belong to the namespace Lorene.
36  *
37  * Revision 1.6 2014/10/06 15:13:21 j_novak
38  * Modified #include directives to use c++ syntax.
39  *
40  * Revision 1.5 2010/10/20 07:58:09 j_novak
41  * Better implementation of the explicit time-integration. Not fully-tested yet.
42  *
43  * Revision 1.4 2008/12/04 19:36:40 j_novak
44  * Removed old tests.
45  *
46  * Revision 1.3 2004/05/12 15:24:20 e_gourgoulhon
47  * Reorganized the #include 's, taking into account that
48  * time_slice.h contains now an #include "metric.h".
49  *
50  * Revision 1.2 2004/05/10 09:11:20 e_gourgoulhon
51  * Corrected bug in Time_slice::adm_mass().
52  * Reorganized outputs.
53  *
54  * Revision 1.1 2004/05/09 20:56:29 e_gourgoulhon
55  * First version.
56  *
57  *
58  * $Header: /cvsroot/Lorene/C++/Source/Time_slice/tslice_adm_mass.C,v 1.7 2014/10/13 08:53:47 j_novak Exp $
59  *
60  */
61 
62 
63 // C headers
64 #include <cmath>
65 
66 // Lorene headers
67 #include "time_slice.h"
68 
69 //--------------------
70 // Time_slice version
71 //--------------------
72 
73 namespace Lorene {
74 double Time_slice::adm_mass() const {
75 
76  if ( !(adm_mass_evol).is_known(jtime) ) { // a new computation is necessary
77 
78  const Map& mp = gam_dd().get_mp() ;
79  Metric_flat ff(mp, *(gam_dd().get_triad())) ;
80  int nz = mp.get_mg()->get_nzone() ;
81  Tbl* tmass = new Tbl(nz) ;
82  tmass->set_etat_qcq() ;
83 
84  Vector ww = gam_dd().derive_con(ff).trace(1,2).up(0,ff)
85  - gam_dd().trace(ff).derive_con(ff) ;
86 
87  for (int l=0; l<nz; l++) {
88  double radius = mp.val_r(l, 1., 0., 0.) ;
89  tmass->set(l) = ww.flux(radius, ff) / (16.* M_PI) ;
90  }
91 
92  adm_mass_evol.update(*tmass, jtime, the_time[jtime]) ;
93 
94  delete tmass ;
95 
96  cout << "Time_slice::adm_mass : " << adm_mass_evol[jtime] << endl ;
97  }
98 
99  const Tbl& tadm = adm_mass_evol[jtime] ;
100  return tadm(tadm.get_taille()-1) ;
101 }
102 
103 
104 //--------------------------
105 // Time_slice_conf version
106 //--------------------------
107 
108 double Time_slice_conf::adm_mass() const {
109 
110  if ( !(adm_mass_evol).is_known(jtime) ) { // a new computation is necessary
111 
112  const Map& mp = psi().get_mp() ;
113  int nz = mp.get_mg()->get_nzone() ;
114  Tbl* tmass = new Tbl(nz) ;
115  tmass->set_etat_qcq() ;
116 
117  Vector ww = psi().derive_con(ff)
118  + 0.125* ( hdirac() - (hh().trace(ff)).derive_con(ff) ) ;
119 
120  for (int l=0; l<nz; l++) {
121  double radius = mp.val_r(l, 1., 0., 0.) ;
122  tmass->set(l) = - ww.flux(radius, ff) / (2.* M_PI) ;
123  }
124 
125  adm_mass_evol.update(*tmass, jtime, the_time[jtime]) ;
126 
127  delete tmass ;
128 
129 #ifndef NDEBUG
130  cout << "Time_slice_conf::adm_mass : " << adm_mass_evol[jtime] << endl ;
131 #endif
132  }
133 
134  const Tbl& tadm = adm_mass_evol[jtime] ;
135  return tadm(tadm.get_taille()-1) ;
136 }
137 
138 
139 //--------------------------
140 // Tslice_dirac_max version
141 //--------------------------
142 
144 
145  if ( !(adm_mass_evol).is_known(jtime) ) { // a new computation is necessary
146 
147  const Map& mp = psi().get_mp() ;
148  int nz = mp.get_mg()->get_nzone() ;
149  Tbl* tmass = new Tbl(nz) ;
150  tmass->set_etat_qcq() ;
151 
152  Vector ww = psi().derive_con(ff)
153  - 0.125* (hh().trace(ff)).derive_con(ff) ;
154  // trh() is not used since it has dzpuis = 4
155 
156  for (int l=0; l<nz; l++) {
157  double radius = mp.val_r(l, 1., 0., 0.) ;
158  tmass->set(l) = - ww.flux(radius, ff) / (2.* M_PI) ;
159  }
160 
161  adm_mass_evol.update(*tmass, jtime, the_time[jtime]) ;
162 
163  delete tmass ;
164 #ifndef NDEBUG
165  cout << "Tslice_dirac_max::adm_mass : " << adm_mass_evol[jtime]
166  << endl ;
167 #endif
168 
169  }
170 
171  const Tbl& tadm = adm_mass_evol[jtime] ;
172  return tadm(tadm.get_taille()-1) ;
173 }
174 
175 
176 
177 
178 }
double flux(double radius, const Metric &met) const
Computes the flux of the vector accross a sphere r = const.
Definition: vector.C:807
virtual double adm_mass() const
Returns the ADM mass at (geometrical units) the current step.
Lorene prototypes.
Definition: app_hor.h:64
Evolution_full< Tbl > adm_mass_evol
ADM mass at each time step, since the creation of the slice.
Definition: time_slice.h:233
Tensor up(int ind, const Metric &gam) const
Computes a new tensor by raising an index of *this.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
double & set(int i)
Read/write of a particular element (index i) (1D case)
Definition: tbl.h:281
Flat metric for tensor calculation.
Definition: metric.h:261
Base class for coordinate mappings.
Definition: map.h:670
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...
int jtime
Time step index of the latest slice.
Definition: time_slice.h:190
Tensor field of valence 1.
Definition: vector.h:188
virtual const Sym_tensor & gam_dd() const
Induced metric (covariant components ) at the current time step (jtime )
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:361
virtual double val_r(int l, double xi, double theta, double pphi) const =0
Returns the value of the radial coordinate r for a given in a given domain.
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
virtual double adm_mass() const
Returns the ADM mass (geometrical units) at the current step.
Evolution_std< double > the_time
Time label of each slice.
Definition: time_slice.h:193
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
int get_taille() const
Gives the total size (ie dim.taille)
Definition: tbl.h:397
Basic array class.
Definition: tbl.h:161
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:861
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.