LORENE
time_slice_access.C
1/*
2 * Methods of class Time_slice to access the various fields
3 *
4 * (see file time_slice.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2004 Eric Gourgoulhon, Jose Luis Jaramillo & 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 time_slice_access_C[] = "$Header: /cvsroot/Lorene/C++/Source/Time_slice/time_slice_access.C,v 1.9 2014/10/13 08:53:47 j_novak Exp $" ;
29
30/*
31 * $Id: time_slice_access.C,v 1.9 2014/10/13 08:53:47 j_novak Exp $
32 * $Log: time_slice_access.C,v $
33 * Revision 1.9 2014/10/13 08:53:47 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.8 2014/10/06 15:13:21 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.7 2008/12/02 15:02:22 j_novak
40 * Implementation of the new constrained formalism, following Cordero et al. 2009
41 * paper. The evolution eqs. are solved as a first-order system. Not tested yet!
42 *
43 * Revision 1.6 2004/05/12 15:24:20 e_gourgoulhon
44 * Reorganized the #include 's, taking into account that
45 * time_slice.h contains now an #include "metric.h".
46 *
47 * Revision 1.5 2004/04/05 11:52:36 j_novak
48 * First operational (but not tested!) version of checks of Eintein equation.
49 *
50 * Revision 1.4 2004/04/01 16:09:02 j_novak
51 * Trace of K_ij is now member of Time_slice (it was member of Time_slice_conf).
52 * Added new methods for checking 3+1 Einstein equations (preliminary).
53 *
54 * Revision 1.3 2004/03/29 12:00:16 e_gourgoulhon
55 * Computation of extrinsic curvature now performed via new methods
56 * Vector::ope_killing.
57 *
58 * Revision 1.2 2004/03/28 21:29:45 e_gourgoulhon
59 * Evolution_std's renamed with suffix "_evol"
60 * Method gam() modified
61 * Added special constructor for derived classes.
62 *
63 * Revision 1.1 2004/03/26 13:33:02 j_novak
64 * New methods for accessing/updating members (nn(), beta(), gam_uu(), k_uu(), ...)
65 *
66 *
67 *
68 * $Header: /cvsroot/Lorene/C++/Source/Time_slice/time_slice_access.C,v 1.9 2014/10/13 08:53:47 j_novak Exp $
69 *
70 */
71
72// C headers
73#include <cassert>
74
75// Lorene headers
76#include "time_slice.h"
77
78namespace Lorene {
79const Scalar& Time_slice::nn() const {
80
82 return n_evol[jtime] ;
83
84}
85
86
87const Vector& Time_slice::beta() const {
88
90 return beta_evol[jtime] ;
91
92
93}
94
95const Metric& Time_slice::gam() const {
96
97 if (p_gamma == 0x0) {
98 gam_dd() ; // may force the computation of p_gamma
99 if (p_gamma == 0x0) p_gamma = new Metric( gam_dd() ) ;
100 }
101
102 return *p_gamma ;
103
104}
105
106
108
109 if (!( gam_dd_evol.is_known(jtime)) ) {
111 if (p_gamma == 0x0) {
112 p_gamma = new Metric( gam_uu_evol[jtime] ) ;
113 }
114
116 }
117
118 return gam_dd_evol[jtime] ;
119
120}
121
123
124 if (!( gam_uu_evol.is_known(jtime)) ) {
126 gam_uu_evol.update(gam().con(), jtime, the_time[jtime] ) ;
127 }
128
129 return gam_uu_evol[jtime] ;
130
131}
132
133
134
136
137 if ( ! (k_dd_evol.is_known(jtime)) ) {
138
139 Vector beta_d = beta().down(0, gam()) ;
140
141 gam_dd() ; // to make sure that gam_dd is up to date before taking its
142 // time derivative
143
144 Sym_tensor resu = beta_d.ope_killing(gam())
146
147 resu = resu / (2*nn()) ;
148
150
151 }
152
153 return k_dd_evol[jtime] ;
154
155}
156
158
159 if ( ! (k_uu_evol.is_known(jtime)) ) {
160
161 gam_uu() ; // to make sure that gam_uu is up to date before taking its
162 // time derivative
163
164 Sym_tensor resu = beta().ope_killing(gam())
166
167 resu = resu / (2*nn()) ;
168
170
171 }
172
173 return k_uu_evol[jtime] ;
174
175}
176
177const Scalar& Time_slice::trk() const {
178
179 if ( ! (trk_evol.is_known(jtime)) ) {
180
181 if ( k_uu_evol.is_known(jtime) )
182 trk_evol.update( k_uu().trace(gam()), jtime, the_time[jtime] ) ;
183 else
184 trk_evol.update( k_dd().trace(gam()), jtime, the_time[jtime] ) ;
185
186 }
187
188 return trk_evol[jtime] ;
189
190}
191
192
193
194
195
196
197
198
199
200
201}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
virtual void update(const TyT &new_value, int j, double time_j)
Sets a new value at a given time step.
TyT time_derive(int j, int n=2) const
Computes the time derivative at time step j by means of a n-th order scheme, from the values at steps...
Definition evolution.C:504
bool is_known(int j) const
Checks whether the value a given time step has been set.
Definition evolution.C:335
Metric for tensor calculation.
Definition metric.h:90
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition metric.C:280
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
virtual const Sym_tensor & k_dd() const
Extrinsic curvature tensor (covariant components ) at the current time step (jtime )
int jtime
Time step index of the latest slice.
Definition time_slice.h:190
Evolution_std< Scalar > trk_evol
Values at successive time steps of the trace K of the extrinsic curvature.
Definition time_slice.h:224
Metric * p_gamma
Pointer on the induced metric at the current time step (jtime)
Definition time_slice.h:239
Evolution_std< Sym_tensor > k_dd_evol
Values at successive time steps of the covariant components of the extrinsic curvature tensor .
Definition time_slice.h:208
virtual const Scalar & trk() const
Trace K of the extrinsic curvature at the current time step (jtime )
virtual const Scalar & nn() const
Lapse function N at the current time step (jtime )
virtual const Sym_tensor & gam_dd() const
Induced metric (covariant components ) at the current time step (jtime )
Evolution_std< Sym_tensor > gam_dd_evol
Values at successive time steps of the covariant components of the induced metric .
Definition time_slice.h:198
Evolution_std< Sym_tensor > gam_uu_evol
Values at successive time steps of the contravariant components of the induced metric .
Definition time_slice.h:203
int scheme_order
Order of the finite-differences scheme for the computation of time derivatives.
Definition time_slice.h:187
virtual const Vector & beta() const
shift vector at the current time step (jtime )
const Metric & gam() const
Induced metric at the current time step (jtime )
virtual const Sym_tensor & k_uu() const
Extrinsic curvature tensor (contravariant components ) at the current time step (jtime )
Evolution_std< Scalar > n_evol
Values at successive time steps of the lapse function N.
Definition time_slice.h:216
virtual const Sym_tensor & gam_uu() const
Induced metric (contravariant components ) at the current time step (jtime )
Evolution_std< double > the_time
Time label of each slice.
Definition time_slice.h:193
Evolution_std< Sym_tensor > k_uu_evol
Values at successive time steps of the contravariant components of the extrinsic curvature tensor .
Definition time_slice.h:213
Evolution_std< Vector > beta_evol
Values at successive time steps of the shift vector .
Definition time_slice.h:219
Tensor field of valence 1.
Definition vector.h:188
Lorene prototypes.
Definition app_hor.h:64