star_QI.C
/*
* Methods of the class Star_QI
*
* (see file compobj.h for documentation).
*
*/
/*
* Copyright (c) 2012 Claire Some, Eric Gourgoulhon
*
* This file is part of LORENE.
*
* LORENE is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License version 2
* as published by the Free Software Foundation.
*
* LORENE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with LORENE; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
*/
char star_QI_C[] = "$Header: /cvsroot/Lorene/C++/Source/Compobj/star_QI.C,v 1.5 2014/10/13 08:52:49 j_novak Exp $" ;
/*
* $Id: star_QI.C,v 1.5 2014/10/13 08:52:49 j_novak Exp $
* $Log: star_QI.C,v $
* Revision 1.5 2014/10/13 08:52:49 j_novak
* Lorene classes and functions now belong to the namespace Lorene.
*
* Revision 1.4 2012/12/03 15:27:30 c_some
* Small changes
*
* Revision 1.3 2012/11/22 16:04:51 c_some
* Minor modifications
*
* Revision 1.2 2012/11/21 14:55:27 c_some
* Added methods fait_shift and fait_nphi
*
* Revision 1.1 2012/11/20 16:29:50 c_some
* New class Star_QI
*
*
* $Header: /cvsroot/Lorene/C++/Source/Compobj/star_QI.C,v 1.5 2014/10/13 08:52:49 j_novak Exp $
*
*/
// C headers
#include <cassert>
// Lorene headers
#include "compobj.h"
#include "unites.h"
//--------------//
// Constructors //
//--------------//
// Standard constructor
// --------------------
namespace Lorene {
Star_QI::Star_QI(Map& mpi) :
Compobj_QI(mpi) ,
logn(mpi),
tnphi(mpi),
nuf(mpi),
nuq(mpi),
dzeta(mpi),
tggg(mpi),
w_shift(mpi, CON, mp.get_bvect_cart()),
khi_shift(mpi),
ssjm1_nuf(mpi),
ssjm1_nuq(mpi),
ssjm1_dzeta(mpi),
ssjm1_tggg(mpi),
ssjm1_khi(mpi),
ssjm1_wshift(mpi, CON, mp.get_bvect_cart())
{
// Pointers of derived quantities initialized to zero :
set_der_0x0() ;
// Initialization to a flat metric :
logn = 0 ;
tnphi = 0 ;
nuf = 0 ;
nuq = 0 ;
dzeta = 0 ;
tggg = 0 ;
w_shift.set_etat_zero() ;
khi_shift = 0 ;
beta.set_etat_zero() ;
ssjm1_nuf = 0 ;
ssjm1_nuq = 0 ;
ssjm1_dzeta = 0 ;
ssjm1_tggg = 0 ;
ssjm1_khi = 0 ;
ssjm1_wshift.set_etat_zero() ;
}
// Copy constructor
// --------------------
Star_QI::Star_QI(const Star_QI& st) :
Compobj_QI(st),
logn(st.logn),
tnphi(st.tnphi),
nuf(st.nuf),
nuq(st.nuq),
dzeta(st.dzeta),
tggg(st.tggg),
w_shift(st.w_shift),
khi_shift(st.khi_shift),
ssjm1_nuf(st.ssjm1_nuf),
ssjm1_nuq(st.ssjm1_nuq),
ssjm1_dzeta(st.ssjm1_dzeta),
ssjm1_tggg(st.ssjm1_tggg),
ssjm1_khi(st.ssjm1_khi),
ssjm1_wshift(st.ssjm1_wshift)
{
// Pointers of derived quantities initialized to zero :
set_der_0x0() ;
}
// Constructor from a file
// -----------------------
Star_QI::Star_QI(Map& mpi, FILE* fich) :
Compobj_QI(mpi) ,
logn(mpi),
tnphi(mpi),
nuf(mpi),
nuq(mpi),
dzeta(mpi),
tggg(mpi),
w_shift(mpi, CON, mp.get_bvect_cart()),
khi_shift(mpi),
ssjm1_nuf(mpi),
ssjm1_nuq(mpi),
ssjm1_dzeta(mpi),
ssjm1_tggg(mpi),
ssjm1_khi(mpi),
ssjm1_wshift(mpi, CON, mp.get_bvect_cart())
{
// Pointers of derived quantities initialized to zero :
set_der_0x0() ;
// Read of the saved fields:
// ------------------------
Scalar nuf_file(mp, *(mp.get_mg()), fich) ;
nuf = nuf_file ;
Scalar nuq_file(mp, *(mp.get_mg()), fich) ;
nuq = nuq_file ;
logn = nuf + nuq ; //## to be checked !
Scalar dzeta_file(mp, *(mp.get_mg()), fich) ;
dzeta = dzeta_file ;
Scalar tggg_file(mp, *(mp.get_mg()), fich) ;
tggg = tggg_file ;
Vector w_shift_file(mp, mp.get_bvect_cart(), fich) ;
w_shift = w_shift_file ;
Scalar khi_shift_file(mp, *(mp.get_mg()), fich) ;
khi_shift = khi_shift_file ;
fait_shift() ; // constructs shift from w_shift and khi_shift
fait_nphi() ; // constructs N^phi from (N^x,N^y,N^z)
Scalar ssjm1_nuf_file(mp, *(mp.get_mg()), fich) ;
ssjm1_nuf = ssjm1_nuf_file ;
Scalar ssjm1_nuq_file(mp, *(mp.get_mg()), fich) ;
ssjm1_nuq = ssjm1_nuq_file ;
Scalar ssjm1_dzeta_file(mp, *(mp.get_mg()), fich) ;
ssjm1_dzeta = ssjm1_dzeta_file ;
Scalar ssjm1_tggg_file(mp, *(mp.get_mg()), fich) ;
ssjm1_tggg = ssjm1_tggg_file ;
Scalar ssjm1_khi_file(mp, *(mp.get_mg()), fich) ;
ssjm1_khi = ssjm1_khi_file ;
Vector ssjm1_wshift_file(mp, mp.get_bvect_cart(), fich) ;
ssjm1_wshift = ssjm1_wshift_file ;
// Initialization of N, A^2, B^2, gamma_ij, tkij and ak_car
update_metric() ;
}
//------------//
// Destructor //
//------------//
Star_QI::~Star_QI(){
del_deriv() ;
}
//----------------------------------//
// Management of derived quantities //
//----------------------------------//
void Star_QI::del_deriv() const {
Compobj_QI::del_deriv() ;
if (p_grv2 != 0x0) delete p_grv2 ;
if (p_grv3 != 0x0) delete p_grv3 ;
if (p_mom_quad != 0x0) delete p_mom_quad ;
if (p_mass_g != 0x0) delete p_mass_g ;
Star_QI::set_der_0x0() ;
}
void Star_QI::set_der_0x0() const {
p_grv2 = 0x0 ;
p_grv3 = 0x0 ;
p_mom_quad = 0x0 ;
p_mass_g = 0x0 ;
}
//--------------//
// Assignment //
//--------------//
// Assignment to another Star_QI
// --------------------------------
void Star_QI::operator=(const Star_QI& st) {
// Assignment of quantities common to all the derived classes of Compobj_QI
Compobj_QI::operator=(st) ;
logn = st.logn ;
tnphi = st.tnphi ;
nuf = st.nuf ;
nuq = st.nuq ;
dzeta = st.dzeta ;
tggg = st.tggg ;
w_shift = st.w_shift ;
khi_shift = st.khi_shift ;
ssjm1_nuf = st.ssjm1_nuf ;
ssjm1_nuq = st.ssjm1_nuq ;
ssjm1_dzeta = st.ssjm1_dzeta ;
ssjm1_tggg = st.ssjm1_tggg ;
ssjm1_khi = st.ssjm1_khi ;
ssjm1_wshift = st.ssjm1_wshift ;
del_deriv() ; // Deletes all derived quantities
}
//--------------//
// Outputs //
//--------------//
// Save in a file
// --------------
void Star_QI::sauve(FILE* fich) const {
nuf.sauve(fich) ;
nuq.sauve(fich) ;
dzeta.sauve(fich) ;
tggg.sauve(fich) ;
w_shift.sauve(fich) ;
khi_shift.sauve(fich) ;
ssjm1_nuf.sauve(fich) ;
ssjm1_nuq.sauve(fich) ;
ssjm1_dzeta.sauve(fich) ;
ssjm1_tggg.sauve(fich) ;
ssjm1_khi.sauve(fich) ;
ssjm1_wshift.sauve(fich) ;
}
// Printing
// --------
ostream& Star_QI::operator>>(ostream& ost) const {
using namespace Unites ;
Compobj_QI::operator>>(ost) ;
ost << endl << "Axisymmetric stationary compact star in quasi-isotropic coordinates (class Star_QI) " << endl ;
ost << "Central values of various fields : " << endl ;
ost << "-------------------------------- " << endl ;
ost << " ln(N) : " << logn.val_grid_point(0,0,0,0) << endl ;
ost << " nuf : " << nuf.val_grid_point(0,0,0,0) << endl ;
ost << " nuq : " << nuq.val_grid_point(0,0,0,0) << endl ;
ost << " zeta = ln(AN): " << dzeta.val_grid_point(0,0,0,0) << endl << endl ;
ost << "Error on the virial identity GRV2 : " << grv2() << endl ;
ost << "Error on the virial identity GRV3 : " << grv3(&ost) << endl ;
double mom_quad_38si = mom_quad() * rho_unit * (pow(r_unit, double(5.))
/ double(1.e38) ) ;
ost << "Quadrupole moment Q : " << mom_quad_38si << " 10^38 kg m^2"
<< endl ;
ost << "c^4 Q / (G^2 M^3) : "
<< mom_quad() / ( pow(qpig/(4*M_PI), 2.) * pow(mass_g(), 3.) )
<< endl ;
ost << "Total angular momentum J : "
<< angu_mom()/( qpig / (4* M_PI) * msol*msol) << " G M_sol^2 / c"
<< endl ;
ost << "c J / (G M^2) : "
<< angu_mom()/( qpig / (4* M_PI) * pow(mass_g(), 2.) ) << endl ;
return ost ;
}
//-------------------------//
// Computational methods //
//-------------------------//
void Star_QI::fait_shift() {
Vector d_khi = khi_shift.derive_con( mp.flat_met_cart() ) ;
d_khi.dec_dzpuis(2) ; // divide by r^2 in the external compactified domain
// x_k dW^k/dx_i
Scalar xx(mp) ;
Scalar yy(mp) ;
Scalar zz(mp) ;
Scalar sintcosp(mp) ;
Scalar sintsinp(mp) ;
Scalar cost(mp) ;
xx = mp.x ;
yy = mp.y ;
zz = mp.z ;
sintcosp = mp.sint * mp.cosp ;
sintsinp = mp.sint * mp.sinp ;
cost = mp.cost ;
int nz = mp.get_mg()->get_nzone() ;
Vector xk(mp, COV, mp.get_bvect_cart()) ;
xk.set(1) = xx ;
xk.set(1).set_domain(nz-1) = sintcosp.domain(nz-1) ;
xk.set(1).set_dzpuis(-1) ;
xk.set(2) = yy ;
xk.set(2).set_domain(nz-1) = sintsinp.domain(nz-1) ;
xk.set(2).set_dzpuis(-1) ;
xk.set(3) = zz ;
xk.set(3).set_domain(nz-1) = cost.domain(nz-1) ;
xk.set(3).set_dzpuis(-1) ;
xk.std_spectral_base() ;
Tensor d_w = w_shift.derive_con( mp.flat_met_cart() ) ;
Vector x_d_w = contract(xk, 0, d_w, 0) ;
x_d_w.dec_dzpuis() ;
double lambda = double(1) / double(3) ;
beta = - (lambda+2)/2./(lambda+1) * w_shift
+ (lambda/2./(lambda+1)) * (d_khi + x_d_w) ;
}
void Star_QI::fait_nphi() {
if ( (beta(1).get_etat() == ETATZERO) && (beta(2).get_etat() == ETATZERO) ) {
tnphi = 0 ;
nphi = 0 ;
return ;
}
// Computation of tnphi
// --------------------
mp.comp_p_from_cartesian( -beta(1), -beta(2), tnphi ) ;
// Computation of nphi
// -------------------
nphi = tnphi ;
nphi.div_rsint() ;
}
// Updates the 3-metric and the shift
void Star_QI::update_metric() {
// Lapse function N
// ----------------
nn = exp( logn ) ;
nn.std_spectral_base() ; // set the bases for spectral expansions
// Metric factor A^2
// -----------------
a_car = exp( 2*( dzeta - logn ) ) ;
a_car.std_spectral_base() ; // set the bases for spectral expansions
// Metric factor B
// ---------------
Scalar tmp = tggg ;
tmp.div_rsint() ; //... Division of tG by r sin(theta)
bbb = (1 + tmp) / nn ;
bbb.std_spectral_base() ; // set the bases for spectral expansions
b_car = bbb * bbb ;
Compobj_QI::update_metric() ; // updates gamma_{ij}
// Tensor B^{-2} K_{ij} and Scalar A^2 K_{ij} K^{ij}
// -------------------------------------------------
//## extrinsic_curvature() ; // should be done by Compobj_QI::update_metric()
// The derived quantities are no longer up to date :
// -----------------------------------------------
del_deriv() ;
}
}
Lorene::Compobj_QI::operator=
void operator=(const Compobj_QI &)
Assignment to another Compobj_QI.
Definition compobj_QI.C:195
Lorene::Compobj_QI::nphi
Scalar nphi
Metric coefficient .
Definition compobj.h:296
Lorene::Compobj_QI::del_deriv
virtual void del_deriv() const
Deletes all the derived quantities.
Definition compobj_QI.C:163
Lorene::Compobj_QI::Compobj_QI
Compobj_QI(Map &map_i)
Standard constructor.
Definition compobj_QI.C:86
Lorene::Compobj_QI::update_metric
virtual void update_metric()
Updates the 3-metric from A and B and the shift vector from .
Definition compobj_QI.C:302
Lorene::Compobj_QI::operator>>
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition compobj_QI.C:264
Lorene::Compobj_QI::b_car
Scalar b_car
Square of the metric factor B.
Definition compobj.h:293
Lorene::Compobj_QI::bbb
Scalar bbb
Metric factor B.
Definition compobj.h:290
Lorene::Compobj_QI::a_car
Scalar a_car
Square of the metric factor A.
Definition compobj.h:287
Lorene::Compobj::nn
Scalar nn
Lapse function N .
Definition compobj.h:135
Lorene::Compobj::beta
Vector beta
Shift vector .
Definition compobj.h:138
Lorene::Compobj::mp
Map & mp
Mapping describing the coordinate system (r,theta,phi).
Definition compobj.h:132
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Lorene::Scalar::set_domain
Tbl & set_domain(int l)
Read/write of the value in a given domain.
Definition scalar.h:615
Lorene::Scalar::set_dzpuis
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition scalar.C:808
Lorene::Scalar::div_rsint
void div_rsint()
Division by everywhere; dzpuis is not changed.
Definition scalar_r_manip.C:348
Lorene::Scalar::domain
const Tbl & domain(int l) const
Read-only of the value in a given domain.
Definition scalar.h:625
Lorene::Star_QI::fait_shift
void fait_shift()
Computes shift from w_shift and khi_shift according to Shibata's prescription [Prog.
Definition star_QI.C:340
Lorene::Star_QI::w_shift
Vector w_shift
Vector used in the decomposition of shift , following Shibata's prescription [Prog.
Definition compobj.h:529
Lorene::Star_QI::grv2
virtual double grv2() const
Error on the virial identity GRV2.
Definition star_QI_global.C:123
Lorene::Star_QI::p_grv2
double * p_grv2
Error on the virial identity GRV2.
Definition compobj.h:585
Lorene::Star_QI::logn
Scalar logn
Logarithm of the lapse N .
Definition compobj.h:495
Lorene::Star_QI::p_mass_g
double * p_mass_g
Gravitational mass (ADM mass as a volume integral).
Definition compobj.h:588
Lorene::Star_QI::nuq
Scalar nuq
Part of the Metric potential = logn generated by the quadratic terms.
Definition compobj.h:510
Lorene::Star_QI::ssjm1_khi
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
Lorene::Star_QI::grv3
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
Definition star_QI_global.C:148
Lorene::Star_QI::Star_QI
Star_QI(Map &mp_i)
Standard constructor.
Definition star_QI.C:68
Lorene::Star_QI::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition star_QI.C:279
Lorene::Star_QI::mom_quad
virtual double mom_quad() const
Quadrupole moment.
Definition star_QI_global.C:250
Lorene::Star_QI::~Star_QI
virtual ~Star_QI()
Destructor.
Definition star_QI.C:211
Lorene::Star_QI::nuf
Scalar nuf
Part of the Metric potential = logn generated by the matter terms.
Definition compobj.h:505
Lorene::Star_QI::ssjm1_nuq
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
Lorene::Star_QI::update_metric
void update_metric()
Computes metric coefficients from known potentials.
Definition star_QI.C:410
Lorene::Star_QI::p_grv3
double * p_grv3
Error on the virial identity GRV3.
Definition compobj.h:586
Lorene::Star_QI::operator>>
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition star_QI.C:300
Lorene::Star_QI::ssjm1_wshift
Vector ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for .
Definition compobj.h:578
Lorene::Star_QI::fait_nphi
void fait_nphi()
Computes tnphi and nphi from the Cartesian components of the shift, stored in shift .
Definition star_QI.C:386
Lorene::Star_QI::angu_mom
virtual double angu_mom() const
Angular momentum.
Definition star_QI_global.C:96
Lorene::Star_QI::khi_shift
Scalar khi_shift
Scalar used in the decomposition of shift , following Shibata's prescription [Prog.
Definition compobj.h:539
Lorene::Star_QI::ssjm1_nuf
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
Lorene::Star_QI::del_deriv
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star_QI.C:222
Lorene::Star_QI::mass_g
virtual double mass_g() const
Gravitational mass.
Definition star_QI_global.C:64
Lorene::Star_QI::operator=
void operator=(const Star_QI &)
Assignment to another Star_QI.
Definition star_QI.C:250
Lorene::Star_QI::tggg
Scalar tggg
Metric potential .
Definition compobj.h:516
Lorene::Star_QI::tnphi
Scalar tnphi
Component of the shift vector.
Definition compobj.h:500
Lorene::Star_QI::ssjm1_tggg
Scalar ssjm1_tggg
Effective source at the previous step for the resolution of the Poisson equation for tggg .
Definition compobj.h:561
Lorene::Star_QI::ssjm1_dzeta
Scalar ssjm1_dzeta
Effective source at the previous step for the resolution of the Poisson equation for dzeta .
Definition compobj.h:556
Lorene::Star_QI::p_mom_quad
double * p_mom_quad
Quadrupole moment.
Definition compobj.h:587
Lorene::Star_QI::dzeta
Scalar dzeta
Metric potential .
Definition compobj.h:513
Lorene::Star_QI::set_der_0x0
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition star_QI.C:235
Lorene::Tensor
Tensor handling.
Definition tensor.h:288
Lorene::Vector
Tensor field of valence 1.
Definition vector.h:188
Lorene::Vector::std_spectral_base
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition vector.C:316
Lorene::Vector::set
Scalar & set(int)
Read/write access to a component.
Definition vector.C:296
Lorene::exp
Cmp exp(const Cmp &)
Exponential.
Definition cmp_math.C:270
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Lorene::Tensor::dec_dzpuis
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition tensor.C:808
Lorene::contract
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Definition tenseur_operateur.C:279
Lorene
Lorene prototypes.
Definition app_hor.h:64
Lorene::get_bvect_cart
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition map.h:791
Lorene::cost
Coord cost
Definition map.h:722
Lorene::Map
Map(const Mg3d &)
Constructor from a multi-domain 3D grid.
Definition map.C:139
Unites
Standard units of space, time and mass.