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.20 2015/11/05 17:31:21 f_vincent Exp $
31  * $Log: compobj.h,v $
32  * Revision 1.20 2015/11/05 17:31:21 f_vincent
33  * Updated class scalarBH.
34  *
35  * Revision 1.19 2015/10/22 09:18:35 f_vincent
36  * New class ScalarBH
37  *
38  * Revision 1.18 2014/10/13 08:52:33 j_novak
39  * Lorene classes and functions now belong to the namespace Lorene.
40  *
41  * Revision 1.17 2014/05/16 11:55:19 o_straub
42  * fixed: GYOTO output from compobj & compobj_QI
43  *
44  * Revision 1.16 2014/01/31 15:34:54 e_gourgoulhon
45  * Added members to class HiggsMonopole
46  *
47  * Revision 1.15 2014/01/29 16:29:16 e_gourgoulhon
48  * Added new class HiggsMonopole
49  *
50  * Revision 1.14 2014/01/14 20:53:39 e_gourgoulhon
51  * Updated documentation of r_isco
52  *
53  * Revision 1.13 2013/07/25 19:44:45 o_straub
54  * calculation of the marginally bound radius
55  *
56  * Revision 1.12 2013/04/17 13:01:50 e_gourgoulhon
57  * Some modifications in the class AltBH_QI
58  *
59  * Revision 1.11 2013/04/16 15:26:45 e_gourgoulhon
60  * Added class AltBH_QI
61  *
62  * Revision 1.10 2013/04/04 15:31:34 e_gourgoulhon
63  * r_isco returns now the coordinate r, not the areal r
64  *
65  * Revision 1.9 2013/04/03 12:08:57 e_gourgoulhon
66  * Added member kk to Compobj; suppressed tkij
67  *
68  * Revision 1.8 2013/04/02 23:17:18 e_gourgoulhon
69  * New class Kerr_QI
70  *
71  * Revision 1.7 2012/12/03 15:26:14 c_some
72  * Added data member m2
73  *
74  * Revision 1.6 2012/11/22 16:02:18 c_some
75  * *** empty log message ***
76  *
77  * Revision 1.5 2012/11/21 14:52:13 c_some
78  * Documentation corrected
79  *
80  * Revision 1.4 2012/11/20 16:21:16 c_some
81  * Added new class Star_QI
82  *
83  * Revision 1.3 2012/11/16 16:13:12 c_some
84  * Added new class Compobj_QI
85  *
86  * Revision 1.2 2012/11/15 20:50:41 e_gourgoulhon
87  * Corrected the documentation
88  *
89  * Revision 1.1 2012/11/15 16:20:51 c_some
90  * New class Compobj
91  *
92  *
93  * $Header: /cvsroot/Lorene/C++/Include/compobj.h,v 1.20 2015/11/05 17:31:21 f_vincent Exp $
94  *
95  */
96 
97 
98 // Headers Lorene
99 #include "tensor.h"
100 #include "metric.h"
101 
102 
103 //---------------------------//
104 // base class Compobj //
105 //---------------------------//
106 
107 namespace Lorene {
126  class Compobj {
127 
128  // Data :
129  // -----
130  protected:
132  Map& mp ;
133 
136 
139 
142 
145 
148 
151 
154 
155  // Derived data :
156  // ------------
157  protected:
158  mutable double* p_adm_mass ;
159 
160  // Constructors - Destructor
161  // -------------------------
162  public:
168  Compobj(Map& map_i) ;
169 
170  Compobj(const Compobj& ) ;
171 
178  Compobj(Map& map_i, FILE* ) ;
179 
180  virtual ~Compobj() ;
181 
182 
183  // Memory management
184  // -----------------
185  protected:
187  virtual void del_deriv() const ;
188 
190  void set_der_0x0() const ;
191 
192 
193  // Mutators / assignment
194  // ---------------------
195  public:
197  void operator=(const Compobj&) ;
198 
200  Map& set_mp() {return mp; } ;
201 
202 
203  // Accessors
204  // ---------
205  public:
207  const Map& get_mp() const {return mp; } ;
208 
210  const Scalar& get_nn() const {return nn;} ;
211 
213  const Vector& get_beta() const {return beta;} ;
214 
216  const Metric& get_gamma() const {return gamma;} ;
217 
219  const Scalar& get_ener_euler() const {return ener_euler;} ;
220 
222  const Vector& get_mom_euler() const {return mom_euler;} ;
223 
225  const Sym_tensor& get_stress_euler() const {return stress_euler;} ;
226 
228  const Sym_tensor& get_kk() const {return kk;} ;
229 
230 
231 
232 
233  // Outputs
234  // -------
235  public:
236  virtual void sauve(FILE *) const ;
237 
238  void gyoto_data(const char* file_name) const ;
239 
240 
241 
243  friend ostream& operator<<(ostream& , const Compobj& ) ;
244 
245  protected:
247  virtual ostream& operator>>(ostream& ) const ;
248 
249  // Computational methods
250  // ---------------------
251  public:
253  virtual void extrinsic_curvature() ;
254 
255 
256  // Global quantities
257  // -----------------
258  public:
260  virtual double adm_mass() const ;
261  };
262 
263 
264  //---------------------------//
265  // base class Compobj_QI //
266  //---------------------------//
267 
280  class Compobj_QI : public Compobj {
281 
282  // Data :
283  // -----
284  protected:
285 
288 
291 
294 
297 
316 
317 
318  // Derived data :
319  // ------------
320  protected:
321  mutable double* p_angu_mom ;
322  mutable double* p_r_isco ;
323  mutable double* p_f_isco ;
325  mutable double* p_espec_isco ;
327  mutable double* p_lspec_isco ;
328  mutable double* p_r_mb ;
329 
330  // Constructors - Destructor
331  // -------------------------
332  public:
338  Compobj_QI(Map& map_i) ;
339 
340  Compobj_QI(const Compobj_QI& ) ;
341 
348  Compobj_QI(Map& map_i, FILE* ) ;
349 
350  virtual ~Compobj_QI() ;
351 
352 
353  // Memory management
354  // -----------------
355  protected:
357  virtual void del_deriv() const ;
358 
360  void set_der_0x0() const ;
361 
362 
363  // Mutators / assignment
364  // ---------------------
365  public:
367  void operator=(const Compobj_QI&) ;
368 
369 
370  // Accessors
371  // ---------
372  public:
373 
375  const Scalar& get_bbb() const {return bbb;} ;
376 
378  const Scalar& get_a_car() const {return a_car;} ;
379 
381  const Scalar& get_b_car() const {return b_car;} ;
382 
384  const Scalar& get_nphi() const {return nphi;} ;
385 
386 
404  const Scalar& get_ak_car() const {return ak_car;} ;
405 
406 
407 
408 
409 
410 
411  // Outputs
412  // -------
413  public:
414  virtual void sauve(FILE *) const ;
415 
416  void gyoto_data(const char* file_name) const ;
417 
418 
419  protected:
421  virtual ostream& operator>>(ostream& ) const ;
422 
423  // Global quantities
424  // -----------------
425  public:
426  virtual double angu_mom() const ;
427 
438  virtual double r_isco(int lmin, ostream* ost = 0x0) const ;
439 
441  virtual double f_isco(int lmin) const ;
442 
444  virtual double espec_isco(int lmin) const ;
445 
447  virtual double lspec_isco(int lmin) const ;
448 
450  virtual double r_mb(int lmin, ostream* ost = 0x0) const ;
451 
452 
453  // Computational routines
454  // ----------------------
455 
460  virtual void update_metric() ;
461 
465  virtual void extrinsic_curvature() ;
466 
467  };
468 
469 
470  //--------------------------//
471  // base class Star_QI //
472  //--------------------------//
473 
487  class Star_QI : public Compobj_QI {
488 
489  // Data :
490  // -----
491  protected:
492 
496 
501 
506 
511 
514 
517 
530 
540 
546 
552 
557 
562 
570 
579 
580 
581  // Derived data :
582  // ------------
583  protected:
584 
585  mutable double* p_grv2 ;
586  mutable double* p_grv3 ;
587  mutable double* p_mom_quad ;
588  mutable double* p_mass_g ;
589 
590 
591  // Constructors - Destructor
592  // -------------------------
593  public:
599  Star_QI(Map& mp_i) ;
600 
601 
602  Star_QI(const Star_QI& ) ;
603 
610  Star_QI(Map& mp_i, FILE* fich) ;
611 
612  virtual ~Star_QI() ;
613 
614 
615  // Memory management
616  // -----------------
617  protected:
619  virtual void del_deriv() const ;
620 
622  virtual void set_der_0x0() const ;
623 
624  // Mutators / assignment
625  // ---------------------
626  public:
628  void operator=(const Star_QI& ) ;
629 
630  // Accessors
631  // ---------
632  public:
633 
636  const Scalar& get_logn() const {return logn;} ;
637 
638 
642  const Scalar& get_tnphi() const {return tnphi;} ;
643 
647  const Scalar& get_nuf() const {return nuf;} ;
648 
652  const Scalar& get_nuq() const {return nuq;} ;
653 
655  const Scalar& get_dzeta() const {return dzeta;} ;
656 
658  const Scalar& get_tggg() const {return tggg;} ;
659 
672  const Vector& get_w_shift() const {return w_shift;} ;
673 
686  const Scalar& get_khi_shift() const {return khi_shift;} ;
687 
688 
689  // Outputs
690  // -------
691  public:
692  virtual void sauve(FILE* ) const ;
693 
694  protected:
696  virtual ostream& operator>>(ostream& ) const ;
697 
698  // Global quantities
699  // -----------------
700  public:
701 
702  virtual double mass_g() const ;
703  virtual double angu_mom() const ;
704 
708  virtual double grv2() const ;
709 
721  virtual double grv3(ostream* ost = 0x0) const ;
722 
732  virtual double mom_quad() const ;
733 
734 
735  // Computational routines
736  // ----------------------
737  public:
738 
749  void update_metric() ;
750 
759  void fait_shift() ;
760 
764  void fait_nphi() ;
765 
795  static double lambda_grv2(const Scalar& sou_m, const Scalar& sou_q) ;
796 
797  };
798 
799 
800  //---------------------//
801  // class Kerr_QI //
802  //---------------------//
803 
816  class Kerr_QI : public Compobj_QI {
817 
818  // Data :
819  // -----
820  protected:
821 
824  double mm ;
825 
828  double aa ;
829 
830 
831  // Derived data :
832  // ------------
833  protected:
834 
835 
836  // Constructors - Destructor
837  // -------------------------
838  public:
846  Kerr_QI(Map& mp_i, double mass, double a_over_m) ;
847 
848 
849  Kerr_QI(const Kerr_QI& ) ;
850 
857  Kerr_QI(Map& mp_i, FILE* fich) ;
858 
859  virtual ~Kerr_QI() ;
860 
861  // Memory management
862  // -----------------
863  protected:
865  virtual void del_deriv() const ;
866 
868  virtual void set_der_0x0() const ;
869 
870  // Mutators / assignment
871  // ---------------------
872  public:
874  void operator=(const Kerr_QI& ) ;
875 
876  // Accessors
877  // ---------
878  public:
879 
880  // Outputs
881  // -------
882  public:
883  virtual void sauve(FILE* ) const ;
884 
885  protected:
887  virtual ostream& operator>>(ostream& ) const ;
888 
889  // Global quantities
890  // -----------------
891  public:
892 
893 
894  // Computational routines
895  // ----------------------
896  public:
897 
898 
899  };
900 
901  //-------------------//
902  // class AltBH_QI //
903  //-------------------//
904 
917  class AltBH_QI : public Compobj_QI {
918 
919  // Data :
920  // -----
921  protected:
922 
923  char description1[256] ;
924  char description2[256] ;
925  double a_spin ;
926 
928 
929  // Derived data :
930  // ------------
931  protected:
932 
933 
934  // Constructors - Destructor
935  // -------------------------
936  public:
944  AltBH_QI(Map& mp_i, const char* file_name, double a_spin_i) ;
945 
946 
947  AltBH_QI(const AltBH_QI& ) ;
948 
955  AltBH_QI(Map& mp_i, FILE* fich) ;
956 
957  virtual ~AltBH_QI() ;
958 
959  // Memory management
960  // -----------------
961  protected:
963  virtual void del_deriv() const ;
964 
966  virtual void set_der_0x0() const ;
967 
968  // Mutators / assignment
969  // ---------------------
970  public:
972  void operator=(const AltBH_QI& ) ;
973 
974  // Accessors
975  // ---------
976  public:
977 
979  const Scalar& get_krphi() const {return krphi;} ;
980 
981  // Outputs
982  // -------
983  public:
984  virtual void sauve(FILE* ) const ;
985 
986  protected:
988  virtual ostream& operator>>(ostream& ) const ;
989 
990  // Global quantities
991  // -----------------
992  public:
993 
994 
995  // Computational routines
996  // ----------------------
997  public:
998 
1000  virtual void extrinsic_curvature() ;
1001 
1002  };
1003 
1004  //-------------------//
1005  // class ScalarBH //
1006  //-------------------//
1007 
1020  class ScalarBH : public Compobj {
1021 
1022  // Data :
1023  // -----
1024  protected:
1025 
1026  //char description1[256] ; ///< String describing the model
1027  // char description2[256] ; ///< String describing the model
1033  double rHor ;
1034 
1035  // Constructors - Destructor
1036  // -------------------------
1037  public:
1045  ScalarBH(Map& mp_i, const char* file_name) ;
1046 
1047  ScalarBH(const ScalarBH& ) ;
1048 
1055  ScalarBH(Map& mp_i, FILE* fich) ;
1056 
1057  virtual ~ScalarBH() ;
1058 
1059  // Memory management
1060  // -----------------
1061  protected:
1063  virtual void del_deriv() const ;
1064 
1066  virtual void set_der_0x0() const ;
1067 
1068  // Mutators / assignment
1069  // ---------------------
1070  public:
1072  void operator=(const ScalarBH& ) ;
1073 
1074  // Accessors
1075  // ---------
1076  public:
1078  const Scalar& get_ff0() const {return ff0; } ;
1079  const Scalar& get_ff1() const {return ff1; } ;
1080  const Scalar& get_ff2() const {return ff2; } ;
1081  const Scalar& get_ww() const {return ww; } ;
1082  const Scalar& get_sfield() const {return sfield; } ;
1083  const double get_rHor() const {return rHor; } ;
1084 
1085  // Outputs
1086  // -------
1087  public:
1088  virtual void sauve(FILE* ) const ;
1089 
1090  protected:
1092  virtual ostream& operator>>(ostream& ) const ;
1093 
1094  // Global quantities
1095  // -----------------
1096  public:
1097 
1098 
1099  // Computational routines
1100  // ----------------------
1101  public:
1102  virtual void update_metric();
1103  };
1104 
1105 
1106  //------------------------//
1107  // class HiggsMonopole //
1108  //------------------------//
1109 
1116  class HiggsMonopole : public Compobj {
1117 
1118  // Data :
1119  // -----
1120  protected:
1121 
1122  char description1[256] ;
1123  char description2[256] ;
1124 
1126 
1128 
1130 
1131  // Derived data :
1132  // ------------
1133  protected:
1134 
1135 
1136  // Constructors - Destructor
1137  // -------------------------
1138  public:
1145  HiggsMonopole(Map& mp_i, const char* file_name) ;
1146 
1147  HiggsMonopole(const HiggsMonopole& ) ;
1148 
1149  virtual ~HiggsMonopole() ;
1150 
1151  // Memory management
1152  // -----------------
1153  protected:
1155  // virtual void del_deriv() const ;
1156 
1158  // virtual void set_der_0x0() const ;
1159 
1160  // Mutators / assignment
1161  // ---------------------
1162  public:
1164  // void operator=(const AltBH_QI& ) ;
1165 
1166  // Accessors
1167  // ---------
1168  public:
1169 
1171  const Scalar& get_higgs() const {return hh;} ;
1172 
1174  const Scalar& get_grr() const {return grr;} ;
1175 
1177  const Scalar& get_press() const {return press;} ;
1178 
1179  // Outputs
1180  // -------
1181  public:
1182  // virtual void sauve(FILE* ) const ; ///< Save in a file
1183 
1184  protected:
1186  virtual ostream& operator>>(ostream& ) const ;
1187 
1188  // Global quantities
1189  // -----------------
1190  public:
1191 
1192 
1193  // Computational routines
1194  // ----------------------
1195  public:
1196 
1198  // virtual void extrinsic_curvature() ;
1199 
1200  };
1201 
1202 
1203 
1204 
1205 }
1206 #endif
Alternative black hole spacetime in Quasi-Isotropic coordinates (under development).
Definition: compobj.h:917
char description2[256]
String describing the model.
Definition: compobj.h:924
void operator=(const AltBH_QI &)
Assignment to another AltBH_QI.
Definition: altBH_QI.C:226
AltBH_QI(Map &mp_i, const char *file_name, double a_spin_i)
Standard constructor.
Definition: altBH_QI.C:69
virtual ~AltBH_QI()
Destructor.
Definition: altBH_QI.C:196
char description1[256]
String describing the model.
Definition: compobj.h:923
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: altBH_QI.C:207
Scalar krphi
K_{(r)(phi)} read in the file.
Definition: compobj.h:927
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: altBH_QI.C:248
double a_spin
Spin parameter of the model.
Definition: compobj.h:925
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: altBH_QI.C:216
const Scalar & get_krphi() const
Returns K_{(r)(phi)}/sin(theta).
Definition: compobj.h:979
virtual void extrinsic_curvature()
Computation of the extrinsic curvature.
Definition: altBH_QI.C:269
virtual void sauve(FILE *) const
Save in a file.
Definition: altBH_QI.C:240
Base class for axisymmetric stationary compact objects in Quasi-Isotropic coordinates (under developm...
Definition: compobj.h:280
const Scalar & get_ak_car() const
Returns the scalar .
Definition: compobj.h:404
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: compobj_QI.C:178
double * p_lspec_isco
Specific angular momentum of a particle at the ISCO.
Definition: compobj.h:327
virtual double f_isco(int lmin) const
Orbital frequency at the innermost stable circular orbit (ISCO).
virtual void extrinsic_curvature()
Computes the extrinsic curvature and ak_car from nphi , nn and b_car .
Definition: compobj_QI.C:339
double * p_r_mb
Coordinate r of the marginally bound orbit.
Definition: compobj.h:328
void operator=(const Compobj_QI &)
Assignment to another Compobj_QI.
Definition: compobj_QI.C:195
double * p_r_isco
Coordinate r of the ISCO.
Definition: compobj.h:322
virtual void sauve(FILE *) const
Save in a file.
Definition: compobj_QI.C:216
const Scalar & get_nphi() const
Returns the metric coefficient .
Definition: compobj.h:384
void gyoto_data(const char *file_name) const
Save in a file for GYOTO.
Definition: compobj_QI.C:232
virtual double lspec_isco(int lmin) const
Angular momentum of a particle at the ISCO.
const Scalar & get_b_car() const
Returns the square of the metric factor B.
Definition: compobj.h:381
virtual double angu_mom() const
Angular momentum.
virtual double r_isco(int lmin, ostream *ost=0x0) const
Coordinate r of the innermost stable circular orbit (ISCO).
Scalar nphi
Metric coefficient .
Definition: compobj.h:296
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: compobj_QI.C:163
Compobj_QI(Map &map_i)
Standard constructor.
Definition: compobj_QI.C:86
const Scalar & get_bbb() const
Returns the metric factor B.
Definition: compobj.h:375
Scalar ak_car
Scalar .
Definition: compobj.h:315
virtual ~Compobj_QI()
Destructor.
Definition: compobj_QI.C:152
virtual void update_metric()
Updates the 3-metric from A and B and the shift vector from .
Definition: compobj_QI.C:302
virtual double espec_isco(int lmin) const
Energy of a particle at the ISCO.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: compobj_QI.C:264
Scalar b_car
Square of the metric factor B.
Definition: compobj.h:293
Scalar bbb
Metric factor B.
Definition: compobj.h:290
double * p_espec_isco
Specific energy of a particle at the ISCO.
Definition: compobj.h:325
Scalar a_car
Square of the metric factor A.
Definition: compobj.h:287
double * p_angu_mom
Angular momentum.
Definition: compobj.h:321
double * p_f_isco
Orbital frequency of the ISCO.
Definition: compobj.h:323
const Scalar & get_a_car() const
Returns the square of the metric factor A.
Definition: compobj.h:378
virtual double r_mb(int lmin, ostream *ost=0x0) const
Coordinate r of the marginally bound circular orbit (R_mb).
Base class for stationary compact objects (under development).
Definition: compobj.h:126
virtual double adm_mass() const
ADM mass (computed as a surface integral at spatial infinity)
Definition: compobj.C:310
const Map & get_mp() const
Returns the mapping.
Definition: compobj.h:207
Sym_tensor kk
Extrinsic curvature tensor
Definition: compobj.h:153
const Metric & get_gamma() const
Returns the 3-metric .
Definition: compobj.h:216
Vector mom_euler
Total 3-momentum density in the Eulerian frame.
Definition: compobj.h:147
const Sym_tensor & get_stress_euler() const
Returns the stress tensor with respect to the Eulerian observer.
Definition: compobj.h:225
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: compobj.C:155
Sym_tensor stress_euler
Stress tensor with respect to the Eulerian observer.
Definition: compobj.h:150
Scalar ener_euler
Total energy density E in the Eulerian frame.
Definition: compobj.h:144
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: compobj.C:163
const Vector & get_beta() const
Returns the shift vector .
Definition: compobj.h:213
void operator=(const Compobj &)
Assignment to another Compobj.
Definition: compobj.C:175
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: compobj.C:239
Metric gamma
3-metric
Definition: compobj.h:141
Compobj(Map &map_i)
Standard constructor.
Definition: compobj.C:82
const Vector & get_mom_euler() const
Returns the total 3-momentum density in the Eulerian frame.
Definition: compobj.h:222
virtual void sauve(FILE *) const
Save in a file.
Definition: compobj.C:196
Scalar nn
Lapse function N .
Definition: compobj.h:135
const Sym_tensor & get_kk() const
Returns the extrinsic curvature tensor .
Definition: compobj.h:228
Vector beta
Shift vector .
Definition: compobj.h:138
const Scalar & get_ener_euler() const
Returns the total energy density E in the Eulerian frame.
Definition: compobj.h:219
void gyoto_data(const char *file_name) const
Save in a file for GYOTO.
Definition: compobj.C:210
const Scalar & get_nn() const
Returns the lapse function N .
Definition: compobj.h:210
friend ostream & operator<<(ostream &, const Compobj &)
Display.
Definition: compobj.C:233
double * p_adm_mass
ADM mass.
Definition: compobj.h:158
Map & set_mp()
Read/write of the mapping.
Definition: compobj.h:200
virtual void extrinsic_curvature()
Computation of the extrinsic curvature.
Definition: compobj.C:290
virtual ~Compobj()
Destructor.
Definition: compobj.C:144
Map & mp
Mapping describing the coordinate system (r,theta,phi)
Definition: compobj.h:132
Higgs monopole (under development).
Definition: compobj.h:1116
HiggsMonopole(Map &mp_i, const char *file_name)
Standard constructor.
Scalar hh
Higgs scalar field.
Definition: compobj.h:1125
virtual ~HiggsMonopole()
Destructor.
const Scalar & get_grr() const
Returns the metric coefficient g_rr.
Definition: compobj.h:1174
Scalar press
Fluid pressure.
Definition: compobj.h:1129
char description2[256]
String describing the model.
Definition: compobj.h:1123
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
char description1[256]
String describing the model.
Definition: compobj.h:1122
Scalar grr
Metric coefficient g_rr.
Definition: compobj.h:1127
const Scalar & get_press() const
Returns the fluid pressure.
Definition: compobj.h:1177
const Scalar & get_higgs() const
Deletes all the derived quantities.
Definition: compobj.h:1171
Kerr spacetime in Quasi-Isotropic coordinates (under development).
Definition: compobj.h:816
virtual ~Kerr_QI()
Destructor.
Definition: kerr_QI.C:148
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: kerr_QI.C:159
double mm
mass parameter
Definition: compobj.h:824
Kerr_QI(Map &mp_i, double mass, double a_over_m)
Standard constructor.
Definition: kerr_QI.C:65
void operator=(const Kerr_QI &)
Assignment to another Kerr_QI.
Definition: kerr_QI.C:178
double aa
angular momentum parameter
Definition: compobj.h:828
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: kerr_QI.C:168
virtual void sauve(FILE *) const
Save in a file.
Definition: kerr_QI.C:195
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: kerr_QI.C:203
Base class for coordinate mappings.
Definition: map.h:670
Metric for tensor calculation.
Definition: metric.h:90
Black hole with scalar hair spacetime (under development).
Definition: compobj.h:1020
Scalar ff0
Metric field F_0 of Herdeiro & Radu (2015)
Definition: compobj.h:1028
Scalar ff2
Metric field F_2 of Herdeiro & Radu (2015)
Definition: compobj.h:1030
ScalarBH(Map &mp_i, const char *file_name)
Standard constructor.
Definition: scalarBH.C:72
Scalar ff1
Metric field F_1 of Herdeiro & Radu (2015)
Definition: compobj.h:1029
double rHor
Event horizon coordinate radius.
Definition: compobj.h:1033
Scalar sfield
Scalar field (modulus of Phi)
Definition: compobj.h:1032
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: scalarBH.C:409
void operator=(const ScalarBH &)
Assignment to another ScalarBH.
Definition: scalarBH.C:419
const Scalar & get_ff0() const
Returns f0.
Definition: compobj.h:1078
virtual ~ScalarBH()
Destructor.
Definition: scalarBH.C:389
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: scalarBH.C:441
virtual void sauve(FILE *) const
Save in a file.
Definition: scalarBH.C:433
Scalar ww
Metric field W of Herdeiro & Radu (2015)
Definition: compobj.h:1031
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: scalarBH.C:400
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Base class for axisymmetric stationary compact stars in Quasi-Isotropic coordinates (under developmen...
Definition: compobj.h:487
void fait_shift()
Computes shift from w_shift and khi_shift according to Shibata's prescription [Prog.
Definition: star_QI.C:340
Vector w_shift
Vector used in the decomposition of shift , following Shibata's prescription [Prog.
Definition: compobj.h:529
const Scalar & get_logn() const
Returns the logarithm of the lapse N.
Definition: compobj.h:636
virtual double grv2() const
Error on the virial identity GRV2.
double * p_grv2
Error on the virial identity GRV2.
Definition: compobj.h:585
Scalar logn
Logarithm of the lapse N .
Definition: compobj.h:495
double * p_mass_g
Gravitational mass (ADM mass as a volume integral)
Definition: compobj.h:588
Scalar nuq
Part of the Metric potential = logn generated by the quadratic terms.
Definition: compobj.h:510
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition: compobj.h:569
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
Star_QI(Map &mp_i)
Standard constructor.
Definition: star_QI.C:68
virtual void sauve(FILE *) const
Save in a file.
Definition: star_QI.C:279
virtual double mom_quad() const
Quadrupole moment.
virtual ~Star_QI()
Destructor.
Definition: star_QI.C:211
Scalar nuf
Part of the Metric potential = logn generated by the matter terms.
Definition: compobj.h:505
Scalar ssjm1_nuq
Effective source at the previous step for the resolution of the Poisson equation for nuq by means of ...
Definition: compobj.h:551
void update_metric()
Computes metric coefficients from known potentials.
Definition: star_QI.C:410
const Scalar & get_nuq() const
Returns the Part of the Metric potential = logn generated by the quadratic terms.
Definition: compobj.h:652
double * p_grv3
Error on the virial identity GRV3.
Definition: compobj.h:586
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: star_QI.C:300
const Scalar & get_tnphi() const
Returns the component of the shift vector.
Definition: compobj.h:642
Vector ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for .
Definition: compobj.h:578
void fait_nphi()
Computes tnphi and nphi from the Cartesian components of the shift, stored in shift .
Definition: star_QI.C:386
static double lambda_grv2(const Scalar &sou_m, const Scalar &sou_q)
Computes the coefficient which ensures that the GRV2 virial identity is satisfied.
virtual double angu_mom() const
Angular momentum.
Scalar khi_shift
Scalar used in the decomposition of shift , following Shibata's prescription [Prog.
Definition: compobj.h:539
const Scalar & get_tggg() const
Returns the Metric potential .
Definition: compobj.h:658
const Scalar & get_khi_shift() const
Returns the scalar used in the decomposition of shift following Shibata's prescription [Prog.
Definition: compobj.h:686
Scalar ssjm1_nuf
Effective source at the previous step for the resolution of the Poisson equation for nuf by means of ...
Definition: compobj.h:545
const Vector & get_w_shift() const
Returns the vector used in the decomposition of shift , following Shibata's prescription [Prog.
Definition: compobj.h:672
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_QI.C:222
virtual double mass_g() const
Gravitational mass.
void operator=(const Star_QI &)
Assignment to another Star_QI.
Definition: star_QI.C:250
const Scalar & get_dzeta() const
Returns the Metric potential .
Definition: compobj.h:655
const Scalar & get_nuf() const
Returns the part of the Metric potential = logn generated by the matter terms.
Definition: compobj.h:647
Scalar tggg
Metric potential .
Definition: compobj.h:516
Scalar tnphi
Component of the shift vector.
Definition: compobj.h:500
Scalar ssjm1_tggg
Effective source at the previous step for the resolution of the Poisson equation for tggg .
Definition: compobj.h:561
Scalar ssjm1_dzeta
Effective source at the previous step for the resolution of the Poisson equation for dzeta .
Definition: compobj.h:556
double * p_mom_quad
Quadrupole moment
Definition: compobj.h:587
Scalar dzeta
Metric potential .
Definition: compobj.h:513
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star_QI.C:235
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:223
Tensor field of valence 1.
Definition: vector.h:188
Lorene prototypes.
Definition: app_hor.h:64