LORENE
|
Surface where boundary conditions for quantities in the bulk will be calculated It relies on geometrical properties of the associated Spheroid() (*** WARNING! under development***) More...
#include <excision_hor.h>
Public Member Functions | |
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) | |
Constructor of an excision surface embedded in a 3-slice (Time_slice ) of 3+1 formalism. More... | |
Excision_hor (const Excision_hor &) | |
Copy constructor. More... | |
Excision_hor (FILE *) | |
Constructor from a file (see sauve(FILE*) ) More... | |
virtual | ~Excision_hor () |
Destructor. More... | |
void | operator= (const Excision_hor &) |
Assignment to another Excision_hor. More... | |
const Spheroid & | get_sph () const |
Returns the spheroid. More... | |
const Scalar & | get_conf_fact () const |
Returns the conformal factor associated with the surface. More... | |
const Scalar & | get_lapse () const |
Returns the lapse function. More... | |
const Vector & | get_shift () const |
Returns the shift vector field. More... | |
const Metric & | get_gamij () const |
Returns the symmetric tensor . More... | |
const Sym_tensor & | get_Kij () const |
returns the 3-d extrinsic curvature More... | |
double | get_delta_t () const |
Returns the timestep used for evolution. More... | |
double | get_no_of_steps () const |
Returns the internal number of timesteps for one iteration. More... | |
const Sym_tensor & | get_Tij () const |
Returns the value of the impulsion-energy tensor. More... | |
Scalar & | set_conf_fact () |
Sets the value of the conformal factor. More... | |
Scalar & | set_lapse () |
Sets the lapse function. More... | |
Vector & | set_shift () |
Sets the shift vector field. More... | |
Metric & | set_gamij () |
Sets the 3d metric of the TimeSlice. More... | |
Sym_tensor & | set_Kij () |
Sets the extrinsic curvature. More... | |
Sym_tensor & | set_Tij () |
Sets the value of the impulsion-energy tensor. More... | |
double & | set_delta_t () |
double & | set_no_of_steps () |
const Scalar & | get_BC_conf_fact () const |
Source of Neumann BC on , derived from the vanishing expansion. More... | |
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 of projected Einstein Equations. More... | |
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 value. More... | |
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 shift (based on a parabolic driver). More... | |
virtual void | sauve (FILE *) const |
Save in a file. More... | |
Protected Member Functions | |
virtual void | del_deriv () const |
Deletes all the derived quantities. More... | |
void | set_der_0x0 () const |
Sets to 0x0 all the pointers on derived quantities. More... | |
Protected Attributes | |
Spheroid | sph |
The associated Spheroid object. More... | |
Scalar | conf_fact |
The value of the conformal factor on the 3-slice. More... | |
Scalar | lapse |
The lapse defined on the 3 slice. More... | |
Vector | shift |
The Shift 3-vector on the slice. More... | |
Metric | gamij |
The 3-d metric on the slice. More... | |
Sym_tensor | Kij |
The 3-d extrinsic curvature on the slice. More... | |
double | delta_t |
The time step for evolution in parabolic drivers. More... | |
double | no_of_steps |
The internal number of timesteps for one iteration. More... | |
Sym_tensor | Tij |
Value of the impulsion-energy tensor on the spheroid. More... | |
Scalar * | p_get_BC_conf_fact |
Source of Neumann BC on , derived from the vanishing expansion. More... | |
Scalar * | p_get_BC_bmN |
Source of Dirichlet BC for (b-N). More... | |
Scalar * | p_get_BC_bpN |
Arbitrary source of Dirichlet BC for (b+N), case 0: based on a parabolic driver towards a constant value, case 1:from a component of projected Einstein Equations. . More... | |
Vector * | p_get_BC_shift |
Source of Dirichlet BC for the shift, issued from BC on bpN and a gauge condition on the tangential shift. More... | |
Friends | |
ostream & | operator<< (ostream &, const Spheroid &) |
Display. More... | |
Surface where boundary conditions for quantities in the bulk will be calculated It relies on geometrical properties of the associated Spheroid() (*** WARNING! under development***)
Definition at line 43 of file excision_hor.h.
Lorene::Excision_hor::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 |
||
) |
Constructor of an excision surface embedded in a 3-slice (Time_slice
) of 3+1 formalism.
This is done from the Time_slice
data.
h_in | : the location of the surface r = h_in (WARNING:must be defined on a mono-domain angular grid) |
gij | : the 3-metric on the 3-slice |
Kij | : the extrinsic curvature of the 3-slice (covariant representation) |
timestep | : time interval associated with the parabolic-driven boundary conditions. |
int_nos | : Number of iterations to be done during timestep. |
Tij | : Value of the impulsion-energy tensor on the spheroid |
Definition at line 46 of file excision_hor.C.
References set_der_0x0().
Lorene::Excision_hor::Excision_hor | ( | const Excision_hor & | exc_in | ) |
Lorene::Excision_hor::Excision_hor | ( | FILE * | ) |
Constructor from a file (see sauve(FILE*)
)
|
virtual |
|
protectedvirtual |
Deletes all the derived quantities.
Definition at line 94 of file excision_hor.C.
References p_get_BC_bmN, p_get_BC_bpN, p_get_BC_conf_fact, p_get_BC_shift, and set_der_0x0().
const Scalar & Lorene::Excision_hor::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 of projected Einstein Equations.
Definition at line 147 of file excision_hor.C.
References Lorene::Scalar::allocate_all(), Lorene::Metric::con(), Lorene::Metric_flat::con(), Lorene::contract(), Lorene::Spheroid::delta(), Lorene::Spheroid::derive_cov2d(), Lorene::Spheroid::derive_cov2dflat(), Lorene::Map::flat_met_spher(), gamij, Lorene::Spheroid::get_hsurf(), Lorene::Spheroid::get_ll(), Lorene::Map::get_mg(), Lorene::Tensor::get_mp(), Lorene::Spheroid::get_qab(), Lorene::Spheroid::get_ricci(), lapse, Lorene::max(), Lorene::maxabs(), p_get_BC_bmN, Lorene::Scalar::poisson_angu(), Lorene::Metric::radial_vect(), Lorene::Tensor::set(), Lorene::Scalar::set_spectral_va(), Lorene::Spheroid::shear(), shift, sph, Lorene::Spheroid::sqrt_q(), Lorene::Tensor::std_spectral_base(), Lorene::Scalar::std_spectral_base(), Lorene::Spheroid::theta_minus(), Tij, Lorene::Tensor::trace(), Lorene::Tensor::up_down(), Lorene::Scalar::val_grid_point(), and Lorene::Valeur::ylm().
const Scalar & Lorene::Excision_hor::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 value.
Case 1: Source of Dirichlet BC for (b+N), from a component of projected Einstein Equations.
Definition at line 306 of file excision_hor.C.
References Lorene::Scalar::allocate_all(), Lorene::contract(), delta_t, Lorene::Spheroid::derive_cov2d(), gamij, get_BC_bmN(), Lorene::Spheroid::get_hsurf(), Lorene::Spheroid::get_ll(), Lorene::Map::get_mg(), Lorene::Tensor::get_mp(), Lorene::Spheroid::get_qab(), Lorene::Spheroid::get_ricci(), Lorene::Scalar::lapang(), lapse, p_get_BC_bpN, Lorene::Metric::radial_vect(), Lorene::Scalar::set_spectral_va(), Lorene::Spheroid::shear(), shift, sph, Lorene::Scalar::std_spectral_base(), Tij, Lorene::Tensor::trace(), Lorene::Tensor::up(), Lorene::Tensor::up_down(), Lorene::Scalar::val_grid_point(), and Lorene::Valeur::ylm().
const Scalar & Lorene::Excision_hor::get_BC_conf_fact | ( | ) | const |
Source of Neumann BC on , derived from the vanishing expansion.
Definition at line 120 of file excision_hor.C.
References conf_fact, Lorene::contract(), Lorene::Metric::cov(), Lorene::Scalar::derive_cov(), Lorene::Vector::divergence(), Lorene::Scalar::dsdr(), gamij, Kij, p_get_BC_conf_fact, Lorene::pow(), Lorene::Metric::radial_vect(), Lorene::Scalar::set_spectral_va(), Lorene::Tensor::std_spectral_base(), Lorene::Scalar::std_spectral_base(), and Lorene::Valeur::ylm().
const Vector & Lorene::Excision_hor::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 shift (based on a parabolic driver).
Definition at line 429 of file excision_hor.C.
References Lorene::Scalar::allocate_all(), Lorene::Tensor::annule_domain(), Lorene::Metric::con(), Lorene::contract(), Lorene::Metric::cov(), delta_t, Lorene::Tensor::derive_con(), Lorene::Tensor::derive_cov(), Lorene::Tensor_sym::derive_cov(), Lorene::Tensor::down(), gamij, Lorene::Map::get_mg(), Lorene::Tensor::get_mp(), Lorene::Spheroid::get_ricci(), lapse, p_get_BC_bmN, p_get_BC_bpN, p_get_BC_shift, Lorene::Metric::radial_vect(), shift, sph, Lorene::Vector::std_spectral_base(), Lorene::Tensor::std_spectral_base(), Lorene::Tensor::up_down(), and Lorene::Scalar::val_grid_point().
|
inline |
Returns the conformal factor associated with the surface.
Definition at line 135 of file excision_hor.h.
References conf_fact.
|
inline |
Returns the timestep used for evolution.
Definition at line 150 of file excision_hor.h.
References delta_t.
|
inline |
|
inline |
|
inline |
|
inline |
Returns the internal number of timesteps for one iteration.
Definition at line 153 of file excision_hor.h.
References no_of_steps.
|
inline |
|
inline |
|
inline |
Returns the value of the impulsion-energy tensor.
Definition at line 156 of file excision_hor.h.
References Tij.
void Lorene::Excision_hor::operator= | ( | const Excision_hor & | ) |
Assignment to another Excision_hor.
|
virtual |
Save in a file.
Definition at line 540 of file excision_hor.C.
|
inline |
Sets the value of the conformal factor.
Definition at line 159 of file excision_hor.h.
References conf_fact, and del_deriv().
|
protected |
Sets to 0x0
all the pointers on derived quantities.
Definition at line 103 of file excision_hor.C.
References p_get_BC_bmN, p_get_BC_bpN, p_get_BC_conf_fact, and p_get_BC_shift.
|
inline |
Sets the 3d metric of the TimeSlice.
Definition at line 168 of file excision_hor.h.
References del_deriv(), and gamij.
|
inline |
Sets the extrinsic curvature.
Definition at line 171 of file excision_hor.h.
References del_deriv(), and Kij.
|
inline |
Sets the lapse function.
Definition at line 162 of file excision_hor.h.
References del_deriv(), and lapse.
|
inline |
Sets the shift vector field.
Definition at line 165 of file excision_hor.h.
References del_deriv(), and shift.
|
inline |
Sets the value of the impulsion-energy tensor.
Definition at line 174 of file excision_hor.h.
References del_deriv(), and Tij.
|
friend |
Display.
|
protected |
The value of the conformal factor on the 3-slice.
Definition at line 53 of file excision_hor.h.
|
protected |
The time step for evolution in parabolic drivers.
Definition at line 68 of file excision_hor.h.
|
protected |
The 3-d metric on the slice.
Definition at line 62 of file excision_hor.h.
|
protected |
The 3-d extrinsic curvature on the slice.
Definition at line 65 of file excision_hor.h.
|
protected |
The lapse defined on the 3 slice.
Definition at line 56 of file excision_hor.h.
|
protected |
The internal number of timesteps for one iteration.
Definition at line 71 of file excision_hor.h.
|
mutableprotected |
Source of Dirichlet BC for (b-N).
Definition at line 81 of file excision_hor.h.
|
mutableprotected |
Arbitrary source of Dirichlet BC for (b+N), case 0: based on a parabolic driver towards a constant value, case 1:from a component of projected Einstein Equations. .
Definition at line 82 of file excision_hor.h.
|
mutableprotected |
Source of Neumann BC on , derived from the vanishing expansion.
Definition at line 80 of file excision_hor.h.
|
mutableprotected |
Source of Dirichlet BC for the shift, issued from BC on bpN and a gauge condition on the tangential shift.
Definition at line 83 of file excision_hor.h.
|
protected |
The Shift 3-vector on the slice.
Definition at line 59 of file excision_hor.h.
|
protected |
The associated Spheroid object.
Definition at line 50 of file excision_hor.h.
|
protected |
Value of the impulsion-energy tensor on the spheroid.
Definition at line 74 of file excision_hor.h.