99 jac2d(map, 2, COV, map.get_bvect_spher()),
100 proj(map, 2, COV, map.get_bvect_spher()),
101 qq(map, COV, map.get_bvect_spher()),
102 ss (map, COV, map.get_bvect_spher()),
103 ephi(map, CON, map.get_bvect_spher()),
104 qab(map.flat_met_spher()),
106 hh(map, COV, map.get_bvect_spher()),
108 ll(map, COV, map.get_bvect_spher()),
109 jj(map, COV, map.get_bvect_spher()),
122 assert(radius > 1.e-15) ;
140 for (
int i=1; i<=3; i++)
141 hh.
set(i,i) = 2./radius ;
155 jac2d(h_in.get_mp(),2, COV, h_in.get_mp().get_bvect_spher()),
156 proj(h_in.get_mp(),2, COV, h_in.get_mp().get_bvect_spher()),
157 qq(h_in.get_mp(), COV, h_in.get_mp().get_bvect_spher()),
158 ss (h_in.get_mp(), COV, h_in.get_mp().get_bvect_spher()),
159 ephi(h_in.get_mp(), CON, h_in.get_mp().get_bvect_spher()),
160 qab(h_in.get_mp().flat_met_spher()),
161 ricci(h_in.get_mp()),
162 hh(h_in.get_mp(), COV, h_in.get_mp().get_bvect_spher()),
164 ll(h_in.get_mp(), COV, h_in.get_mp().get_bvect_spher()),
165 jj(h_in.get_mp(), COV, h_in.get_mp().get_bvect_spher()),
178 assert(mp_rad != 0x0) ;
183 int np = gri2d.
get_np(0) ;
184 int nt = gri2d.
get_nt(0) ;
186 int nr3 = gri3d.
get_nr(0) ;
187 int nt3 = gri3d.
get_nt(0) ;
188 int np3 = gri3d.
get_np(0) ;
190 assert( nt == nt3 ) ;
191 assert ( np == np3 );
211 for (
int f= 0; f<nz; f++)
212 for (
int k=0; k<np3; k++)
213 for (
int j=0; j<nt3; j++) {
214 for (
int l=0; l<nr3; l++) {
233 jac.std_spectral_base();
237 jac.set(1,2) = -h_surf3.
srdsdt() ;
244 for (
int l=1; l<4; l++)
245 for (
int m=1; m<4; m++)
246 for (
int k=0; k<np; k++)
247 for (
int j=0; j<nt; j++) {
249 j, k, pipo, lz, xi) ;
251 jac(l,m).get_spectral_va().val_point_jk(lz, xi, j, k) ;
256 jac_inv.
set(1,2) = - jac_inv(1,2);
257 jac_inv.
set(1,3) = - jac_inv(1,3) ;
267 ss3d = ss3d /
sqrt (ssnorm) ;
268 ss3dcon = ss3dcon /
sqrt (ssnorm) ;
283 for (
int ii=1; ii <=3; ii++){
287 for (
int ii=1; ii <=3; ii++)
288 for (
int jjy = 1; jjy <=3; jjy ++){
289 jac_inv.
set(ii, jjy).
dec_dzpuis(jac_inv(ii, jjy).get_dzpuis());
290 jac.set(ii, jjy).dec_dzpuis(jac(ii, jjy).get_dzpuis());
298 ss3 =
contract(jac_inv, 0, ss3d, 0);
302 for (
int l=1; l<4; l++)
303 for (
int k=0; k<np; k++)
304 for (
int j=0; j<nt; j++) {
308 ss3(l).get_spectral_va().val_point_jk(lz, xi, j, k) ;
315 ephi3.set(2).annule_hard();
317 ephi33.mult_r(); ephi33.mult_sint();
318 ephi3.set(3) = ephi33;
319 ephi3.std_spectral_base();
324 for (
int l=1; l<4; l++)
325 for (
int k=0; k<np; k++)
326 for (
int j=0; j<nt; j++) {
328 j, k, pipo, lz, xi) ;
330 ephi3(l).get_spectral_va().val_point_jk(lz, xi, j, k) ;
335 Tensor proj_prov = gamij.
con().
down(1, gamij) - ss3dcon*ss3d;
338 proj3d.allocate_all();
340 proj3d.std_spectral_base();
342 for (
int l=1; l<4; l++)
343 for (
int m=1; m<4; m++)
344 for (
int k=0; k<np; k++)
345 for (
int j=0; j<nt; j++) {
347 j, k, pipo, lz, xi) ;
349 proj3d(l,m).get_spectral_va().val_point_jk(lz, xi, j, k) ;
362 qq3d.set(2,2)= gamij.
cov()(1,1) * (h_surf3.
srdsdt())* (h_surf3.
srdsdt())
363 + 2* gamij.
cov()(1,2)* (h_surf3.
srdsdt()) + gamij.
cov()(2,2);
373 qab2.allocate_all() ;
374 qab2.std_spectral_base();
375 for (
int l=1; l<4; l++)
376 for (
int m=1; m<4; m++)
377 for (
int k=0; k<np; k++)
378 for (
int j=0; j<nt; j++) {
380 j, k, pipo, lz, xi) ;
381 qab2.set(l,m).set_grid_point(0, k, j, 0) =
382 qq3d(l,m).get_spectral_va().val_point_jk(lz, xi, j, k) ;
388 - ss3d * ss3d) , 0), 1) ;
391 for (
int l=1; l<4; l++)
392 for (
int m=1; m<4; m++)
393 for (
int k=0; k<np; k++)
394 for (
int j=0; j<nt; j++) {
396 j, k, pipo, lz, xi) ;
398 qqq(l,m).get_spectral_va().val_point_jk(lz, xi, j, k) ;
404 for (
int k=0; k<np; k++)
405 for (
int j=0; j<nt; j++) {
407 j, k, pipo, lz, xi) ;
417 for (
int k=0; k<np; k++)
418 for (
int j=0; j<nt; j++) {
430 for (
int k=0; k<np; k++)
431 for (
int j=0; j<nt; j++) {
440 double rayon =
sqrt(
area()/(4.*M_PI));
441 ftilde = -ftilde/(rayon*rayon);
449 double *a_tilde =
new double[nombre];
453 for (
int k=0; k<np+1; k++)
454 for (
int j=0; j<nt; j++) {
455 int l_q, m_q, base_r ;
457 if (nullite_plm(j, nt, k, np, base) == 1) {
458 donne_lm(1, lz, j, k, base, m_q, l_q,base_r) ;
460 a_tilde[l_q] = (*coefftilde)(0, k, j, 0);
470 for (
int k=0; k<np+1; k++)
471 for (
int j=0; j<nt; j++) {
472 int l_q2, m_q2, base_r2 ;
474 if (nullite_plm(j, nt, k, np, base2) == 1) {
475 donne_lm(1, lz, j, k, base2, m_q2, l_q2,base_r2) ;
478 (*dzeta).set(0,k,j,0) = a_tilde[0] + (a_tilde[1]/3.)*
sqrt(2.*1. + 1.);
481 (*dzeta).set(0,k,j,0) =
482 (a_tilde[l_q2 +1]/(2.*l_q2 + 3.))*
sqrt((2.*(l_q2 +1.)+1.)/(2.*l_q2 + 1.))
483 - (a_tilde[l_q2 -1]/(2.*l_q2 - 1.))
484 *
sqrt((2.*(l_q2 - 1.)+1.)/(2.*l_q2 + 1.));
501 for (
int l=1; l<4; l++)
502 for (
int k=0; k<np; k++)
503 for (
int j=0; j<nt; j++) {
507 ll3(l).get_spectral_va().val_point_jk(lz, xi, j, k) ;
516 Tensor jj3d = Kij - ss3d*ll3d - ll3d*ss3d
517 -
contract(Kij, 0 , 1, sxss3dcon , 0, 1)* sxss3d ;
521 for (
int l=1; l<4; l++)
522 for (
int m=1; m<4; m++)
523 for (
int k=0; k<np; k++)
524 for (
int j=0; j<nt; j++) {
526 j, k, pipo, lz, xi) ;
528 jj3(l,m).get_spectral_va().val_point_jk(lz, xi, j, k) ;
534 contract(proj_prov, 0,ss3d.derive_cov(gamij),1), 1 ) ;
538 for (
int l=1; l<4; l++)
539 for (
int m=1; m<4; m++)
540 for (
int k=0; k<np; k++)
541 for (
int j=0; j<nt; j++) {
543 j, k, pipo, lz, xi) ;
545 hh3(l,m).get_spectral_va().val_point_jk(lz, xi, j, k) ;
550 Tensor hh3dupdown = hh3d.
up(0, gamij);
556 Tensor hh3dupup = hh3dupdown.
up(1,gamij);
561 Scalar ricci22 = ricciscal3
567 ricci22 += (hh3d.
trace(gamij)*hh3d.
trace(gamij))
575 for (
int k=0; k<np; k++)
576 for (
int j=0; j<nt; j++) {
609 issphere(sph_in.issphere)
660 if (p_epsilon_A_minus_one != 0x0)
delete p_epsilon_A_minus_one ;
665 if (p_delta != 0x0)
delete p_delta ;
676 p_epsilon_A_minus_one = 0x0;
741 double rayon =
sqrt(
area()/(4.*M_PI));
753 double rayon =
sqrt(
area()/(4.*M_PI));
755 double factor =
mass()/(8. * M_PI);
757 {
for (
int compte=0; compte <=order -1; compte++)
758 factor = factor*rayon;
768 for (
int nn=1; nn<order; nn++){
770 Scalar Pnnew = (2.*nn +1.)*Pn;
772 Pnnew = Pnnew - nn*Pnold;
773 Pnnew = Pnnew/(double(nn) + 1.);
810 double rayon =
sqrt(
area()/(4.*M_PI));
812 double factor = 1./(8. * M_PI);
814 {
for (
int compte=0; compte <=order -2; compte++)
815 factor = factor*rayon;
828 for (
int nn=1; nn<order; nn++){
830 Scalar Pnnew = (2.*nn +1.)*Pn;
832 Pnnew = Pnnew - nn*Pnold;
833 Pnnew = Pnnew/(double(nn) + 1.);
842 dPn = Pn* zeta; dPn = dPn - Pnold; dPn = double(order)*dPn;
845 quotient = quotient*zeta*zeta; quotient = quotient -1.;
868 if (p_epsilon_A_minus_one == 0x0) {
873 return *p_epsilon_A_minus_one;
963 int valence1 = valence0 + 1 ;
964 int valence1m1 = valence1 - 1 ;
985 Itbl tipe(valence1) ;
987 for (
int id = 0;
id<valence0;
id++) {
988 tipe.
set(
id) = tipeuu(
id) ;
990 tipe.
set(valence1m1) = COV ;
1007 Itbl ind1(valence1) ;
1008 Itbl ind0(valence0) ;
1009 Itbl ind(valence0) ;
1029 for (
int ic=0; ic<ncomp1; ic++) {
1038 for (
int id = 0;
id < valence0;
id++) {
1039 ind0.
set(
id) = ind1(
id) ;
1043 int k = ind1(valence1m1) ;
1059 cresu = (uu(ind0)).srdsdt() ;
1062 for (
int id=0;
id<valence0;
id++) {
1064 switch ( ind0(
id) ) {
1096 cerr <<
"Connection_fspher::p_derive_cov : index problem ! " 1110 cresu = (uu(ind0)).srstdsdp() ;
1113 for (
int id=0;
id<valence0;
id++) {
1115 switch ( ind0(
id) ) {
1133 tmp.div_r_dzpuis(dz_resu) ;
1153 tmp.div_r_dzpuis(dz_resu) ;
1162 cerr <<
"Connection_fspher::p_derive_cov : index problem ! \n" 1174 cerr <<
"Connection_fspher::p_derive_cov : index problem ! \n" ;
1191 cout <<
"c'est pas fait!" << endl ;
1208 if (p_delta == 0x0) {
1212 christoflat.set_etat_zero() ;
1217 for (
int k=1; k<=3; k++) {
1218 for (
int i=1; i<=3; i++) {
1219 for (
int j=1; j<=i; j++) {
1221 for (
int l=1; l<=3; l++) {
1222 cc += qab.
con()(k,l) * (
1223 dgam(l,j,i) + dgam(i,l,j) - dgam(i,j,l) ) ;
1231 p_delta =
new Tensor (christoflat) ;
1248 for (
int y=1; y<=nbboucle; y++){
Sym_tensor * p_shear
The shear tensor.
Scalar h_surf
The location of the 2-surface as r = h_surf .
const Base_val & get_spectral_base() const
Returns the spectral bases of the Valeur va.
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
const Scalar & theta_plus() const
Computes the outgoing null expansion .
Metric for tensor calculation.
double * p_epsilon_P_minus_one
Refined Penrose parameter, difference wrt one.
double multipole_mass(const int order) const
Computes the mass multipole of a given order for the spheroid, assumed to be spherical.
double * p_area
The area of the 2-surface.
double angu_mom() const
Computes the angular momentum with respect to a divergence-free vector field tangent to the 2-surface...
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Tensor jac2d
The jacobian of the adaptation of coordinates (contravariant/covariant representation) ...
const Scalar & sqrt_q() const
Computes the normal vector field to the 2-surface.
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
int & set(int i)
Read/write of a particular element (index i ) (1D case)
double integrale_surface(const Cmp &ci, double rayon) const
Performs the surface integration of ci on the sphere of radius rayon .
virtual void del_deriv() const
Deletes all the derived quantities.
double * p_angu_mom
The angular momentum.
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
void coef() const
Computes the coeffcients of *this.
Cmp sqrt(const Cmp &)
Square root.
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
double multipole_angu(const int order) const
Computes the angular multipole of a given order for the spheroid, assumed to be spherical.
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
int sym_index1() const
Number of the first symmetric index (0<= id_sym1 < valence )
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Tensor up(int ind, const Metric &gam) const
Computes a new tensor by raising an index of *this.
void ylm()
Computes the coefficients of *this.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Tensor field of valence 0 (or component of a tensorial field).
Base class for coordinate mappings.
Vector ephi
The conformal Killing vector field on the 2-surface (set to by default to the axial vector associated...
int get_n_comp() const
Returns the number of stored components.
int sym_index2() const
Number of the second symmetric index (id_sym1 < id_sym2 < valence )
Basic integer array class.
const Sym_tensor & shear() const
Computes the shear of the 2-surface .
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
const Sym_tensor & ricci() const
Returns the Ricci tensor (given by the Connection p_connect )
virtual void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Tensor field of valence 1.
Scalar * p_theta_minus
Null ingoing expansion.
const Map & get_mp() const
Returns the mapping.
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
void annule_hard()
Sets the Scalar to zero in a hard way.
Spheroidal 2-surfaces embedded in a time-slice of the 3+1 formalism.
Spheroid(const Map_af &map, double radius)
The delta tensorial fields linked to Christoffel symbols.
const Vector & get_ephi() const
Returns the conformal Killing symmetry vector on the 2-surface.
Tensor derive_cov2dflat(const Tensor &uu) const
Computes the round covariant derivative on the spheroid.
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Tensor derive_cov2d(const Tensor &uu) const
Computes the total covariant derivative on the spheroid.
int get_index_type(int i) const
Gives the type (covariant or contravariant) of the index number i .
virtual void allocate_all()
Performs the memory allocation of all the elements, down to the double arrays of the Tbl s...
bool issphere
Flag to know whether the horizon is geometrically round or distorted.
Scalar trk
Trace K of the extrinsic curvature of the 3-slice.
Base class for pure radial mappings.
Vector ll
Normal-tangent component of the extrinsic curvature of the 3-slice.
double * p_multipole_angu
Angular momentum multipole for the spheroid.
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
int get_nzone() const
Returns the number of domains.
void operator=(const Spheroid &)
Assignment to another Spheroid.
virtual void sauve(FILE *) const
Save in a file.
double epsilon_P_minus_one() const
Computation of the classical Penrose parameter, and its difference wrt one.
double val_point_jk(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
void mult_cost()
Multiplication by .
const Scalar & get_ricci() const
Returns the 2-ricci scalar .
double area() const
Computes the area of the 2-surface.
Cmp pow(const Cmp &, int)
Power .
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Sym_tensor jj
Tangent components of the extrinsic curvature of the 3-slice.
int get_valence() const
Returns the valence.
int & set_index_type(int i)
Sets the type of the index number i .
virtual void val_lx_jk(double rr, int j, int k, const Param &par, int &l, double &xi) const =0
Computes the domain index l and the value of corresponding to a point of arbitrary r but collocation...
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Scalar ggg
Normalization function for the ingoing null vector k.
Bases of the spectral expansions.
double & set_grid_point(int l, int k, int j, int i)
Setting the value of the field at a given grid point.
const Scalar & srdsdt() const
Returns of *this .
const Scalar & theta_minus() const
Computes the ingoing null expansion .
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Sym_tensor hh
The ricci scalar on the 2-surface.
const Sym_tensor & get_qq() const
returns the 3-d degenerate 2-metric
Scalar fff
Normalization function for the outgoing null vector l.
Coefficients storage for the multi-domain spectral method.
double * p_mass
Mass defined from angular momentum.
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Valeur & set_spectral_va()
Returns va (read/write version)
Scalar & set(int)
Read/write access to a component.
Symmetric tensors (with respect to two of their arguments).
Vector ss
The adapted normal vector field to spheroid in the 3-slice.
double * p_multipole_mass
Mass multipole for the spheroid.
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Tensor proj
The 3-d projector on the 2-surface (contravariant-covariant form).
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
const Map & get_mp() const
Returns the mapping.
const Mg3d * get_angu_1dom() const
Returns the pointer on the associated mono-domain angular grid.
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Scalar * p_theta_plus
Classical Penrose parameter, difference wrt one.
Scalar ricci
Induced metric on the 2-surface .
Class intended to describe valence-2 symmetric tensors.
double epsilon_A_minus_one() const
Computation of the refined Penrose parameter for axisymmetric spacetimes, and its difference wrt one...
Scalar * p_sqrt_q
Surface element .
const Valeur & get_spectral_va() const
Returns va (read only version)
Sym_tensor qq
The 3-d covariant degenerated 2-metric on the surface.
virtual Itbl indices(int pos) const
Returns the indices of a component given by its position in the array cmp .
const Scalar & srstdsdp() const
Returns of *this .
double mass() const
Computes the mass as defined from the calculus of angular momentum, done with respect to a divergence...
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.
Coord r
r coordinate centered on the grid
virtual ~Spheroid()
Destructor.
const Tensor & delta() const
Computes the delta coefficients for covariant derivative.
Tensor down(int ind, const Metric &gam) const
Computes a new tensor by lowering an index of *this.