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