LORENE
star_bhns_upmetr.C
1/*
2 * Method of class Star_bhns to compute metric quantities from
3 * the companion black-hole and total metric quantities
4 *
5 * (see file star_bhns.h for documentation).
6 *
7 */
8
9/*
10 * Copyright (c) 2005,2007 Keisuke Taniguchi
11 *
12 * This file is part of LORENE.
13 *
14 * LORENE is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU General Public License version 2
16 * as published by the Free Software Foundation.
17 *
18 * LORENE is distributed in the hope that it will be useful,
19 * but WITHOUT ANY WARRANTY; without even the implied warranty of
20 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 * GNU General Public License for more details.
22 *
23 * You should have received a copy of the GNU General Public License
24 * along with LORENE; if not, write to the Free Software
25 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26 *
27 */
28
29char star_bhns_upmetr_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns_upmetr.C,v 1.3 2014/10/13 08:53:41 j_novak Exp $" ;
30
31/*
32 * $Id: star_bhns_upmetr.C,v 1.3 2014/10/13 08:53:41 j_novak Exp $
33 * $Log: star_bhns_upmetr.C,v $
34 * Revision 1.3 2014/10/13 08:53:41 j_novak
35 * Lorene classes and functions now belong to the namespace Lorene.
36 *
37 * Revision 1.2 2008/05/15 19:17:39 k_taniguchi
38 * Change of some parameters.
39 *
40 * Revision 1.1 2007/06/22 01:32:37 k_taniguchi
41 * *** empty log message ***
42 *
43 *
44 * $Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns_upmetr.C,v 1.3 2014/10/13 08:53:41 j_novak Exp $
45 *
46 */
47
48// C++ headers
49//#include <>
50
51// C headers
52//#include <>
53
54// Lorene headers
55#include "star_bhns.h"
56#include "hole_bhns.h"
57#include "utilitaires.h"
58
59namespace Lorene {
61 const Star_bhns& star_prev,
62 double relax_met_comp) {
63
64 //-----------------------------------------------------
65 // Computation of quantities coming from the companion
66 //-----------------------------------------------------
67
68 const Map& mp_bh (hole.get_mp()) ;
69
70 // Lapconf function
71 // ----------------
72
73 if ( (hole.get_lapconf_auto()).get_etat() == ETATZERO ) {
75 }
76 else {
80 }
81
82 // Shift vector
83 // ------------
84
85 if ( (hole.get_shift_auto())(2).get_etat() == ETATZERO ) {
86 assert( (hole.get_shift_auto())(1).get_etat() == ETATZERO ) ;
87 assert( (hole.get_shift_auto())(3).get_etat() == ETATZERO ) ;
88
90 }
91 else {
94
96 comp_shift.change_triad(mp_bh.get_bvect_cart()) ;
97 comp_shift.change_triad(mp.get_bvect_cart()) ;
98
99 assert( *(shift_comp.get_triad()) == *(comp_shift.get_triad()) ) ;
100
101 (shift_comp.set(1)).import( comp_shift(1) ) ;
102 (shift_comp.set(2)).import( comp_shift(2) ) ;
103 (shift_comp.set(3)).import( comp_shift(3) ) ;
104
106 }
107
108 // Conformal factor
109 // ----------------
110
111 if ( (hole.get_confo_auto()).get_etat() == ETATZERO ) {
113 }
114 else {
118 }
119
120 //----------------------------------------------------
121 // Relaxation on lapconf_comp, shift_comp, confo_comp
122 //----------------------------------------------------
123
124 double relax_met_comp_jm1 = 1. - relax_met_comp ;
125
127 + relax_met_comp_jm1 * (star_prev.lapconf_comp) ;
128
130 + relax_met_comp_jm1 * (star_prev.shift_comp) ;
131
133 + relax_met_comp_jm1 * (star_prev.confo_comp) ;
134
135 //-------------------------
136 // Total metric quantities
137 //-------------------------
138
141
144
147
148 psi4 = pow(confo_tot, 4.) ;
150
153
156
157 //---------------------------------
158 // Derivative of metric quantities
159 //---------------------------------
160
162 for (int i=1; i<=3; i++)
164
166
168 for (int i=1; i<=3; i++) {
169 for (int j=1; j<=3; j++) {
170 d_shift_auto.set(i,j) = shift_auto(j).deriv(i) ;
171 }
172 }
173
175
177 for (int i=1; i<=3; i++)
179
181
182 // ... extrinsic curvature (taij_auto and taij_quad_auto)
183 // extr_curv_bhns() ;
184
185 // The derived quantities are obsolete
186 // -----------------------------------
187
188 del_deriv() ;
189
190}
191}
const Map & get_mp() const
Returns the mapping.
Definition blackhole.h:213
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Class for black holes in black hole-neutron star binaries.
Definition hole_bhns.h:65
const Scalar & get_lapconf_auto() const
Returns the part of the lapconf function generated by the black hole.
Definition hole_bhns.h:370
const Scalar & get_confo_auto() const
Returns the part of the conformal factor generated by the black hole.
Definition hole_bhns.h:439
const Vector & get_shift_auto() const
Returns the part of the shift vector generated by the black hole.
Definition hole_bhns.h:408
Base class for coordinate mappings.
Definition map.h:670
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
const Scalar & deriv(int i) const
Returns of *this , where .
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition scalar.C:353
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition scalar.C:784
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
void import(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
Class for stars in black hole-neutron star binaries.
Definition star_bhns.h:59
Vector shift_tot
Total shift vector.
Definition star_bhns.h:144
Scalar lapconf_auto
Lapconf function generated by the star.
Definition star_bhns.h:113
Scalar confo_comp
Conformal factor generated by the companion black hole.
Definition star_bhns.h:160
Tensor d_shift_auto
Derivative of the shift vector generated by the star .
Definition star_bhns.h:149
Scalar lapse_tot
Total lapse function.
Definition star_bhns.h:125
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star_bhns.C:341
Scalar confo_tot
Total conformal factor.
Definition star_bhns.h:163
Scalar psi4
Fourth power of the total conformal factor.
Definition star_bhns.h:176
Vector d_lapconf_auto
Derivative of the lapconf function generated by the star .
Definition star_bhns.h:130
Vector d_confo_auto
Derivative of the conformal factor generated by the star .
Definition star_bhns.h:168
Vector shift_comp
Shift vector generated by the companion black hole.
Definition star_bhns.h:141
Scalar lapconf_tot
Total lapconf function.
Definition star_bhns.h:119
Scalar lapse_auto
Lapse function generated by the "star".
Definition star_bhns.h:122
Scalar lapconf_comp
Lapconf function generated by the companion black hole.
Definition star_bhns.h:116
void update_metric_bhns(const Hole_bhns &hole, const Star_bhns &star_prev, double relax)
Computes metric coefficients from known potentials with relaxation when the companion is a black hole...
Vector shift_auto
Shift vector generated by the star.
Definition star_bhns.h:138
Scalar confo_auto
Conformal factor generated by the star.
Definition star_bhns.h:157
Map & mp
Mapping associated with the star.
Definition star.h:180
Tensor field of valence 1.
Definition vector.h:188
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition vector.C:316
Scalar & set(int)
Read/write access to a component.
Definition vector.C:296
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition tensor.C:519
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition tensor.h:866
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition tensor.C:481
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition tensor.C:497
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition tensor.C:654
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition tensor.C:926
Lorene prototypes.
Definition app_hor.h:64