24 char binhor_hh_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Bin_hor/binhor_hh.C,v 1.5 2014/10/13 08:52:42 j_novak Exp $" ;
57 #include "graphique.h" 58 #include "utilitaires.h" 95 rr1.std_spectral_base() ;
112 for (
int ind=1; ind<=3; ind++){
130 rr2_2.
set(1) = xx_2 ;
131 rr2_2.
set(2) = yy_2 ;
132 rr2_2.
set(3) = zz_2 ;
133 rr2_2.std_spectral_base() ;
150 for (
int ind=1; ind<=3; ind++){
151 nn2_2.
set(ind).
set_grid_point(nz2-1,k,j,i) = nn2_2(ind).val_grid_point(1,k,j,0) ;
192 for (
int i=1 ; i<=3 ; i++){
193 nn2.set(i).import(nn2_2(i)) ;
194 nn2.set(i).set_spectral_va().set_base(nn2_2(i).get_spectral_va().get_base()) ;
212 Scalar r2sr1 (1./unsr2*unsr1) ;
244 rr12.std_spectral_base() ;
248 r12 =
sqrt( rr12(1)*rr12(1) + rr12(2)*rr12(2) + rr12(3)*rr12(3)) ;
278 double sigma = 1.*r0 ;
284 fexp =
exp( -(rr - r0)*(rr - r0)/sigma/sigma ) ;
285 for (
int ii=0; ii<nz1-1; ii++)
301 f_delta = -5.*r1/(8.*r12*r12*r12) - 15./(8.*r1*r12) +
302 5.*r1*r1*unsr2/(8.*r12*r12*r12) + 1./(r1+1./unsr2+r12)/(r1+1./unsr2+r12)*
303 (1 + r1/r12 + r12/r1 - r1sr2 - r1*r1sr2/r12 + r12*r12*unsr2/(2*r1)) +
304 1./(r1+1./unsr2+r12)*(-7./r1 + 2./r12) ;
308 f_delta_zec = - 15./(8.*r1*r12) + 1./(r1+1./unsr2+r12)/(r1+1./unsr2+r12)*
309 (1 + r1/r12 + r12/r1 - r1sr2 - r1*r1sr2/r12 + r12*r12*unsr2/(2*r1)) +
310 1./(r1+1./unsr2+r12)*(-7./r1 + 2./r12) ;
311 f_delta_zec += fexp*(-5.*r1/(8.*r12*r12*r12)+5.*r1*r1*unsr2/(8.*r12*r12*r12)) ;
314 for (
int i=0 ;i<nz1-1 ; i++){
318 f_delta = f_delta + f_delta_zec ;
334 f_1_1 = r1/(8.*r12*r12*r12) + 11./(8.*r1*r12) -
335 1./(8.*r1*unsr2*unsr2*r12*r12*r12) + 7./(r1+1./unsr2+r12)/(r1+1./unsr2+r12) +
336 7./r1/(r1+1./unsr2+r12) ;
339 f_1_1_zec = 11./(8.*r1*r12) + 7./(r1+1./unsr2+r12)/(r1+1./unsr2+r12) +
340 7./r1/(r1+1./unsr2+r12) ;
341 f_1_1_zec += fexp*(r1/(8.*r12*r12*r12)-1./(8.*r1*unsr2*unsr2*r12*r12*r12)) ;
344 for (
int i=0 ; i<nz1-1 ; i++){
348 f_1_1 = f_1_1 + f_1_1_zec ;
364 f_1_12 = - 7./(2*r12*r12) + 8./(r1+1./unsr2+r12)/(r1+1./unsr2+r12) ;
367 f_1_12_zec = 8./(r1+1./unsr2+r12)/(r1+1./unsr2+r12) ;
368 f_1_12_zec += fexp*(- 7./(2*r12*r12)) ;
371 for (
int i=0 ; i<nz1-1 ; i++){
375 f_1_12 = f_1_12 + f_1_12_zec ;
391 f_12_12 = (-4./(r1+1./unsr2+r12)/(r1+1./unsr2+r12) -
392 4./r12/(r1+1./unsr2+r12)) ;
409 f_1_2 = 11./(r1+1./unsr2+r12)/(r1+1./unsr2+r12);
428 for (
int i=1 ; i<= 3 ; i++){
429 for (
int j=i ; j<= 3 ; j++){
430 hh_temp.
set(i,j) = f_delta *
hole1.
ff.
con()(i,j) + f_1_1 * nn1(i)*nn1(j)
431 + f_1_12 * 0.5 *(nn1(i) * nn12(j) + nn1(j) * nn12(i))
432 + f_12_12 * nn12(i)*nn12(j)
433 + f_1_2 * 0.5*(nn1(i)*nn2(j) + nn1(j)*nn2(i) ) ;
487 rr2.std_spectral_base() ;
504 for (
int ind=1; ind<=3; ind++){
511 rr1_1.
set(1) = xx_1 ;
512 rr1_1.
set(2) = yy_1 ;
513 rr1_1.
set(3) = zz_1 ;
514 rr1_1.std_spectral_base() ;
531 for (
int ind=1; ind<=3; ind++){
532 nn1_1.
set(ind).
set_grid_point(nz1-1,k,j,i) = nn1_1(ind).val_grid_point(1,k,j,0) ;
545 for (
int i=1 ; i<=3 ; i++){
546 nn1.set(i).import(nn1_1(i)) ;
547 nn1.set(i).set_spectral_va().set_base(nn1_1(i).get_spectral_va().get_base()) ;
565 Scalar r1sr2 (1./unsr1*unsr2) ;
576 rr21.std_spectral_base() ;
580 r21 =
sqrt( rr21(1)*rr21(1) + rr21(2)*rr21(2) + rr21(3)*rr21(3)) ;
610 double sigma = 1.*r0 ;
616 fexp =
exp( -(rr - r0)*(rr - r0)/sigma/sigma ) ;
617 for (
int ii=0; ii<nz2-1; ii++)
633 f_delta = -5.*r2/(8.*r21*r21*r21) - 15./(8.*r2*r21) +
634 5.*r2*r2*unsr1/(8.*r21*r21*r21) + 1./(r2+1./unsr1+r21)/(r2+1./unsr1+r21)*
635 (1 + r2/r21 + r21/r2 - r2sr1 - r2*r2sr1/r21 + r21*r21*unsr1/(2*r2)) +
636 1./(r2+1./unsr1+r21)*(-7./r2 + 2./r21) ;
640 f_delta_zec = - 15./(8.*r2*r21) + 1./(r2+1./unsr1+r21)/(r2+1./unsr1+r21)*
641 (1 + r2/r21 + r21/r2 - r2sr1 - r2*r2sr1/r21 + r21*r21*unsr1/(2*r2)) +
642 1./(r2+1./unsr1+r21)*(-7./r2 + 2./r21) ;
643 f_delta_zec += fexp*(-5.*r2/(8.*r21*r21*r21)+5.*r2*r2*unsr1/(8.*r21*r21*r21)) ;
646 for (
int i=0 ;i<nz2-1 ; i++){
650 f_delta = f_delta + f_delta_zec ;
666 f_2_2 = r2/(8.*r21*r21*r21) + 11./(8.*r2*r21) -
667 1./(8.*r2*unsr1*unsr1*r21*r21*r21) + 7./(r2+1./unsr1+r21)/(r2+1./unsr1+r21) +
668 7./r2/(r2+1./unsr1+r21) ;
671 f_2_2_zec = 11./(8.*r2*r21) + 7./(r2+1./unsr1+r21)/(r2+1./unsr1+r21) +
672 7./r2/(r2+1./unsr1+r21) ;
673 f_2_2_zec += fexp*(r2/(8.*r21*r21*r21)-1./(8.*r2*unsr1*unsr1*r21*r21*r21)) ;
676 for (
int i=0 ; i<nz2-1 ; i++){
680 f_2_2 = f_2_2 + f_2_2_zec ;
696 f_2_21 = - 7./(2*r21*r21) + 8./(r2+1./unsr1+r21)/(r2+1./unsr1+r21) ;
699 f_2_21_zec = 8./(r2+1./unsr1+r21)/(r2+1./unsr1+r21) ;
700 f_2_21_zec += fexp*(- 7./(2*r21*r21)) ;
703 for (
int i=0 ; i<nz2-1 ; i++){
707 f_2_21 = f_2_21 + f_2_21_zec ;
723 f_21_21 = (-4./(r2+1./unsr1+r21)/(r2+1./unsr1+r21) -
724 4./r21/(r2+1./unsr1+r21)) ;
741 f_2_1 = 11./(r2+1./unsr1+r21)/(r2+1./unsr1+r21);
760 for (
int i=1 ; i<= 3 ; i++){
761 for (
int j=i ; j<= 3 ; j++){
762 hh_temp.
set(i,j) = f_delta *
hole2.
ff.
con()(i,j) + f_2_2 * nn2(i)*nn2(j)
763 - f_2_21 * 0.5 *(nn2(i) * nn21(j) + nn2(j) * nn21(i))
764 + f_21_21 * nn21(i)*nn21(j)
765 + f_2_1 * 0.5*(nn2(i)*nn1(j) + nn2(j)*nn1(i) ) ;
850 for (
int i=1 ; i<=3 ; i++)
851 for (
int j=i ; j<=3 ; j++){
866 for (
int i=1 ; i<=3 ; i++){
867 for (
int j=i ; j<=3 ; j++){
868 hh2_1.set(i,j).import(hh2(i,j)) ;
869 hh2_1.set(i,j).set_spectral_va().set_base(hh2(i,j).get_spectral_va().get_base()) ;
875 for (
int i=1 ; i<=3 ; i++){
876 for (
int j=i ; j<=3 ; j++){
877 hh1_2.set(i,j).import(hh1(i,j)) ;
878 hh1_2.set(i,j).set_spectral_va().set_base(hh1(i,j).get_spectral_va().get_base()) ;
897 for (
int i=1 ; i<= 3 ; i++)
898 for (
int j=i ; j<= 3 ; j++){
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Metric for tensor calculation.
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
void des_coupe_z(const Scalar &uu, double z0, int nzdes, const char *title=0x0, const Scalar *defsurf=0x0, double zoom=1.2, bool draw_bound=true, int ncour=15, int nx=100, int ny=100)
Draws isocontour lines of a Scalar in a plane Z=constant.
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Cmp exp(const Cmp &)
Exponential.
Metric tgam
3 metric tilde
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Cmp sqrt(const Cmp &)
Square root.
double area_hor() const
Area of the horizon.
void annule(int l)
Sets the Cmp to zero in a given domain.
const Tbl & domain(int l) const
Read-only of the value in a given domain.
double get_ori_y() const
Returns the y coordinate of the origin.
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.
Tensor field of valence 0 (or component of a tensorial field).
double get_ori_x() const
Returns the x coordinate of the origin.
virtual double val_r(int l, double xi, double theta, double pphi) const
Returns the value of the radial coordinate r for a given in a given domain.
Single_hor hole1
Black hole one.
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
void set(const Map *mp, Mtbl *(*construct)(const Map *))
Semi-constructor from a mapping and a method.
Tensor field of valence 1.
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Tbl & set_domain(int l)
Read/write of the value in a given domain.
Sym_tensor hh_Samaya_hole2()
Calculation of the hole2 part of the Post-Newtonian correction to .
Coord tet
coordinate centered on the grid
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Coord phi
coordinate centered on the grid
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Sym_tensor hh
Deviation metric.
void set_outer_boundary(int l, double x)
Sets the value of the Scalar at the outer boundary of a given domain.
Sym_tensor hh_Samaya_hole1()
Calculation of the hole1 part of the Post-Newtonian correction to .
double get_radius() const
Returns the radius of the horizon.
int get_nzone() const
Returns the number of domains.
Cmp pow(const Cmp &, int)
Power .
Active physical coordinates and mapping derivatives.
Single_hor hole2
Black hole two.
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
void import(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
void des_meridian(const Scalar &uu, double r_min, double r_max, const char *nomy, int ngraph, const char *device=0x0, bool closeit=false, bool draw_bound=true)
Draws 5 profiles of a scalar field along various radial axes in two meridional planes and ...
double & set_grid_point(int l, int k, int j, int i)
Setting the value of the field at a given grid point.
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Map_af & mp
Affine mapping.
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Coord y
y coordinate centered on the grid
void set_hh_Samaya()
Calculation of the Post-Newtonian correction to .
Coord x
x coordinate centered on the grid
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
double get_ori_z() const
Returns the z coordinate of the origin.
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.
double ang_mom_hor() const
Angular momentum (modulo)
double radius
Radius of the horizon in LORENE's units.
Metric_flat ff
3 metric flat
Coord z
z coordinate centered on the grid
Class intended to describe valence-2 symmetric tensors.
const Valeur & get_spectral_va() const
Returns va (read only version)
Coord r
r coordinate centered on the grid
const Base_val & get_base() const
Return the bases for spectral expansions (member base )