3#include "param_elliptic.h"
4#include "excised_slice.h"
27 const Scalar& psi4 = conf_fact*conf_fact*conf_fact*conf_fact;
41 const Metric tgam(tgam_uu);
48 for (
int ii=1; ii<=3; ii++)
49 for(
int jj=1; jj<=3; jj++)
58 for (
int ii=1; ii<=3; ii++)
74 for (
int ii=1; ii<=3; ii++)
75 for(
int jj=1; jj<=3; jj++)
118 for (
int i=1; i<=3; i++) {
119 for (
int j=1; j<=i; j++) {
121 for (
int k=1; k<=3; k++) {
122 for (
int l=1; l<=3; l++) {
123 tmp += dhh(i,k,l) * dhh(j,l,k) ;
126 sym_tmp.
set(i,j) = tmp ;
129 ricci_star -= sym_tmp ;
130 for (
int i=1; i<=3; i++) {
131 for (
int j=1; j<=i; j++) {
133 for (
int k=1; k<=3; k++) {
134 for (
int l=1; l<=3; l++) {
135 for (
int m=1; m<=3; m++) {
136 for (
int n=1; n<=3; n++) {
138 tmp += 0.5 * tgam_uu(i,k)* tgam_uu(j,l)
139 * dhh(m,n,k) * dtgam(m,n,l)
140 + tgam_dd(n,l) * dhh(m,n,k)
141 * (tgam_uu(i,k) * dhh(j,l,m) + tgam_uu(j,k) * dhh(i,l,m) )
142 - tgam_dd(k,l) *tgam_uu(m,n) * dhh(i,k,m) * dhh(j,l,n) ;
147 sym_tmp.
set(i,j) = tmp ;
150 ricci_star += sym_tmp ;
152 ricci_star = 0.5 * ricci_star ;
159 contract(dhh, 0, 1, dtgam, 0, 1), 0, 1 )
161 contract(dhh, 0, 1, dtgam, 0, 2), 0, 1 ) ;
168 sym_tmp =
lapse * (ricci_star + 8.* tdln_psi_u * tdln_psi_u)
169 + 4.*( tdln_psi_u * tdnn_u + tdnn_u * tdln_psi_u )
170 - 0.3333333333333333 *
171 (
lapse * (tricci_scal + 8.*
contract(dln_psi, 0, tdln_psi_u, 0) )
172 + 8.*
contract(dln_psi, 0, tdnn_u, 0) ) *tgam_uu ;
174 ss = sym_tmp / psi4 ;
183 for (
int i=1; i<=3; i++) {
184 for (
int j=1; j<=i; j++) {
186 for (
int k=1; k<=3; k++) {
187 for (
int l=1; l<=3; l++) {
188 tmp += ( hh(i,k)*dhh(l,j,k) + hh(k,j)*dhh(i,l,k)
189 - hh(k,l)*dhh(i,j,k) ) * dqq(l) ;
192 sym_tmp.
set(i,j) += 0.5 * tmp ;
201 sym_tmp -= 0.3333333333333333 * tmp *tgam_uu ;
203 ss -= sym_tmp / (psi4*conf_fact*conf_fact) ;
205 for (
int i=1; i<=3; i++) {
206 for (
int j=1; j<=i; j++) {
208 for (
int k=1; k<=3; k++) {
209 for (
int l=1; l<=3; l++) {
210 tmp += tgam_dd(k,l) * aa(i,k) * aa(j,l) ;
213 sym_tmp.
set(i,j) = tmp ;
217 tmp = psi4 * strain_tens.
trace(tgam) ;
219 ss += (2.*
lapse) * ( sym_tmp);
221 Sym_tensor ss2 =2.*
lapse*( qpig*(psi4*strain_tens - 0.33333333333333 * tmp * tgam_uu));
248 source_hh += 2.*
lapse * ss ;
252 Vector vtmp = beta_point ;
261 + 1.3333333333333333 * div_beta* (hh_point - lbh)
265 for (
int i=1; i<=3; i++) {
266 for (
int j=1; j<=i; j++) {
268 for (
int k=1; k<=3; k++) {
273 sym_tmp.
set(i,j) = tmp ;
277 source_hh -=
lapse / (psi4*conf_fact*conf_fact) * sym_tmp ;
284 source_hh += 0.6666666666666666*
285 ( tmp - 0.6666666666666666* div_beta * div_beta ) * hh ;
298 source_hh += 0.6666666666666666* div_beta * l_beta - sym_tmp ;
302 source_hh = -source_hh;
Spherical orthonormal vectorial bases (triads).
Sym_tensor hatA
Rescaled tracefree extrinsic curvature tensor: rescaled same way as Cordero et al.
void secmembre_kerr(Sym_tensor &source_hh)
If hor_type=1, computes the rhs of hyperbolic equation for conformal metric assuming stationarity; WA...
Sym_tensor hij
Deviation tensor( non-flat part of the conformal 3-metric on the slice; see Bonazzola et al.
Scalar lapse
Lapse function.
Vector shift
Shift vector.
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Flat metric for tensor calculation.
Metric for tensor calculation.
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Tensor field of valence 0 (or component of a tensorial field).
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the co...
void annule_hard()
Sets the Scalar to zero in a hard way.
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Scalar derive_lie(const Vector &v) const
Computes the derivative of this along a vector field v.
Class intended to describe valence-2 symmetric tensors.
Sym_tensor derive_lie(const Vector &v) const
Computes the Lie derivative of this with respect to some vector field v.
Symmetric tensors (with respect to two of their arguments).
Tensor field of valence 1.
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
const Scalar & divergence(const Metric &) const
The divergence of this with respect to a Metric .
Sym_tensor ope_killing_conf(const Metric &gam) const
Computes the conformal Killing operator associated with a given metric.
Scalar & set(int)
Read/write access to a component.
Cmp sqrt(const Cmp &)
Square root.
Cmp log(const Cmp &)
Neperian logarithm.
const Map & get_mp() const
Returns the mapping.
const Tensor_sym & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
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...
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
const Tensor & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Standard units of space, time and mass.