4#include "utilitaires.h"
9#include "param_elliptic.h"
33 const Scalar& psi4 = conf_fact*conf_fact*conf_fact*conf_fact;
47 const Metric tgam(tgam_uu);
54 for (
int ii=1; ii<=3; ii++)
55 for(
int jj=1; jj<=3; jj++)
64 for (
int ii=1; ii<=3; ii++)
80 for (
int ii=1; ii<=3; ii++)
81 for(
int jj=1; jj<=3; jj++)
124 for (
int i=1; i<=3; i++) {
125 for (
int j=1; j<=i; j++) {
127 for (
int k=1; k<=3; k++) {
128 for (
int l=1; l<=3; l++) {
129 tmp += dhh(i,k,l) * dhh(j,l,k) ;
132 sym_tmp.
set(i,j) = tmp ;
135 ricci_star -= sym_tmp ;
136 for (
int i=1; i<=3; i++) {
137 for (
int j=1; j<=i; j++) {
139 for (
int k=1; k<=3; k++) {
140 for (
int l=1; l<=3; l++) {
141 for (
int m=1; m<=3; m++) {
142 for (
int n=1; n<=3; n++) {
144 tmp += 0.5 * tgam_uu(i,k)* tgam_uu(j,l)
145 * dhh(m,n,k) * dtgam(m,n,l)
146 + tgam_dd(n,l) * dhh(m,n,k)
147 * (tgam_uu(i,k) * dhh(j,l,m) + tgam_uu(j,k) * dhh(i,l,m) )
148 - tgam_dd(k,l) *tgam_uu(m,n) * dhh(i,k,m) * dhh(j,l,n) ;
153 sym_tmp.
set(i,j) = tmp ;
156 ricci_star += sym_tmp ;
158 ricci_star = 0.5 * ricci_star ;
165 contract(dhh, 0, 1, dtgam, 0, 1), 0, 1 )
167 contract(dhh, 0, 1, dtgam, 0, 2), 0, 1 ) ;
174 sym_tmp =
lapse * (ricci_star + 8.* tdln_psi_u * tdln_psi_u)
175 + 4.*( tdln_psi_u * tdnn_u + tdnn_u * tdln_psi_u )
176 - 0.3333333333333333 *
177 (
lapse * (tricci_scal + 8.*
contract(dln_psi, 0, tdln_psi_u, 0) )
178 + 8.*
contract(dln_psi, 0, tdnn_u, 0) ) *tgam_uu ;
180 ss = sym_tmp / psi4 ;
189 for (
int i=1; i<=3; i++) {
190 for (
int j=1; j<=i; j++) {
192 for (
int k=1; k<=3; k++) {
193 for (
int l=1; l<=3; l++) {
194 tmp += ( hh(i,k)*dhh(l,j,k) + hh(k,j)*dhh(i,l,k)
195 - hh(k,l)*dhh(i,j,k) ) * dqq(l) ;
198 sym_tmp.
set(i,j) += 0.5 * tmp ;
207 sym_tmp -= 0.3333333333333333 * tmp *tgam_uu ;
209 ss -= sym_tmp / (psi4*conf_fact*conf_fact) ;
211 for (
int i=1; i<=3; i++) {
212 for (
int j=1; j<=i; j++) {
214 for (
int k=1; k<=3; k++) {
215 for (
int l=1; l<=3; l++) {
216 tmp += tgam_dd(k,l) * aa(i,k) * aa(j,l) ;
219 sym_tmp.
set(i,j) = tmp ;
223 tmp = psi4 * strain_tens.
trace(tgam) ;
225 ss += (2.*
lapse) * ( sym_tmp);
227 Sym_tensor ss2 =2.*
lapse*( qpig*(psi4*strain_tens - 0.33333333333333 * tmp * tgam_uu));
254 source_hh += 2.*
lapse * ss ;
258 Vector vtmp = beta_point ;
267 + 1.3333333333333333 * div_beta* (hh_point - lbh)
271 for (
int i=1; i<=3; i++) {
272 for (
int j=1; j<=i; j++) {
274 for (
int k=1; k<=3; k++) {
279 sym_tmp.
set(i,j) = tmp ;
283 source_hh -=
lapse / (psi4*conf_fact*conf_fact) * sym_tmp ;
290 source_hh += 0.6666666666666666*
291 ( tmp - 0.6666666666666666* div_beta * div_beta ) * hh ;
304 source_hh += 0.6666666666666666* div_beta * l_beta - sym_tmp ;
308 source_hh = -source_hh;
Spherical orthonormal vectorial bases (triads).
Scalar lapse
Lapse function.
Sym_tensor hij
Deviation tensor( non-flat part of the conformal 3-metric on the slice; see Bonazzola et al.
void secmembre_kerr(Sym_tensor &source_hh)
Computes the rhs of hyperbolic equation for conformal metric assuming statioarity; WARNING; up to now...
Vector shift
Shift vector.
Sym_tensor hatA
Rescaled tracefree extrinsic curvature tensor: see Cordero et al.
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.