LORENE
excision_hor.h
1/*
2 * Definition of Lorene class Excision_hor, friend class of Spheroid.
3 *
4 */
5
6/*
7 * Copyright (c) 2009 Jose-Luis Jaramillo & Nicolas Vasset
8 *
9 * This file is part of LORENE.
10 *
11 * LORENE is free software; you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License version 2
13 * as published by the Free Software Foundation.
14 *
15 * LORENE is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with LORENE; if not, write to the Free Software
22 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 *
24 */
25
26#ifndef __EXCISIONHOR_H_
27#define __EXCISIONHOR_H_
28
29/*
30 * $Header: /cvsroot/Lorene/C++/Include/excision_hor.h,v 1.3 2014/10/13 08:52:34 j_novak Exp $
31 *
32 */
33#include "metric.h"
34#include "spheroid.h"
35
36namespace Lorene {
44
45
46 // Data :
47 // -----
48 protected:
51
54
57
60
63
66
68 double delta_t;
69
72
75
76
77 // Derived data :
78 // ------------
79 protected:
81 mutable Scalar* p_get_BC_bmN ;
82 mutable Scalar* p_get_BC_bpN ;
84
85
86 // Constructors - Destructor
87 // -------------------------
88 public:
89
102 Excision_hor(const Scalar& h_in, const Metric& gij, const Sym_tensor& Kij2, const Scalar& ppsi, const Scalar& nn, const Vector& beta,const Sym_tensor& Tij2, double timestep, int int_nos = 1) ;
103 Excision_hor(const Excision_hor& ) ;
104
106 Excision_hor(FILE* ) ;
107
108 virtual ~Excision_hor() ;
109
110
111 // Memory management
112 // -----------------
113 protected:
115 virtual void del_deriv() const ;
116
118 void set_der_0x0() const ;
119
120
121 // Mutators / assignment
122 // ---------------------
123 public:
125 void operator=(const Excision_hor&) ;
126
127
128 // Accessors
129 // ---------
130 public:
132 const Spheroid& get_sph() const {return sph; };
133
135 const Scalar& get_conf_fact() const {return conf_fact; } ;
136
138 const Scalar& get_lapse() const {return lapse ; } ;
139
141 const Vector& get_shift() const {return shift ; } ;
142
144 const Metric& get_gamij() const {return gamij ; } ;
145
147 const Sym_tensor& get_Kij() const {return Kij ; } ;
148
150 double get_delta_t() const {return delta_t ;};
151
153 double get_no_of_steps() const {return no_of_steps ;};
154
156 const Sym_tensor& get_Tij() const{return Tij; } ;
157
160
162 Scalar& set_lapse() {del_deriv() ; return lapse ; } ;
163
165 Vector& set_shift() {del_deriv() ; return shift ; } ;
166
168 Metric& set_gamij() {del_deriv() ; return gamij ; } ;
169
171 Sym_tensor& set_Kij() {del_deriv() ; return Kij ; } ;
172
174 Sym_tensor& set_Tij() {del_deriv() ; return Tij; } ;
175
176
177 double& set_delta_t() {del_deriv() ; return delta_t ; } ;
178
179 double& set_no_of_steps() {del_deriv() ; return no_of_steps ; } ;
180
181
182 // Computational functions
183 // -----------------------
184 public:
185
187 const Scalar& get_BC_conf_fact() const ;
189 const Scalar& get_BC_bmN(int choice_bmN, double value = 1.) const ;
192 const Scalar& get_BC_bpN(int choice_bpN, double c_bpn_lap = 1., double c_bpn_fin = 1., Scalar *bpN_fin = 0x0) const ;
194 const Vector& get_BC_shift(double c_V_lap) const ;
195
196 // Outputs
197 // -------
198 public:
199 virtual void sauve(FILE *) const ;
200
202 friend ostream& operator<<(ostream& , const Spheroid& ) ;
203
204};
205
206ostream& operator<<(ostream& , const Spheroid& ) ;
207
208
209}
210#endif
Surface where boundary conditions for quantities in the bulk will be calculated It relies on geometri...
const Vector & get_BC_shift(double c_V_lap) const
Source of Dirichlet BC for the shift, issued from BC on bpN and a gauge condition on the tangential s...
Vector shift
The Shift 3-vector on the slice.
double get_no_of_steps() const
Returns the internal number of timesteps for one iteration.
Scalar & set_lapse()
Sets the lapse function.
Sym_tensor Tij
Value of the impulsion-energy tensor on the spheroid.
Sym_tensor & set_Tij()
Sets the value of the impulsion-energy tensor.
const Sym_tensor & get_Tij() const
Returns the value of the impulsion-energy tensor.
Scalar * p_get_BC_conf_fact
Source of Neumann BC on , derived from the vanishing expansion.
Excision_hor(FILE *)
Constructor from a file (see sauve(FILE*) )
double no_of_steps
The internal number of timesteps for one iteration.
Vector * p_get_BC_shift
Source of Dirichlet BC for the shift, issued from BC on bpN and a gauge condition on the tangential s...
Scalar * p_get_BC_bpN
Arbitrary source of Dirichlet BC for (b+N), case 0: based on a parabolic driver towards a constant va...
Metric & set_gamij()
Sets the 3d metric of the TimeSlice.
Scalar & set_conf_fact()
Sets the value of the conformal factor.
Metric gamij
The 3-d metric on the slice.
friend ostream & operator<<(ostream &, const Spheroid &)
Display.
Sym_tensor Kij
The 3-d extrinsic curvature on the slice.
virtual ~Excision_hor()
Destructor.
const Sym_tensor & get_Kij() const
returns the 3-d extrinsic curvature
void operator=(const Excision_hor &)
Assignment to another Excision_hor.
const Scalar & get_BC_bpN(int choice_bpN, double c_bpn_lap=1., double c_bpn_fin=1., Scalar *bpN_fin=0x0) const
Case 0: Arbitrary source of Dirichlet BC for (b+N), based on a parabolic driver towards a constant va...
Vector & set_shift()
Sets the shift vector field.
double get_delta_t() const
Returns the timestep used for evolution.
const Scalar & get_lapse() const
Returns the lapse function.
virtual void sauve(FILE *) const
Save in a file.
Spheroid sph
The associated Spheroid object.
Scalar * p_get_BC_bmN
Source of Dirichlet BC for (b-N).
const Scalar & get_conf_fact() const
Returns the conformal factor associated with the surface.
Scalar lapse
The lapse defined on the 3 slice.
Scalar conf_fact
The value of the conformal factor on the 3-slice.
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Sym_tensor & set_Kij()
Sets the extrinsic curvature.
const Scalar & get_BC_conf_fact() const
Source of Neumann BC on , derived from the vanishing expansion.
double delta_t
The time step for evolution in parabolic drivers.
const Metric & get_gamij() const
Returns the symmetric tensor .
const Spheroid & get_sph() const
Returns the spheroid.
const Scalar & get_BC_bmN(int choice_bmN, double value=1.) const
Source of Dirichlet BC for (b-N): case 0: based on an entropy prescription, case 1: from a component ...
virtual void del_deriv() const
Deletes all the derived quantities.
const Vector & get_shift() const
Returns the shift vector field.
Metric for tensor calculation.
Definition metric.h:90
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Spheroidal 2-surfaces embedded in a time-slice of the 3+1 formalism.
Definition spheroid.h:84
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
Tensor field of valence 1.
Definition vector.h:188
Lorene prototypes.
Definition app_hor.h:64