LORENE
compobj.h
1 /*
2  * Definition of Lorene class Compobj, Compobj_QI, Star_QI, Kerr_QI, AltBH_QI, HiggsMonopole, ScalarBH
3  *
4  */
5 
6 /*
7  * Copyright (c) 2012, 2013 Claire Some, Eric Gourgoulhon
8  *
9  * This file is part of LORENE.
10  *
11  * LORENE is free software; you can redistribute it and/or modify
12  * it under the terms of the GNU General Public License version 2
13  * as published by the Free Software Foundation.
14  *
15  * LORENE is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  * GNU General Public License for more details.
19  *
20  * You should have received a copy of the GNU General Public License
21  * along with LORENE; if not, write to the Free Software
22  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23  *
24  */
25 
26 #ifndef __COMPOBJ_H_
27 #define __COMPOBJ_H_
28 
29 /*
30  * $Id: compobj.h,v 1.21 2018/11/16 14:34:34 j_novak Exp $
31  * $Log: compobj.h,v $
32  * Revision 1.21 2018/11/16 14:34:34 j_novak
33  * Changed minor points to avoid some compilation warnings.
34  *
35  * Revision 1.20 2015/11/05 17:31:21 f_vincent
36  * Updated class scalarBH.
37  *
38  * Revision 1.19 2015/10/22 09:18:35 f_vincent
39  * New class ScalarBH
40  *
41  * Revision 1.18 2014/10/13 08:52:33 j_novak
42  * Lorene classes and functions now belong to the namespace Lorene.
43  *
44  * Revision 1.17 2014/05/16 11:55:19 o_straub
45  * fixed: GYOTO output from compobj & compobj_QI
46  *
47  * Revision 1.16 2014/01/31 15:34:54 e_gourgoulhon
48  * Added members to class HiggsMonopole
49  *
50  * Revision 1.15 2014/01/29 16:29:16 e_gourgoulhon
51  * Added new class HiggsMonopole
52  *
53  * Revision 1.14 2014/01/14 20:53:39 e_gourgoulhon
54  * Updated documentation of r_isco
55  *
56  * Revision 1.13 2013/07/25 19:44:45 o_straub
57  * calculation of the marginally bound radius
58  *
59  * Revision 1.12 2013/04/17 13:01:50 e_gourgoulhon
60  * Some modifications in the class AltBH_QI
61  *
62  * Revision 1.11 2013/04/16 15:26:45 e_gourgoulhon
63  * Added class AltBH_QI
64  *
65  * Revision 1.10 2013/04/04 15:31:34 e_gourgoulhon
66  * r_isco returns now the coordinate r, not the areal r
67  *
68  * Revision 1.9 2013/04/03 12:08:57 e_gourgoulhon
69  * Added member kk to Compobj; suppressed tkij
70  *
71  * Revision 1.8 2013/04/02 23:17:18 e_gourgoulhon
72  * New class Kerr_QI
73  *
74  * Revision 1.7 2012/12/03 15:26:14 c_some
75  * Added data member m2
76  *
77  * Revision 1.6 2012/11/22 16:02:18 c_some
78  * *** empty log message ***
79  *
80  * Revision 1.5 2012/11/21 14:52:13 c_some
81  * Documentation corrected
82  *
83  * Revision 1.4 2012/11/20 16:21:16 c_some
84  * Added new class Star_QI
85  *
86  * Revision 1.3 2012/11/16 16:13:12 c_some
87  * Added new class Compobj_QI
88  *
89  * Revision 1.2 2012/11/15 20:50:41 e_gourgoulhon
90  * Corrected the documentation
91  *
92  * Revision 1.1 2012/11/15 16:20:51 c_some
93  * New class Compobj
94  *
95  *
96  * $Header: /cvsroot/Lorene/C++/Include/compobj.h,v 1.21 2018/11/16 14:34:34 j_novak Exp $
97  *
98  */
99 
100 
101 // Headers Lorene
102 #include "tensor.h"
103 #include "metric.h"
104 
105 
106 //---------------------------//
107 // base class Compobj //
108 //---------------------------//
109 
110 namespace Lorene {
129  class Compobj {
130 
131  // Data :
132  // -----
133  protected:
135  Map& mp ;
136 
139 
142 
145 
148 
151 
154 
157 
158  // Derived data :
159  // ------------
160  protected:
161  mutable double* p_adm_mass ;
162 
163  // Constructors - Destructor
164  // -------------------------
165  public:
171  Compobj(Map& map_i) ;
172 
173  Compobj(const Compobj& ) ;
174 
181  Compobj(Map& map_i, FILE* ) ;
182 
183  virtual ~Compobj() ;
184 
185 
186  // Memory management
187  // -----------------
188  protected:
190  virtual void del_deriv() const ;
191 
193  void set_der_0x0() const ;
194 
195 
196  // Mutators / assignment
197  // ---------------------
198  public:
200  void operator=(const Compobj&) ;
201 
203  Map& set_mp() {return mp; } ;
204 
205 
206  // Accessors
207  // ---------
208  public:
210  const Map& get_mp() const {return mp; } ;
211 
213  const Scalar& get_nn() const {return nn;} ;
214 
216  const Vector& get_beta() const {return beta;} ;
217 
219  const Metric& get_gamma() const {return gamma;} ;
220 
222  const Scalar& get_ener_euler() const {return ener_euler;} ;
223 
225  const Vector& get_mom_euler() const {return mom_euler;} ;
226 
228  const Sym_tensor& get_stress_euler() const {return stress_euler;} ;
229 
231  const Sym_tensor& get_kk() const {return kk;} ;
232 
233 
234 
235 
236  // Outputs
237  // -------
238  public:
239  virtual void sauve(FILE *) const ;
240 
241  void gyoto_data(const char* file_name) const ;
242 
243 
244 
246  friend ostream& operator<<(ostream& , const Compobj& ) ;
247 
248  protected:
250  virtual ostream& operator>>(ostream& ) const ;
251 
252  // Computational methods
253  // ---------------------
254  public:
256  virtual void extrinsic_curvature() ;
257 
258 
259  // Global quantities
260  // -----------------
261  public:
263  virtual double adm_mass() const ;
264  };
265 
266 
267  //---------------------------//
268  // base class Compobj_QI //
269  //---------------------------//
270 
283  class Compobj_QI : public Compobj {
284 
285  // Data :
286  // -----
287  protected:
288 
291 
294 
297 
300 
319 
320 
321  // Derived data :
322  // ------------
323  protected:
324  mutable double* p_angu_mom ;
325  mutable double* p_r_isco ;
326  mutable double* p_f_isco ;
327  mutable double* p_espec_isco ;
330  mutable double* p_lspec_isco ;
331  mutable double* p_r_mb ;
332 
333  // Constructors - Destructor
334  // -------------------------
335  public:
341  Compobj_QI(Map& map_i) ;
342 
343  Compobj_QI(const Compobj_QI& ) ;
344 
351  Compobj_QI(Map& map_i, FILE* ) ;
352 
353  virtual ~Compobj_QI() ;
354 
355 
356  // Memory management
357  // -----------------
358  protected:
360  virtual void del_deriv() const ;
361 
363  void set_der_0x0() const ;
364 
365 
366  // Mutators / assignment
367  // ---------------------
368  public:
370  void operator=(const Compobj_QI&) ;
371 
372 
373  // Accessors
374  // ---------
375  public:
376 
378  const Scalar& get_bbb() const {return bbb;} ;
379 
381  const Scalar& get_a_car() const {return a_car;} ;
382 
384  const Scalar& get_b_car() const {return b_car;} ;
385 
387  const Scalar& get_nphi() const {return nphi;} ;
388 
389 
407  const Scalar& get_ak_car() const {return ak_car;} ;
408 
409 
410 
411 
412 
413 
414  // Outputs
415  // -------
416  public:
417  virtual void sauve(FILE *) const ;
418 
419  void gyoto_data(const char* file_name) const ;
420 
421 
422  protected:
424  virtual ostream& operator>>(ostream& ) const ;
425 
426  // Global quantities
427  // -----------------
428  public:
429  virtual double angu_mom() const ;
430 
441  virtual double r_isco(int lmin, ostream* ost = 0x0) const ;
442 
444  virtual double f_isco(int lmin) const ;
445 
447  virtual double espec_isco(int lmin) const ;
448 
450  virtual double lspec_isco(int lmin) const ;
451 
453  virtual double r_mb(int lmin, ostream* ost = 0x0) const ;
454 
455 
456  // Computational routines
457  // ----------------------
458 
463  virtual void update_metric() ;
464 
468  virtual void extrinsic_curvature() ;
469 
470  };
471 
472 
473  //--------------------------//
474  // base class Star_QI //
475  //--------------------------//
476 
490  class Star_QI : public Compobj_QI {
491 
492  // Data :
493  // -----
494  protected:
495 
499 
504 
509 
514 
517 
520 
533 
543 
549 
555 
560 
565 
573 
582 
583 
584  // Derived data :
585  // ------------
586  protected:
587 
588  mutable double* p_grv2 ;
589  mutable double* p_grv3 ;
590  mutable double* p_mom_quad ;
591  mutable double* p_mass_g ;
592 
593 
594  // Constructors - Destructor
595  // -------------------------
596  public:
602  Star_QI(Map& mp_i) ;
603 
604 
605  Star_QI(const Star_QI& ) ;
606 
613  Star_QI(Map& mp_i, FILE* fich) ;
614 
615  virtual ~Star_QI() ;
616 
617 
618  // Memory management
619  // -----------------
620  protected:
622  virtual void del_deriv() const ;
623 
625  virtual void set_der_0x0() const ;
626 
627  // Mutators / assignment
628  // ---------------------
629  public:
631  void operator=(const Star_QI& ) ;
632 
633  // Accessors
634  // ---------
635  public:
636 
639  const Scalar& get_logn() const {return logn;} ;
640 
641 
645  const Scalar& get_tnphi() const {return tnphi;} ;
646 
650  const Scalar& get_nuf() const {return nuf;} ;
651 
655  const Scalar& get_nuq() const {return nuq;} ;
656 
658  const Scalar& get_dzeta() const {return dzeta;} ;
659 
661  const Scalar& get_tggg() const {return tggg;} ;
662 
675  const Vector& get_w_shift() const {return w_shift;} ;
676 
689  const Scalar& get_khi_shift() const {return khi_shift;} ;
690 
691 
692  // Outputs
693  // -------
694  public:
695  virtual void sauve(FILE* ) const ;
696 
697  protected:
699  virtual ostream& operator>>(ostream& ) const ;
700 
701  // Global quantities
702  // -----------------
703  public:
704 
705  virtual double mass_g() const ;
706  virtual double angu_mom() const ;
707 
711  virtual double grv2() const ;
712 
724  virtual double grv3(ostream* ost = 0x0) const ;
725 
735  virtual double mom_quad() const ;
736 
737 
738  // Computational routines
739  // ----------------------
740  public:
741 
752  void update_metric() ;
753 
762  void fait_shift() ;
763 
767  void fait_nphi() ;
768 
798  static double lambda_grv2(const Scalar& sou_m, const Scalar& sou_q) ;
799 
800  };
801 
802 
803  //---------------------//
804  // class Kerr_QI //
805  //---------------------//
806 
819  class Kerr_QI : public Compobj_QI {
820 
821  // Data :
822  // -----
823  protected:
824 
827  double mm ;
828 
831  double aa ;
832 
833 
834  // Derived data :
835  // ------------
836  protected:
837 
838 
839  // Constructors - Destructor
840  // -------------------------
841  public:
849  Kerr_QI(Map& mp_i, double mass, double a_over_m) ;
850 
851 
852  Kerr_QI(const Kerr_QI& ) ;
853 
860  Kerr_QI(Map& mp_i, FILE* fich) ;
861 
862  virtual ~Kerr_QI() ;
863 
864  // Memory management
865  // -----------------
866  protected:
868  virtual void del_deriv() const ;
869 
871  virtual void set_der_0x0() const ;
872 
873  // Mutators / assignment
874  // ---------------------
875  public:
877  void operator=(const Kerr_QI& ) ;
878 
879  // Accessors
880  // ---------
881  public:
882 
883  // Outputs
884  // -------
885  public:
886  virtual void sauve(FILE* ) const ;
887 
888  protected:
890  virtual ostream& operator>>(ostream& ) const ;
891 
892  // Global quantities
893  // -----------------
894  public:
895 
896 
897  // Computational routines
898  // ----------------------
899  public:
900 
901 
902  };
903 
904  //-------------------//
905  // class AltBH_QI //
906  //-------------------//
907 
920  class AltBH_QI : public Compobj_QI {
921 
922  // Data :
923  // -----
924  protected:
925 
926  char description1[256] ;
927  char description2[256] ;
928  double a_spin ;
929 
931 
932  // Derived data :
933  // ------------
934  protected:
935 
936 
937  // Constructors - Destructor
938  // -------------------------
939  public:
947  AltBH_QI(Map& mp_i, const char* file_name, double a_spin_i) ;
948 
949 
950  AltBH_QI(const AltBH_QI& ) ;
951 
958  AltBH_QI(Map& mp_i, FILE* fich) ;
959 
960  virtual ~AltBH_QI() ;
961 
962  // Memory management
963  // -----------------
964  protected:
966  virtual void del_deriv() const ;
967 
969  virtual void set_der_0x0() const ;
970 
971  // Mutators / assignment
972  // ---------------------
973  public:
975  void operator=(const AltBH_QI& ) ;
976 
977  // Accessors
978  // ---------
979  public:
980 
982  const Scalar& get_krphi() const {return krphi;} ;
983 
984  // Outputs
985  // -------
986  public:
987  virtual void sauve(FILE* ) const ;
988 
989  protected:
991  virtual ostream& operator>>(ostream& ) const ;
992 
993  // Global quantities
994  // -----------------
995  public:
996 
997 
998  // Computational routines
999  // ----------------------
1000  public:
1001 
1003  virtual void extrinsic_curvature() ;
1004 
1005  };
1006 
1007  //-------------------//
1008  // class ScalarBH //
1009  //-------------------//
1010 
1023  class ScalarBH : public Compobj {
1024 
1025  // Data :
1026  // -----
1027  protected:
1028 
1029  //char description1[256] ; ///< String describing the model
1030  // char description2[256] ; ///< String describing the model
1036  double rHor ;
1037 
1038  // Constructors - Destructor
1039  // -------------------------
1040  public:
1048  ScalarBH(Map& mp_i, const char* file_name) ;
1049 
1050  ScalarBH(const ScalarBH& ) ;
1051 
1058  ScalarBH(Map& mp_i, FILE* fich) ;
1059 
1060  virtual ~ScalarBH() ;
1061 
1062  // Memory management
1063  // -----------------
1064  protected:
1066  virtual void del_deriv() const ;
1067 
1069  virtual void set_der_0x0() const ;
1070 
1071  // Mutators / assignment
1072  // ---------------------
1073  public:
1075  void operator=(const ScalarBH& ) ;
1076 
1077  // Accessors
1078  // ---------
1079  public:
1081  const Scalar& get_ff0() const {return ff0; } ;
1082  const Scalar& get_ff1() const {return ff1; } ;
1083  const Scalar& get_ff2() const {return ff2; } ;
1084  const Scalar& get_ww() const {return ww; } ;
1085  const Scalar& get_sfield() const {return sfield; } ;
1086  double get_rHor() const {return rHor; } ;
1087 
1088  // Outputs
1089  // -------
1090  public:
1091  virtual void sauve(FILE* ) const ;
1092 
1093  protected:
1095  virtual ostream& operator>>(ostream& ) const ;
1096 
1097  // Global quantities
1098  // -----------------
1099  public:
1100 
1101 
1102  // Computational routines
1103  // ----------------------
1104  public:
1105  virtual void update_metric();
1106  };
1107 
1108 
1109  //------------------------//
1110  // class HiggsMonopole //
1111  //------------------------//
1112 
1119  class HiggsMonopole : public Compobj {
1120 
1121  // Data :
1122  // -----
1123  protected:
1124 
1125  char description1[256] ;
1126  char description2[256] ;
1127 
1129 
1131 
1133 
1134  // Derived data :
1135  // ------------
1136  protected:
1137 
1138 
1139  // Constructors - Destructor
1140  // -------------------------
1141  public:
1148  HiggsMonopole(Map& mp_i, const char* file_name) ;
1149 
1150  HiggsMonopole(const HiggsMonopole& ) ;
1151 
1152  virtual ~HiggsMonopole() ;
1153 
1154  // Memory management
1155  // -----------------
1156  protected:
1158  // virtual void del_deriv() const ;
1159 
1161  // virtual void set_der_0x0() const ;
1162 
1163  // Mutators / assignment
1164  // ---------------------
1165  public:
1167  // void operator=(const AltBH_QI& ) ;
1168 
1169  // Accessors
1170  // ---------
1171  public:
1172 
1174  const Scalar& get_higgs() const {return hh;} ;
1175 
1177  const Scalar& get_grr() const {return grr;} ;
1178 
1180  const Scalar& get_press() const {return press;} ;
1181 
1182  // Outputs
1183  // -------
1184  public:
1185  // virtual void sauve(FILE* ) const ; ///< Save in a file
1186 
1187  protected:
1189  virtual ostream& operator>>(ostream& ) const ;
1190 
1191  // Global quantities
1192  // -----------------
1193  public:
1194 
1195 
1196  // Computational routines
1197  // ----------------------
1198  public:
1199 
1201  // virtual void extrinsic_curvature() ;
1202 
1203  };
1204 
1205 
1206 
1207 
1208 }
1209 #endif
Base class for axisymmetric stationary compact objects in Quasi-Isotropic coordinates (under developm...
Definition: compobj.h:283
AltBH_QI(Map &mp_i, const char *file_name, double a_spin_i)
Standard constructor.
Definition: altBH_QI.C:72
Scalar logn
Logarithm of the lapse N .
Definition: compobj.h:498
virtual void extrinsic_curvature()
Computation of the extrinsic curvature.
Definition: altBH_QI.C:272
Scalar ssjm1_dzeta
Effective source at the previous step for the resolution of the Poisson equation for dzeta ...
Definition: compobj.h:559
Metric for tensor calculation.
Definition: metric.h:90
Scalar khi_shift
Scalar used in the decomposition of shift , following Shibata&#39;s prescription [Prog.
Definition: compobj.h:542
virtual ~ScalarBH()
Destructor.
Definition: scalarBH.C:392
Scalar ff2
Metric field F_2 of Herdeiro & Radu (2015)
Definition: compobj.h:1033
Vector mom_euler
Total 3-momentum density in the Eulerian frame.
Definition: compobj.h:150
virtual ~HiggsMonopole()
Destructor.
virtual void sauve(FILE *) const
Save in a file.
Definition: scalarBH.C:436
const Vector & get_w_shift() const
Returns the vector used in the decomposition of shift , following Shibata&#39;s prescription [Prog...
Definition: compobj.h:675
Scalar tnphi
Component of the shift vector.
Definition: compobj.h:503
Scalar press
Fluid pressure.
Definition: compobj.h:1132
virtual void sauve(FILE *) const
Save in a file.
Definition: kerr_QI.C:201
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: altBH_QI.C:210
virtual double adm_mass() const
ADM mass (computed as a surface integral at spatial infinity)
Definition: compobj.C:313
Compobj_QI(Map &map_i)
Standard constructor.
Definition: compobj_QI.C:89
Kerr spacetime in Quasi-Isotropic coordinates (under development).
Definition: compobj.h:819
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: compobj_QI.C:181
double a_spin
Spin parameter of the model.
Definition: compobj.h:928
Scalar ak_car
Scalar .
Definition: compobj.h:318
friend ostream & operator<<(ostream &, const Compobj &)
Display.
Definition: compobj.C:236
Lorene prototypes.
Definition: app_hor.h:67
virtual ~Compobj()
Destructor.
Definition: compobj.C:147
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
virtual void extrinsic_curvature()
Computes the extrinsic curvature and ak_car from nphi , nn and b_car .
Definition: compobj_QI.C:342
double * p_adm_mass
ADM mass.
Definition: compobj.h:161
virtual double mass_g() const
Gravitational mass.
Base class for coordinate mappings.
Definition: map.h:688
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: altBH_QI.C:251
virtual double lspec_isco(int lmin) const
Angular momentum of a particle at the ISCO.
Scalar nuq
Part of the Metric potential = logn generated by the quadratic terms.
Definition: compobj.h:513
const Metric & get_gamma() const
Returns the 3-metric .
Definition: compobj.h:219
double * p_lspec_isco
Specific angular momentum of a particle at the ISCO.
Definition: compobj.h:330
virtual double mom_quad() const
Quadrupole moment.
Scalar ww
Metric field W of Herdeiro & Radu (2015)
Definition: compobj.h:1034
virtual double espec_isco(int lmin) const
Energy of a particle at the ISCO.
Scalar ssjm1_nuq
Effective source at the previous step for the resolution of the Poisson equation for nuq by means of ...
Definition: compobj.h:554
void gyoto_data(const char *file_name) const
Save in a file for GYOTO.
Definition: compobj_QI.C:235
const Scalar & get_nuf() const
Returns the part of the Metric potential = logn generated by the matter terms.
Definition: compobj.h:650
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition: compobj.h:572
void operator=(const Kerr_QI &)
Assignment to another Kerr_QI.
Definition: kerr_QI.C:184
Sym_tensor kk
Extrinsic curvature tensor .
Definition: compobj.h:156
Tensor field of valence 1.
Definition: vector.h:188
virtual void extrinsic_curvature()
Computation of the extrinsic curvature.
Definition: compobj.C:293
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: scalarBH.C:403
Scalar a_car
Square of the metric factor A.
Definition: compobj.h:290
Map & set_mp()
Read/write of the mapping.
Definition: compobj.h:203
double * p_grv2
Error on the virial identity GRV2.
Definition: compobj.h:588
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: scalarBH.C:412
const Scalar & get_nuq() const
Returns the Part of the Metric potential = logn generated by the quadratic terms.
Definition: compobj.h:655
virtual double angu_mom() const
Angular momentum.
Scalar nuf
Part of the Metric potential = logn generated by the matter terms.
Definition: compobj.h:508
void operator=(const ScalarBH &)
Assignment to another ScalarBH.
Definition: scalarBH.C:422
virtual double f_isco(int lmin) const
Orbital frequency at the innermost stable circular orbit (ISCO).
void gyoto_data(const char *file_name) const
Save in a file for GYOTO.
Definition: compobj.C:213
Base class for axisymmetric stationary compact stars in Quasi-Isotropic coordinates (under developmen...
Definition: compobj.h:490
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
const Scalar & get_khi_shift() const
Returns the scalar used in the decomposition of shift following Shibata&#39;s prescription [Prog...
Definition: compobj.h:689
Scalar nphi
Metric coefficient .
Definition: compobj.h:299
Scalar b_car
Square of the metric factor B.
Definition: compobj.h:296
Scalar ff1
Metric field F_1 of Herdeiro & Radu (2015)
Definition: compobj.h:1032
Base class for stationary compact objects (under development).
Definition: compobj.h:129
const Scalar & get_tggg() const
Returns the Metric potential .
Definition: compobj.h:661
double * p_angu_mom
Angular momentum.
Definition: compobj.h:324
const Scalar & get_bbb() const
Returns the metric factor B.
Definition: compobj.h:378
const Vector & get_mom_euler() const
Returns the total 3-momentum density in the Eulerian frame.
Definition: compobj.h:225
Scalar tggg
Metric potential .
Definition: compobj.h:519
virtual ~AltBH_QI()
Destructor.
Definition: altBH_QI.C:199
void operator=(const Compobj &)
Assignment to another Compobj.
Definition: compobj.C:178
virtual double r_isco(int lmin, ostream *ost=0x0) const
Coordinate r of the innermost stable circular orbit (ISCO).
Scalar bbb
Metric factor B.
Definition: compobj.h:293
virtual void sauve(FILE *) const
Save in a file.
Definition: altBH_QI.C:243
void fait_shift()
Computes shift from w_shift and khi_shift according to Shibata&#39;s prescription [Prog.
Definition: star_QI.C:343
static double lambda_grv2(const Scalar &sou_m, const Scalar &sou_q)
Computes the coefficient which ensures that the GRV2 virial identity is satisfied.
double * p_espec_isco
Specific energy of a particle at the ISCO.
Definition: compobj.h:328
void operator=(const AltBH_QI &)
Assignment to another AltBH_QI.
Definition: altBH_QI.C:229
Scalar grr
Metric coefficient g_rr.
Definition: compobj.h:1130
Scalar ssjm1_tggg
Effective source at the previous step for the resolution of the Poisson equation for tggg ...
Definition: compobj.h:564
const Scalar & get_ener_euler() const
Returns the total energy density E in the Eulerian frame.
Definition: compobj.h:222
Vector beta
Shift vector .
Definition: compobj.h:141
Star_QI(Map &mp_i)
Standard constructor.
Definition: star_QI.C:71
ScalarBH(Map &mp_i, const char *file_name)
Standard constructor.
Definition: scalarBH.C:75
void operator=(const Compobj_QI &)
Assignment to another Compobj_QI.
Definition: compobj_QI.C:198
const Scalar & get_ff0() const
Returns f0.
Definition: compobj.h:1081
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: altBH_QI.C:219
char description2[256]
String describing the model.
Definition: compobj.h:1126
const Scalar & get_higgs() const
Deletes all the derived quantities.
Definition: compobj.h:1174
Scalar ff0
Metric field F_0 of Herdeiro & Radu (2015)
Definition: compobj.h:1031
double rHor
Event horizon coordinate radius.
Definition: compobj.h:1036
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: compobj.C:158
Vector ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for ...
Definition: compobj.h:581
virtual ~Compobj_QI()
Destructor.
Definition: compobj_QI.C:155
Scalar dzeta
Metric potential .
Definition: compobj.h:516
const Scalar & get_a_car() const
Returns the square of the metric factor A.
Definition: compobj.h:381
const Scalar & get_press() const
Returns the fluid pressure.
Definition: compobj.h:1180
double * p_r_isco
Coordinate r of the ISCO.
Definition: compobj.h:325
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star_QI.C:238
double * p_grv3
Error on the virial identity GRV3.
Definition: compobj.h:589
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: compobj.C:242
char description1[256]
String describing the model.
Definition: compobj.h:1125
const Sym_tensor & get_stress_euler() const
Returns the stress tensor with respect to the Eulerian observer.
Definition: compobj.h:228
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: scalarBH.C:444
double * p_r_mb
Coordinate r of the marginally bound orbit.
Definition: compobj.h:331
double * p_mass_g
Gravitational mass (ADM mass as a volume integral)
Definition: compobj.h:591
double * p_mom_quad
Quadrupole moment.
Definition: compobj.h:590
const Scalar & get_dzeta() const
Returns the Metric potential .
Definition: compobj.h:658
Compobj(Map &map_i)
Standard constructor.
Definition: compobj.C:85
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: kerr_QI.C:174
const Scalar & get_grr() const
Returns the metric coefficient g_rr.
Definition: compobj.h:1177
const Scalar & get_logn() const
Returns the logarithm of the lapse N.
Definition: compobj.h:639
Kerr_QI(Map &mp_i, double mass, double a_over_m)
Standard constructor.
Definition: kerr_QI.C:71
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: kerr_QI.C:165
Alternative black hole spacetime in Quasi-Isotropic coordinates (under development).
Definition: compobj.h:920
double * p_f_isco
Orbital frequency of the ISCO.
Definition: compobj.h:326
const Scalar & get_tnphi() const
Returns the component of the shift vector.
Definition: compobj.h:645
Scalar ener_euler
Total energy density E in the Eulerian frame.
Definition: compobj.h:147
Black hole with scalar hair spacetime (under development).
Definition: compobj.h:1023
virtual void update_metric()
Updates the 3-metric from A and B and the shift vector from .
Definition: compobj_QI.C:305
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: compobj_QI.C:267
Scalar ssjm1_nuf
Effective source at the previous step for the resolution of the Poisson equation for nuf by means of ...
Definition: compobj.h:548
Scalar krphi
K_{(r)(phi)} read in the file.
Definition: compobj.h:930
Scalar nn
Lapse function N .
Definition: compobj.h:138
const Scalar & get_krphi() const
Returns K_{(r)(phi)}/sin(theta).
Definition: compobj.h:982
const Scalar & get_nn() const
Returns the lapse function N .
Definition: compobj.h:213
Scalar sfield
Scalar field (modulus of Phi)
Definition: compobj.h:1035
virtual double r_mb(int lmin, ostream *ost=0x0) const
Coordinate r of the marginally bound circular orbit (R_mb).
Metric gamma
3-metric
Definition: compobj.h:144
virtual double grv2() const
Error on the virial identity GRV2.
Scalar hh
Higgs scalar field.
Definition: compobj.h:1128
const Sym_tensor & get_kk() const
Returns the extrinsic curvature tensor .
Definition: compobj.h:231
virtual double angu_mom() const
Angular momentum.
virtual void sauve(FILE *) const
Save in a file.
Definition: compobj_QI.C:219
Higgs monopole (under development).
Definition: compobj.h:1119
double mm
mass parameter
Definition: compobj.h:827
HiggsMonopole(Map &mp_i, const char *file_name)
Standard constructor.
virtual ~Kerr_QI()
Destructor.
Definition: kerr_QI.C:154
char description1[256]
String describing the model.
Definition: compobj.h:926
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_QI.C:225
virtual void sauve(FILE *) const
Save in a file.
Definition: compobj.C:199
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: compobj.C:166
const Scalar & get_b_car() const
Returns the square of the metric factor B.
Definition: compobj.h:384
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: kerr_QI.C:209
char description2[256]
String describing the model.
Definition: compobj.h:927
double aa
angular momentum parameter
Definition: compobj.h:831
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: compobj_QI.C:166
void operator=(const Star_QI &)
Assignment to another Star_QI.
Definition: star_QI.C:253
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:226
const Scalar & get_ak_car() const
Returns the scalar .
Definition: compobj.h:407
void fait_nphi()
Computes tnphi and nphi from the Cartesian components of the shift, stored in shift ...
Definition: star_QI.C:389
void update_metric()
Computes metric coefficients from known potentials.
Definition: star_QI.C:413
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: star_QI.C:303
const Vector & get_beta() const
Returns the shift vector .
Definition: compobj.h:216
const Scalar & get_nphi() const
Returns the metric coefficient .
Definition: compobj.h:387
const Map & get_mp() const
Returns the mapping.
Definition: compobj.h:210
Vector w_shift
Vector used in the decomposition of shift , following Shibata&#39;s prescription [Prog.
Definition: compobj.h:532
Map & mp
Mapping describing the coordinate system (r,theta,phi)
Definition: compobj.h:135
Sym_tensor stress_euler
Stress tensor with respect to the Eulerian observer.
Definition: compobj.h:153
virtual void sauve(FILE *) const
Save in a file.
Definition: star_QI.C:282
virtual ~Star_QI()
Destructor.
Definition: star_QI.C:214