LORENE
star_bin_upmetr_der.C
1 /*
2  * Methods Star_bin::update_metric_der_comp
3  *
4  * (see file star.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 2004 Francois Limousin
10  *
11  * This file is part of LORENE.
12  *
13  * LORENE is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License as published by
15  * the Free Software Foundation; either version 2 of the License, or
16  * (at your option) any later version.
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 
29 
30 
31 
32 /*
33  * $Id: star_bin_upmetr_der.C,v 1.18 2016/12/05 16:18:15 j_novak Exp $
34  * $Log: star_bin_upmetr_der.C,v $
35  * Revision 1.18 2016/12/05 16:18:15 j_novak
36  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
37  *
38  * Revision 1.17 2014/10/13 08:53:38 j_novak
39  * Lorene classes and functions now belong to the namespace Lorene.
40  *
41  * Revision 1.16 2014/10/06 15:13:17 j_novak
42  * Modified #include directives to use c++ syntax.
43  *
44  * Revision 1.15 2006/05/24 16:52:50 f_limousin
45  * New computation of tkij_comp
46  *
47  * Revision 1.14 2005/09/13 19:38:31 f_limousin
48  * Reintroduction of the resolution of the equations in cartesian coordinates.
49  *
50  * Revision 1.13 2005/02/24 16:26:46 f_limousin
51  * Add the name of the variable for the companion which is now used
52  * to compute dlogn_comp, dlnq_comp...
53  *
54  * Revision 1.12 2005/02/24 16:07:23 f_limousin
55  * Improve the computation of dlogn, dlnq and dlnpsi. We import the part
56  * coming from the companion and add the auto part.
57  *
58  * Revision 1.11 2005/02/18 13:14:18 j_novak
59  * Changing of malloc/free to new/delete + suppression of some unused variables
60  * (trying to avoid compilation warnings).
61  *
62  * Revision 1.10 2005/02/17 17:34:28 f_limousin
63  * Change the name of some quantities to be consistent with other classes
64  * (for instance nnn is changed to nn, shift to beta, beta to lnq...)
65  *
66  * Revision 1.9 2004/06/22 12:52:47 f_limousin
67  * Change qq, qq_auto and qq_comp to beta, beta_auto and beta_comp.
68  *
69  * Revision 1.8 2004/06/07 16:25:14 f_limousin
70  * Minor modif.
71  *
72  * Revision 1.7 2004/04/08 16:33:32 f_limousin
73  * The new variable is ln(Q) instead of Q=psi^2*N. It improves the
74  * convergence of the code.
75  *
76  * Revision 1.6 2004/03/23 10:00:09 f_limousin
77  * We now make the derivation with respect to the metric tilde
78  * instead of the flat metric for the computation of dshift_comp.
79  *
80  * Revision 1.5 2004/02/27 09:53:14 f_limousin
81  * Correction of an error on the computation of kcar_comp.
82  *
83  * Revision 1.4 2004/02/18 18:47:01 e_gourgoulhon
84  * divshift_comp now computed via Tensor::divergence, the
85  * method Tensor::scontract having disappeared.
86  *
87  * Revision 1.3 2004/01/20 15:20:23 f_limousin
88  * First version
89  *
90  *
91  * $Header: /cvsroot/Lorene/C++/Source/Star/star_bin_upmetr_der.C,v 1.18 2016/12/05 16:18:15 j_novak Exp $
92  *
93  */
94 
95 // C headers
96 #include <cmath>
97 
98 // Headers Lorene
99 #include "star.h"
100 #include "utilitaires.h"
101 #include "graphique.h"
102 
103 namespace Lorene {
104 void Star_bin::update_metric_der_comp(const Star_bin& comp, double om) {
105 
106  // Derivatives of metric coefficients
107  // ----------------------------------
108 
109  // dcov_logn
110  Vector temp ((comp.logn_auto).derive_cov(comp.get_flat())) ;
111  temp.dec_dzpuis(2) ;
112  temp.change_triad(temp.get_mp().get_bvect_cart()) ;
113  temp.change_triad(mp.get_bvect_cart()) ;
114  Base_val sauve_base1 (temp(1).get_spectral_va().get_base()) ;
115  Base_val sauve_base2 (temp(2).get_spectral_va().get_base()) ;
116  Base_val sauve_base3 (temp(3).get_spectral_va().get_base()) ;
117  assert ( *(temp.get_triad()) == *(dcov_logn.get_triad())) ;
118 
119  dcov_logn.set(1).import(temp(1)) ;
120  dcov_logn.set(2).import(temp(2)) ;
121  dcov_logn.set(3).import(temp(3)) ;
122 
123  dcov_logn.set(1).set_spectral_va().set_base(sauve_base1) ;
124  dcov_logn.set(2).set_spectral_va().set_base(sauve_base2) ;
125  dcov_logn.set(3).set_spectral_va().set_base(sauve_base3) ;
126  dcov_logn.inc_dzpuis(2) ;
127 
129 
130  // dcon_logn
131  Vector temp_con((comp.logn_auto).derive_con(comp.get_flat())) ;
132  temp_con.dec_dzpuis(2) ;
133  temp_con.change_triad(temp_con.get_mp().get_bvect_cart()) ;
134  temp_con.change_triad(mp.get_bvect_cart()) ;
135  sauve_base1 = temp_con(1).get_spectral_va().get_base() ;
136  sauve_base2 = temp_con(2).get_spectral_va().get_base() ;
137  sauve_base3 = temp_con(3).get_spectral_va().get_base() ;
138  assert ( *(temp_con.get_triad()) == *(dcon_logn.get_triad())) ;
139 
140  dcon_logn.set(1).import(temp_con(1)) ;
141  dcon_logn.set(2).import(temp_con(2)) ;
142  dcon_logn.set(3).import(temp_con(3)) ;
143 
144  dcon_logn.set(1).set_spectral_va().set_base(sauve_base1) ;
145  dcon_logn.set(2).set_spectral_va().set_base(sauve_base2) ;
146  dcon_logn.set(3).set_spectral_va().set_base(sauve_base3) ;
147  dcon_logn.inc_dzpuis(2) ;
148 
150 
151  // dcov_phi
152  temp = (0.5*(comp.lnq_auto-comp.logn_auto)).derive_cov(comp.get_flat()) ;
153  temp.dec_dzpuis(2) ;
154  temp.change_triad(temp.get_mp().get_bvect_cart()) ;
155  temp.change_triad(mp.get_bvect_cart()) ;
156  sauve_base1 = temp(1).get_spectral_va().get_base() ;
157  sauve_base2 = temp(2).get_spectral_va().get_base() ;
158  sauve_base3 = temp(3).get_spectral_va().get_base() ;
159  assert ( *(temp.get_triad()) == *(dcov_phi.get_triad())) ;
160 
161  dcov_phi.set(1).import(temp(1)) ;
162  dcov_phi.set(2).import(temp(2)) ;
163  dcov_phi.set(3).import(temp(3)) ;
164 
165  dcov_phi.set(1).set_spectral_va().set_base(sauve_base1) ;
166  dcov_phi.set(2).set_spectral_va().set_base(sauve_base2) ;
167  dcov_phi.set(3).set_spectral_va().set_base(sauve_base3) ;
168  dcov_phi.inc_dzpuis(2) ;
169 
170  dcov_phi = dcov_phi + 0.5*(lnq_auto - logn_auto).derive_cov(flat) ;
171 
172  // dcon_phi
173  temp_con = (0.5*(comp.lnq_auto-comp.logn_auto))
174  .derive_con(comp.get_flat()) ;
175  temp_con.dec_dzpuis(2) ;
176  temp_con.change_triad(temp_con.get_mp().get_bvect_cart()) ;
177  temp_con.change_triad(mp.get_bvect_cart()) ;
178  sauve_base1 = temp_con(1).get_spectral_va().get_base() ;
179  sauve_base2 = temp_con(2).get_spectral_va().get_base() ;
180  sauve_base3 = temp_con(3).get_spectral_va().get_base() ;
181  assert ( *(temp_con.get_triad()) == *(dcon_phi.get_triad())) ;
182 
183  dcon_phi.set(1).import(temp_con(1)) ;
184  dcon_phi.set(2).import(temp_con(2)) ;
185  dcon_phi.set(3).import(temp_con(3)) ;
186 
187  dcon_phi.set(1).set_spectral_va().set_base(sauve_base1) ;
188  dcon_phi.set(2).set_spectral_va().set_base(sauve_base2) ;
189  dcon_phi.set(3).set_spectral_va().set_base(sauve_base3) ;
190  dcon_phi.inc_dzpuis(2) ;
191 
192  dcon_phi = dcon_phi + 0.5*(lnq_auto - logn_auto).derive_con(flat) ;
193 
194 
195  // Construction of Omega d/dphi
196  // ----------------------------
197 
198  const Mg3d* mg = mp.get_mg() ;
199  int nz = mg->get_nzone() ; // total number of domains
200  Vector omdsdp (mp, CON, mp.get_bvect_cart()) ;
201  Scalar yya (mp) ;
202  yya = mp.ya ;
203  Scalar xxa (mp) ;
204  xxa = mp.xa ;
205 
206  if (fabs(mp.get_rot_phi()) < 1e-10){
207  omdsdp.set(1) = - om * yya ;
208  omdsdp.set(2) = om * xxa ;
209  omdsdp.set(3).annule_hard() ;
210  }
211  else{
212  omdsdp.set(1) = om * yya ;
213  omdsdp.set(2) = - om * xxa ;
214  omdsdp.set(3).annule_hard() ;
215  }
216 
217  omdsdp.set(1).set_spectral_va()
218  .set_base(*(mp.get_mg()->std_base_vect_cart()[0])) ;
219  omdsdp.set(2).set_spectral_va()
220  .set_base(*(mp.get_mg()->std_base_vect_cart()[1])) ;
221  omdsdp.set(3).set_spectral_va()
222  .set_base(*(mp.get_mg()->std_base_vect_cart()[2])) ;
223 
224  omdsdp.annule_domain(nz-1) ;
225 
226  // Computation of tkij_comp
227  // ------------------------
228 
229  // Gradient tilde (partial derivatives with respect to
230  // the cartesian coordinates of the mapping)
231  // D~^j beta^i
232 
233  const Tensor& dbeta_comp = beta_comp.derive_con(gtilde) ;
234 
235  // Trace of D~_j beta^i :
236  Scalar divbeta_comp = beta_comp.divergence(gtilde) ;
237 
238  // Computation of K^{ij}
239  // -------------------------
240 
241  for (int i=1; i<=3; i++)
242  for (int j=i; j<=3; j++) {
243 
244  tkij_comp.set(i, j) = dbeta_comp(i, j) + dbeta_comp(j, i) -
245  double(2) /double(3) * divbeta_comp * (gtilde.con())(i,j) ;
246  }
247 
248  // Addition (or not !) of u^{ij}
249  tkij_comp = tkij_comp - 0*hij_comp.derive_lie(omdsdp) ;
250 
251  tkij_comp = 0.5 * tkij_comp / nn ;
252 
253  /*
254  Sym_tensor aa_comp (mp, CON, mp.get_bvect_cart()) ;
255  aa_comp.set_etat_qcq() ;
256  Sym_tensor comp_aa (comp.get_tkij_auto()) ;
257  comp_aa.dec_dzpuis(2) ;
258  comp_aa.change_triad(mp.get_bvect_cart()) ;
259 
260 
261  assert(*(aa_comp.get_triad()) == *(comp_aa.get_triad())) ;
262  // importations :
263  for (int i=1 ; i<=3 ; i++)
264  for (int j=i ; j<=3 ; j++) {
265  aa_comp.set(i, j).import(comp_aa(i, j)) ;
266  aa_comp.set(i, j).set_spectral_va().set_base(comp_aa(i, j).
267  get_spectral_va().get_base()) ;
268  }
269  aa_comp.inc_dzpuis(2) ;
270 
271  for (int i=1 ; i<=3 ; i++)
272  for (int j=i ; j<=3 ; j++)
273  for (int l=0 ; l<=nz-2 ; l++)
274  tkij_comp.set(i,j).set_domain(l) = aa_comp(i,j).domain(l) ;
275  */
276 
277 
278  // Computation of kcar_comp
279  // ------------------------
280 
281  Sym_tensor tkij_auto_cov = tkij_auto.up_down(gtilde) ;
282 
283  kcar_comp = contract(tkij_auto_cov, 0, 1, tkij_comp, 0, 1,true) ;
284 
285  // The derived quantities are obsolete
286  // -----------------------------------
287 
288  del_deriv() ;
289 
290 }
291 
292 }
Sym_tensor hij_comp
Deviation of the inverse conformal metric from the inverse flat metric generated principally by the ...
Definition: star.h:594
Vector dcov_phi
Covariant derivative of the logarithm of the conformal factor.
Definition: star.h:555
Coord xa
Absolute x coordinate.
Definition: map.h:748
Base_val ** std_base_vect_cart() const
Returns the standard spectral bases for the Cartesian components of a vector.
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Definition: metric.C:293
Map & mp
Mapping associated with the star.
Definition: star.h:180
Lorene prototypes.
Definition: app_hor.h:67
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
Tensor up_down(const Metric &gam) const
Computes a new tensor by raising or lowering all the indices of *this .
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition: map.h:793
void set(const Map *mp, Mtbl *(*construct)(const Map *))
Semi-constructor from a mapping and a method.
Definition: coord.C:137
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Definition: scalar_deriv.C:402
Tensor field of valence 1.
Definition: vector.h:188
Scalar logn_auto
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
Definition: star.h:527
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:813
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:817
void annule_hard()
Sets the Scalar to zero in a hard way.
Definition: scalar.C:386
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tensor.h:879
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) ...
Definition: star.h:562
Scalar lnq_auto
Scalar field generated principally by the star.
Definition: star.h:543
const Metric & get_flat() const
Return the flat metric defined on the mapping (Spherical components with respect to the mapping of th...
Definition: star.h:859
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_bin.C:372
const Scalar & divergence(const Metric &) const
The divergence of this with respect to a Metric .
Definition: vector.C:387
Sym_tensor tkij_comp
Part of the extrinsic curvature tensor generated by beta_comp.
Definition: star.h:606
Tensor handling.
Definition: tensor.h:294
Vector dcov_logn
Covariant derivative of the total logarithm of the lapse.
Definition: star.h:535
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
const Tensor & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
Definition: tensor.C:1011
Class for stars in binary system.
Definition: star.h:483
Metric gtilde
Conformal metric .
Definition: star.h:565
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition: tensor.C:825
Vector dcon_logn
Contravariant derivative of the total logarithm of the lapse.
Definition: star.h:538
void import(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
Definition: scalar_import.C:71
Coord ya
Absolute y coordinate.
Definition: map.h:749
const Tensor & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
Definition: tensor.C:1023
Multi-domain grid.
Definition: grilles.h:279
Bases of the spectral expansions.
Definition: base_val.h:325
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:809
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition: tensor.C:663
Scalar nn
Lapse function N .
Definition: star.h:225
Vector beta_comp
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition: star.h:575
Sym_tensor tkij_auto
Part of the extrinsic curvature tensor generated by beta_auto.
Definition: star.h:600
Scalar kcar_comp
Part of the scalar generated by beta_auto and beta_comp, i.e.
Definition: star.h:618
Valeur & set_spectral_va()
Returns va (read/write version)
Definition: scalar.h:610
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:302
void update_metric_der_comp(const Star_bin &comp, double omega)
Computes the derivative of metric functions related to the companion star.
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
Definition: scalar_deriv.C:390
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:874
Sym_tensor derive_lie(const Vector &v) const
Computes the Lie derivative of this with respect to some vector field v.
Definition: sym_tensor.C:363
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:226
Vector dcon_phi
Contravariant derivative of the logarithm of the conformal factor.
Definition: star.h:557