LORENE
vector.C
1 /*
2  * Methods of class Vector
3  *
4  * (see file vector.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
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 vector_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/vector.C,v 1.30 2014/10/13 08:53:44 j_novak Exp $" ;
31 
32 /*
33  * $Id: vector.C,v 1.30 2014/10/13 08:53:44 j_novak Exp $
34  * $Log: vector.C,v $
35  * Revision 1.30 2014/10/13 08:53:44 j_novak
36  * Lorene classes and functions now belong to the namespace Lorene.
37  *
38  * Revision 1.29 2014/10/06 15:13:20 j_novak
39  * Modified #include directives to use c++ syntax.
40  *
41  * Revision 1.28 2008/10/29 14:09:14 jl_cornou
42  * Spectral bases for pseudo vectors and curl added
43  *
44  * Revision 1.27 2008/08/27 08:52:23 jl_cornou
45  * Added fonctions for angular potential A
46  *
47  * Revision 1.26 2007/12/21 16:07:08 j_novak
48  * Methods to filter Tensor, Vector and Sym_tensor objects.
49  *
50  * Revision 1.25 2005/02/14 13:01:50 j_novak
51  * p_eta and p_mu are members of the class Vector. Most of associated functions
52  * have been moved from the class Vector_divfree to the class Vector.
53  *
54  * Revision 1.24 2005/01/25 15:37:35 j_novak
55  * Solved some dzpuis problem...
56  *
57  * Revision 1.23 2005/01/12 16:48:23 j_novak
58  * Better treatment of the case where all vector components are null in
59  * decompose_div .
60  *
61  * Revision 1.22 2004/10/12 09:58:25 j_novak
62  * Better memory management.
63  *
64  * Revision 1.21 2004/10/11 09:46:31 j_novak
65  * Speed improvements.
66  *
67  * Revision 1.20 2004/05/09 20:55:05 e_gourgoulhon
68  * Added method flux.
69  *
70  * Revision 1.19 2004/03/29 11:57:45 e_gourgoulhon
71  * Added methods ope_killing and ope_killing_conf.
72  *
73  * Revision 1.18 2004/02/26 22:48:50 e_gourgoulhon
74  * -- Method divergence: call to Tensor::divergence and cast of the
75  * result.
76  * -- Added method derive_lie.
77  *
78  * Revision 1.17 2004/02/24 09:46:20 j_novak
79  * Correction to cope with SGI compiler's warnings.
80  *
81  * Revision 1.16 2004/02/20 10:53:04 j_novak
82  * Added the matching of the potential adapted to the case where the
83  * vector is the source of a Poisson equation (dzpuis = 4).
84  *
85  * Revision 1.15 2004/01/30 10:30:17 j_novak
86  * Changed dzpuis handling in Vector::decompose_div (this may be temporary).
87  *
88  * Revision 1.14 2003/12/30 23:09:47 e_gourgoulhon
89  * Change in methods derive_cov() and divergence() to take into account
90  * the change of name: Metric::get_connect() --> Metric::connect().
91  *
92  * Revision 1.13 2003/12/19 15:18:16 j_novak
93  * Shadow variables hunt
94  *
95  * Revision 1.12 2003/10/29 11:04:34 e_gourgoulhon
96  * dec2_dpzuis() replaced by dec_dzpuis(2).
97  * inc2_dpzuis() replaced by inc_dzpuis(2).
98  *
99  * Revision 1.11 2003/10/22 14:24:19 j_novak
100  * *** empty log message ***
101  *
102  * Revision 1.9 2003/10/20 13:00:38 j_novak
103  * Memory error corrected
104  *
105  * Revision 1.8 2003/10/20 09:32:12 j_novak
106  * Members p_potential and p_div_free of the Helmholtz decomposition
107  * + the method decompose_div(Metric).
108  *
109  * Revision 1.7 2003/10/16 14:21:37 j_novak
110  * The calculation of the divergence of a Tensor is now possible.
111  *
112  * Revision 1.6 2003/10/13 13:52:40 j_novak
113  * Better managment of derived quantities.
114  *
115  * Revision 1.5 2003/10/06 13:58:48 j_novak
116  * The memory management has been improved.
117  * Implementation of the covariant derivative with respect to the exact Tensor
118  * type.
119  *
120  * Revision 1.4 2003/10/05 21:14:20 e_gourgoulhon
121  * Added method std_spectral_base().
122  *
123  * Revision 1.3 2003/10/03 14:10:32 e_gourgoulhon
124  * Added constructor from Tensor.
125  *
126  * Revision 1.2 2003/10/03 14:08:46 j_novak
127  * Removed old change_trid...
128  *
129  * Revision 1.1 2003/09/26 08:05:31 j_novak
130  * New class Vector.
131  *
132  *
133  * $Header: /cvsroot/Lorene/C++/Source/Tensor/vector.C,v 1.30 2014/10/13 08:53:44 j_novak Exp $
134  *
135  */
136 
137 // Headers C
138 #include <cstdlib>
139 #include <cassert>
140 #include <cmath>
141 
142 // Headers Lorene
143 #include "metric.h"
144 #include "proto.h"
145 #include "matrice.h"
146 #include "nbr_spx.h"
147 
148 
149  //--------------//
150  // Constructors //
151  //--------------//
152 
153 // Standard constructor
154 // --------------------
155 namespace Lorene {
156 Vector::Vector(const Map& map, int tipe, const Base_vect& triad_i)
157  : Tensor(map, 1, tipe, triad_i) {
158 
159  set_der_0x0() ;
160 
161 }
162 
163 // Standard constructor with the triad passed as a pointer
164 // -------------------------------------------------------
165 Vector::Vector(const Map& map, int tipe, const Base_vect* triad_i)
166  : Tensor(map, 1, tipe, *triad_i) {
167 
168  set_der_0x0() ;
169 }
170 
171 // Copy constructor
172 // ----------------
173 Vector::Vector (const Vector& source) :
174  Tensor(source) {
175 
176  assert(valence == 1) ;
177  set_der_0x0() ;
178 
179 }
180 
181 
182 // Constructor from a {\tt Tensor}.
183 //--------------------------------
184 Vector::Vector(const Tensor& uu) : Tensor(uu) {
185 
186  assert(valence == 1) ;
187  set_der_0x0() ;
188 
189 }
190 
191 
192 // Constructor from a file
193 // -----------------------
194 Vector::Vector(const Map& mapping, const Base_vect& triad_i, FILE* fd) :
195  Tensor(mapping, triad_i, fd) {
196 
197  assert ( (valence == 1) && (n_comp == 3) ) ;
198  set_der_0x0() ;
199 
200 }
201 
202 
203  //--------------//
204  // Destructor //
205  //--------------//
206 
207 
209 
211 
212 }
213 
214 
215  //-------------------//
216  // Memory managment //
217  //-------------------//
218 
219 void Vector::del_deriv() const {
220 
221  for (int i=0; i<N_MET_MAX; i++)
222  del_derive_met(i) ;
223 
224  if (p_A != 0x0) delete p_A ;
225  if (p_eta != 0x0) delete p_eta ;
226  if (p_mu != 0x0) delete p_mu ;
227  set_der_0x0() ;
229 
230 }
231 
232 void Vector::set_der_0x0() const {
233 
234  for (int i=0; i<N_MET_MAX; i++)
235  set_der_met_0x0(i) ;
236  p_A = 0x0 ;
237  p_eta = 0x0 ;
238  p_mu = 0x0 ;
239 
240 }
241 
242 void Vector::del_derive_met(int j) const {
243 
244  assert( (j>=0) && (j<N_MET_MAX) ) ;
245 
246  if (met_depend[j] != 0x0) {
247  if (p_potential[j] != 0x0)
248  delete p_potential[j] ;
249  if (p_div_free[j] != 0x0)
250  delete p_div_free[j] ;
251 
252  set_der_met_0x0(j) ;
253 
255  }
256 }
257 
258 void Vector::set_der_met_0x0(int i) const {
259 
260  assert( (i>=0) && (i<N_MET_MAX) ) ;
261 
262  p_potential[i] = 0x0 ;
263  p_div_free[i] = 0x0 ;
264 
265 }
266 
267 void Vector::operator=(const Vector& t) {
268 
269  triad = t.triad ;
270 
271  assert(t.type_indice(0) == type_indice(0)) ;
272 
273  for (int i=0 ; i<3 ; i++) {
274  *cmp[i] = *t.cmp[i] ;
275  }
276  del_deriv() ;
277 }
278 
279 void Vector::operator=(const Tensor& t) {
280 
281  assert (t.valence == 1) ;
282 
283  triad = t.triad ;
284 
285  assert(t.type_indice(0) == type_indice(0)) ;
286 
287  for (int i=0 ; i<3 ; i++) {
288  *cmp[i] = *t.cmp[i] ;
289  }
290  del_deriv() ;
291 }
292 
293 
294 
295 // Affectation d'une composante :
296 Scalar& Vector::set(int index) {
297 
298  assert ( (index>=1) && (index<=3) ) ;
299 
300  del_deriv() ;
301 
302  return *cmp[index - 1] ;
303 }
304 
305 const Scalar& Vector::operator()(int index) const {
306 
307  assert ( (index>=1) && (index<=3) ) ;
308 
309  return *cmp[index - 1] ;
310 
311 }
312 
313 
314 // Sets the standard spectal bases of decomposition for each component
315 
317 
318  Base_val** bases = 0x0 ;
319 
320  if ( triad->identify() == (mp->get_bvect_cart()).identify() ) {
321 
322  // Cartesian case
323  bases = mp->get_mg()->std_base_vect_cart() ;
324 
325  }
326  else {
327  // Spherical case
328  assert( triad->identify() == (mp->get_bvect_spher()).identify()) ;
329  bases = mp->get_mg()->std_base_vect_spher() ;
330  }
331 
332  for (int i=0 ; i<3 ; i++) {
333  cmp[i]->set_spectral_base( *bases[i] ) ;
334  }
335 
336  for (int i=0 ; i<3 ; i++) {
337  delete bases[i] ;
338  }
339  delete [] bases ;
340 
341 
342 }
343 
344 // Sets the standard spectral bases of decomposition for each component for a pseudo vector
345 
347 
348  Base_val** bases = 0x0 ;
349 
350  if ( triad->identify() == (mp->get_bvect_cart()).identify() ) {
351 
352  // Cartesian case
353  bases = mp->get_mg()->pseudo_base_vect_cart() ;
354 
355  }
356  else {
357  // Spherical case
358  assert( triad->identify() == (mp->get_bvect_spher()).identify()) ;
359  bases = mp->get_mg()->pseudo_base_vect_spher() ;
360  }
361 
362  for (int i=0 ; i<3 ; i++) {
363  cmp[i]->set_spectral_base( *bases[i] ) ;
364  }
365 
366  for (int i=0 ; i<3 ; i++) {
367  delete bases[i] ;
368  }
369  delete [] bases ;
370 
371 
372 }
373 
374 
375 
376  //-------------------------------//
377  // Computational methods //
378  //-------------------------------//
379 
380 
381 const Scalar& Vector::divergence(const Metric& metre) const {
382 
383  const Scalar* pscal =
384  dynamic_cast<const Scalar*>( &(Tensor::divergence(metre)) ) ;
385 
386  assert(pscal != 0x0) ;
387 
388  return *pscal ;
389 }
390 
391 
392 Vector Vector::derive_lie(const Vector& vv) const {
393 
394  Vector resu(*mp, type_indice(0), triad) ;
395 
396  compute_derive_lie(vv, resu) ;
397 
398  return resu ;
399 
400 }
401 
403 
404  if ( triad->identify() == (mp->get_bvect_cart()).identify() ) {
405  const Metric_flat& metc = mp->flat_met_cart() ;
406  Vector_divfree resu(*mp, mp->get_bvect_cart(), metc) ;
407  resu.set(1)= cmp[3]->dsdy() - cmp[2]->dsdz();
408  resu.set(2)= cmp[1]->dsdz() - cmp[3]->dsdx();
409  resu.set(3)= cmp[2]->dsdx() - cmp[1]->dsdy();
410  resu.pseudo_spectral_base();
411  return resu ;
412  }
413  else {
414  assert( triad->identify() == (mp->get_bvect_spher()).identify()) ;
415  const Metric_flat& mets = mp->flat_met_spher() ;
416  Vector_divfree resu(*mp, mp->get_bvect_spher(), mets);
417  Scalar tmp = *cmp[3] ;
418  tmp.div_tant();
419  tmp += cmp[3]->dsdt();
420  tmp.div_r();
421  resu.set(1) = tmp - cmp[2]->srstdsdp() ;
422  tmp = *cmp[3] ;
423  tmp.mult_r();
424  tmp = tmp.dsdr();
425  tmp.div_r();
426  resu.set(2) = cmp[1]->srstdsdp() - tmp ;
427  tmp = *cmp[2];
428  tmp.mult_r();
429  resu.set(3) = tmp.dsdr() - cmp[1]->dsdt() ;
430  resu.set(3).div_r();
431  resu.pseudo_spectral_base();
432  return resu ;
433  }
434 
435 }
436 
437 
439 
440  Sym_tensor resu(*mp, type_indice(0), *triad) ;
441 
442  if (type_indice(0) == CON ) {
443  for (int i=1; i<=3; i++) {
444  for (int j=i; j<=3; j++) {
445  resu.set(i,j) = derive_con(gam)(i,j) + derive_con(gam)(j,i) ;
446  }
447  }
448  }
449  else {
450  for (int i=1; i<=3; i++) {
451  for (int j=i; j<=3; j++) {
452  resu.set(i,j) = derive_cov(gam)(i,j) + derive_cov(gam)(j,i) ;
453  }
454  }
455  }
456 
457  return resu ;
458 
459 }
460 
461 
463 
464  Sym_tensor resu(*mp, type_indice(0), *triad) ;
465 
466  if (type_indice(0) == CON ) {
467  for (int i=1; i<=3; i++) {
468  for (int j=i; j<=3; j++) {
469  resu.set(i,j) = derive_con(gam)(i,j) + derive_con(gam)(j,i)
470  - 0.6666666666666666* divergence(gam)
471  * gam.con()(i,j) ;
472  }
473  }
474  }
475  else {
476  for (int i=1; i<=3; i++) {
477  for (int j=i; j<=3; j++) {
478  resu.set(i,j) = derive_cov(gam)(i,j) + derive_cov(gam)(j,i)
479  - 0.6666666666666666* derive_con(gam).trace()
480  * gam.cov()(i,j) ;
481  }
482  }
483  }
484 
485  return resu ;
486 
487 }
488 
489 
490 
491 
492 const Scalar& Vector::potential(const Metric& metre) const {
493 
494  set_dependance(metre) ;
495  int j = get_place_met(metre) ;
496  assert ((j>=0) && (j<N_MET_MAX)) ;
497  if (p_potential[j] == 0x0) {
498  decompose_div(metre) ;
499  }
500 
501  return *p_potential[j] ;
502 }
503 
504 const Vector_divfree& Vector::div_free(const Metric& metre) const {
505 
506  set_dependance(metre) ;
507  int j = get_place_met(metre) ;
508  assert ((j>=0) && (j<N_MET_MAX)) ;
509  if (p_div_free[j] == 0x0) {
510  decompose_div(metre) ;
511  }
512 
513  return *p_div_free[j] ;
514 }
515 
516 void Vector::decompose_div(const Metric& metre) const {
517 
518  assert( type_indice(0) == CON ) ; //Only for contravariant vectors...
519 
520  set_dependance(metre) ;
521  int ind = get_place_met(metre) ;
522  assert ((ind>=0) && (ind<N_MET_MAX)) ;
523 
524  if ( (p_potential[ind] != 0x0) && (p_div_free[ind] != 0x0) )
525  return ; // Nothing to do ...
526 
527  else {
528  if (p_div_free[ind] != 0x0)
529  delete p_div_free[ind] ;
530  if (p_potential[ind] != 0x0)
531  delete p_potential[ind] ;
532 
533  const Mg3d* mmg = mp->get_mg() ;
534  int nz = mmg->get_nzone() ;
535 
536  int dzp = cmp[0]->get_dzpuis() ;
537 #ifndef NDEBUG
538  bool dz_zero = cmp[0]->check_dzpuis(0) ;
539  bool dz_four = cmp[0]->check_dzpuis(4) ;
540 #endif
541  assert( dz_zero || dz_four) ;
542  assert (cmp[1]->check_dzpuis(dzp)) ;
543  assert (cmp[2]->check_dzpuis(dzp)) ;
544 
545  Scalar dive = divergence(metre) ;
546 
547  if (dive.get_etat() == ETATZERO) {
548  p_potential[ind] = new Scalar(*mp) ;
549  p_potential[ind]->set_etat_zero() ;
550  p_potential[ind]->set_dzpuis(dzp) ;
551  }
552  else {
553  //No matching of the solution at this point
554  p_potential[ind] = new Scalar(dive.poisson()) ;
555 
556  if (dzp == 4) {
557  assert (mmg->get_type_r(nz-1) == UNSURR) ;
558  // Let's now do the matching ... in the case dzpuis = 4
559  const Map_af* mapping = dynamic_cast<const Map_af*>(mp) ;
560  assert (mapping != 0x0) ;
561  Valeur& val = p_potential[ind]->set_spectral_va() ;
562  val.ylm() ;
563  Mtbl_cf& mtc = *val.c_cf ;
564  if (nz > 1) {
565  int np = mmg->get_np(0) ;
566  int nt = mmg->get_nt(0) ;
567 #ifndef NDEBUG
568  for (int lz=0; lz<nz; lz++) {
569  assert (mmg->get_np(lz) == np) ;
570  assert (mmg->get_nt(lz) == nt) ;
571  }
572 #endif
573  int lmax = 0 ;
574  for (int k=0; k<np+1; k++)
575  for (int j=0; j<nt; j++)
576  if ( nullite_plm(j, nt, k, np, val.base)) {
577  int m_quant, l_quant, base_r ;
578  donne_lm (nz, 0, j, k, val.base, m_quant,
579  l_quant, base_r) ;
580  lmax = (l_quant > lmax ? l_quant : lmax) ;
581  }
582  Scalar** ri = new Scalar*[lmax+1] ;
583  Scalar** rmi = new Scalar*[lmax+1] ;
584  Scalar erre(*mp) ;
585  erre = mp->r ;
586  for (int l=0; l<=lmax; l++) {
587  ri[l] = new Scalar(*mp) ;
588  rmi[l] = new Scalar(*mp) ;
589  if (l == 0) *(ri[l]) = 1. ;
590  else *(ri[l]) = pow(erre, l) ;
591  ri[l]->annule_domain(nz - 1) ;
592  ri[l]->std_spectral_base() ; //exact base in r will be set later
593  if (l==2) *(rmi[l]) = 1 ;
594  else *(rmi[l]) = pow(erre, 2-l) ;
595  rmi[l]->annule(0,nz-2) ;
596  Scalar tmp = pow(erre, -l-1) ;
597  tmp.annule_domain(nz-1) ;
598  tmp.annule_domain(0) ;
599  *(rmi[l]) += tmp ;
600  rmi[l]->set_dzpuis(3) ;
601  rmi[l]->std_spectral_base() ;//exact base in r will be set later
602  }
603 
604  for (int k=0; k<np+1; k++) {
605  for (int j=0; j<nt; j++) {
606  if ( nullite_plm(j, nt, k, np, val.base)) {
607  int m_quant, l_quant, base_r, l, c ;
608  donne_lm (nz, 0, j, k, val.base, m_quant, l_quant, base_r) ;
609  bool lzero = (l_quant == 0) || (l_quant == 1) ;
610  int nb_hom_ced = (l_quant < 2 ? 0 : 1) ;
611 
612  int taille = 2*(nz-1) - 1 + nb_hom_ced ;
613  Tbl deuz(taille) ;
614  deuz.set_etat_qcq() ;
615  Matrice systeme(taille,taille) ;
616  systeme.set_etat_qcq() ;
617  for (l=0; l<taille; l++)
618  for (c=0; c<taille; c++) systeme.set(l,c) = 0 ;
619  for (l=0; l<taille; l++) deuz.set(l) = 0 ;
620 
621  //---------
622  // Nucleus
623  //---------
624  assert(mmg->get_type_r(0) == RARE) ;
625  assert ((base_r == R_CHEBP)||(base_r == R_CHEBI)) ;
626  ri[l_quant]->set_spectral_va().set_base_r(0, base_r) ;
627 
628  double alpha = mapping->get_alpha()[0] ;
629  int nr = mmg->get_nr(0) ;
630  Tbl partn(nr) ;
631  partn.set_etat_qcq() ;
632  for (int i=0; i<nr; i++)
633  partn.set(i) = mtc(0, k, j, i) ;
634  l=0 ; c=0 ;
635 
636  systeme.set(l,c) += pow(alpha, double(l_quant)) ;
637 
638  deuz.set(l) -= val1_dern_1d(0, partn, base_r) ;
639  l++ ;
640 
641  if (!lzero) {
642  systeme.set(l,c) += l_quant * pow(alpha, double(l_quant - 1)) ;
643 
644  deuz.set(l) -= val1_dern_1d(1, partn, base_r) / alpha ;
645  }
646 
647  //----------
648  // Shells
649  //----------
650 
651  for (int lz=1; lz<nz-1; lz++) {
652  alpha = mapping->get_alpha()[lz] ;
653  double beta = mapping->get_beta()[lz] ;
654  double xm1 = beta - alpha ;
655  double xp1 = alpha + beta ;
656  nr = mmg->get_nr(lz) ;
657  Tbl parts(nr) ;
658  parts.set_etat_qcq() ;
659  for (int i=0; i<nr; i++)
660  parts.set(i) = mtc(lz, k, j, i) ;
661 
662  //Function at x = -1
663  l-- ;
664  c = 2*lz - 1 ;
665  systeme.set(l,c) -= pow(xm1, double(l_quant)) ;
666  c++;
667  systeme.set(l,c) -= pow(xm1, double(-l_quant-1)) ;
668 
669  deuz.set(l) += valm1_dern_1d(0, parts, R_CHEB) ;
670 
671  if ((lz>1) || (!lzero)) {
672  //First derivative at x=-1
673  l++ ;
674  c-- ;
675  systeme.set(l,c) -= l_quant*pow(xm1, double(l_quant - 1)) ;
676  c++;
677  systeme.set(l,c) -= (-l_quant - 1)*
678  pow(xm1, double(-l_quant - 2)) ;
679 
680  deuz.set(l) += valm1_dern_1d(1, parts, R_CHEB) / alpha ;
681  }
682 
683  //Function at x = 1
684  l++ ;
685  c-- ;
686  systeme.set(l,c) += pow(xp1, double(l_quant)) ;
687  c++;
688  systeme.set(l,c) += pow(xp1, double(-l_quant-1)) ;
689 
690  deuz.set(l) -= val1_dern_1d(0, parts, R_CHEB) ;
691 
692  //First derivative at x = 1
693  l++ ;
694  c-- ;
695  systeme.set(l,c) += l_quant*pow(xp1, double(l_quant - 1)) ;
696  c++;
697  systeme.set(l,c) += (-l_quant - 1)*
698  pow(xp1, double(-l_quant - 2)) ;
699 
700  deuz.set(l) -= val1_dern_1d(1, parts, R_CHEB) / alpha ;
701 
702  } //End of the loop on different shells
703 
704  //-------------------------------
705  // Compactified external domain
706  //-------------------------------
707  assert(mmg->get_type_r(nz-1) == UNSURR) ;
708  nr = mmg->get_nr(nz-1) ;
709  Tbl partc(nr) ;
710  partc.set_etat_qcq() ;
711  for (int i=0; i<nr; i++)
712  partc.set(i) = mtc(nz-1, k, j, i) ;
713 
714  alpha = mapping->get_alpha()[nz-1] ;
715  double beta = mapping->get_beta()[nz-1] ;
716  double xm1 = beta - alpha ; // 1 / r_left
717  double part0, part1 ;
718  part0 = valm1_dern_1d(0, partc, R_CHEB) ;
719  part1 = xm1*valm1_dern_1d(1, partc, R_CHEB) / alpha ;
720  assert (p_potential[ind]->get_dzpuis() == 3) ;
721 
722  //Function at x = -1
723  l--;
724  if (!lzero) {
725  c++;
726  systeme.set(l,c) -= pow(xm1, double(l_quant+1)) ;
727  }
728  deuz.set(l) += part0*xm1*xm1*xm1 ;
729 
730  // First derivative at x = -1
731  l++ ;
732  if (!lzero) {
733  systeme.set(l,c) -= (-l_quant - 1)*
734  pow(xm1, double(l_quant + 2)) ;
735  }
736  if ( (nz > 2) || (!lzero))
737  deuz.set(l) += -xm1*xm1*xm1*xm1*(3*part0 + part1) ;
738 
739  //--------------------------------------
740  // Solution of the linear system
741  //--------------------------------------
742 
743  int inf = 1 + nb_hom_ced;
744  int sup = 3 - nb_hom_ced;
745  systeme.set_band(sup, inf) ;
746  systeme.set_lu() ;
747  Tbl facteur(systeme.inverse(deuz)) ;
748  ri[l_quant]->set_spectral_va().coef() ;
749  rmi[l_quant]->set_spectral_va().coef() ;
750 
751  //Linear combination in the nucleus
752  nr = mmg->get_nr(0) ;
753  for (int i=0; i<nr; i++)
754  mtc.set(0, k, j, i) +=
755  facteur(0)*(*(ri[l_quant]->get_spectral_va().c_cf))(0, 0, 0, i) ;
756 
757  //Linear combination in the shells
758  for (int lz=1; lz<nz-1; lz++) {
759  nr = mmg->get_nr(lz) ;
760  for (int i=0; i<nr; i++)
761  mtc.set(lz, k, j, i) += facteur(2*lz - 1)*
762  (*(ri[l_quant]->get_spectral_va().c_cf))(lz, 0, 0, i) ;
763  for (int i=0; i<nr; i++)
764  mtc.set(lz, k, j, i) += facteur(2*lz)*
765  (*(rmi[l_quant]->get_spectral_va().c_cf))(lz, 0, 0, i) ;
766  }
767 
768  //Linear combination in the CED
769  nr = mmg->get_nr(nz-1) ;
770  if (!lzero) {
771  for (int i=0; i<nr; i++)
772  mtc.set(nz-1, k, j, i) +=
773  facteur(taille - 1)*
774  (*(rmi[l_quant]->get_spectral_va().c_cf))(nz-1, 0, 0, i) ;
775  }
776  } //End of nullite_plm ...
777 
778  } //End of j/theta loop
779  } //End of k/phi loop
780 
781  for (int l=0; l<=lmax; l++) {
782  delete ri[l] ;
783  delete rmi[l] ;
784  }
785  delete [] ri ;
786  delete [] rmi ;
787 
788  } //End of the case of more than one domain
789 
790  val.ylm_i() ;
791 
792  } //End of the case dzp = 4
793  }
794  p_div_free[ind] = new Vector_divfree(*mp, *triad, metre) ;
795 
796  Vector gradient = p_potential[ind]->derive_con(metre) ;
797  if (dzp != 4) gradient.dec_dzpuis(2) ;
798 
799  *p_div_free[ind] = ( *this - gradient ) ;
800 
801  }
802 
803 }
804 
805 
806 
807 double Vector::flux(double radius, const Metric& met) const {
808 
809  assert(type_indice(0) == CON) ;
810 
811  const Map_af* mp_af = dynamic_cast<const Map_af*>(mp) ;
812  if (mp_af == 0x0) {
813  cerr <<
814  "Vector::flux : the case of a mapping outside the class Map_af\n"
815  << " is not implemented yet !" << endl ;
816  abort() ;
817  }
818 
819  const Metric_flat* ff = dynamic_cast<const Metric_flat*>(&met) ;
820  if (ff == 0x0) {
821  cerr <<
822  "Vector::flux : the case of a non flat metric is not implemented yet !"
823  << endl ;
824  abort() ;
825  }
826 
827  const Base_vect_cart* bcart = dynamic_cast<const Base_vect_cart*>(triad) ;
828  Vector* vspher = 0x0 ;
829  if (bcart != 0x0) { // switch to spherical components:
830  vspher = new Vector(*this) ;
831  vspher->change_triad(mp->get_bvect_spher()) ;
832  }
833  else assert( dynamic_cast<const Base_vect_spher*>(triad) != 0x0 ) ;
834 
835  const Vector* vv = (bcart != 0x0) ? vspher : this ;
836 
837  const Scalar& vr = vv->operator()(1) ;
838 
839  double resu ;
840  if (radius == __infinity) {
841  resu = mp_af->integrale_surface_infini(vr) ;
842  }
843  else {
844  resu = mp_af->integrale_surface(vr, radius) ;
845  }
846 
847  return resu ;
848 }
849 
850 void Vector::exponential_filter_r(int lzmin, int lzmax, int p,
851  double alpha) {
852  if( triad->identify() == (mp->get_bvect_cart()).identify() )
853  for (int i=0; i<n_comp; i++)
854  cmp[i]->exponential_filter_r(lzmin, lzmax, p, alpha) ;
855  else {
856  assert( triad->identify() == (mp->get_bvect_spher()).identify()) ;
857  Scalar vr_tmp = operator()(1) ;
858  vr_tmp.exponential_filter_r(lzmin, lzmax, p, alpha) ;
859  Scalar eta_tmp = eta() ;
860  eta_tmp.exponential_filter_r(lzmin, lzmax, p, alpha) ;
861  Scalar mu_tmp = mu() ;
862  mu_tmp.exponential_filter_r(lzmin, lzmax, p, alpha) ;
863  set_vr_eta_mu(vr_tmp, eta_tmp, mu_tmp) ;
864  }
865 }
866 
867 void Vector::exponential_filter_ylm(int lzmin, int lzmax, int p,
868  double alpha) {
869  if( triad->identify() == (mp->get_bvect_cart()).identify() )
870  for (int i=0; i<n_comp; i++)
871  cmp[i]->exponential_filter_ylm(lzmin, lzmax, p, alpha) ;
872  else {
873  assert( triad->identify() == (mp->get_bvect_spher()).identify()) ;
874  Scalar vr_tmp = operator()(1) ;
875  vr_tmp.exponential_filter_ylm(lzmin, lzmax, p, alpha) ;
876  Scalar eta_tmp = eta() ;
877  eta_tmp.exponential_filter_ylm(lzmin, lzmax, p, alpha) ;
878  Scalar mu_tmp = mu() ;
879  mu_tmp.exponential_filter_ylm(lzmin, lzmax, p, alpha) ;
880  set_vr_eta_mu(vr_tmp, eta_tmp, mu_tmp) ;
881  }
882 }
883 }
Lorene::Vector::exponential_filter_r
virtual void exponential_filter_r(int lzmin, int lzmax, int p, double alpha=-16.)
Applies exponential filters to all components (see Scalar::exponential_filter_r ).
Definition: vector.C:850
Lorene::Tensor::set_dependance
void set_dependance(const Metric &) const
To be used to describe the fact that the derivatives members have been calculated with met .
Definition: tensor.C:453
R_CHEB
#define R_CHEB
base de Chebychev ordinaire (fin)
Definition: type_parite.h:166
Lorene::Vector::ope_killing
Sym_tensor ope_killing(const Metric &gam) const
Computes the Killing operator associated with a given metric.
Definition: vector.C:438
Lorene::Tbl::set
double & set(int i)
Read/write of a particular element (index i) (1D case)
Definition: tbl.h:281
Lorene::Base_val
Bases of the spectral expansions.
Definition: base_val.h:322
Lorene::Mg3d::get_np
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
Lorene::Matrice::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: matrice.C:175
Lorene::Valeur
Values and coefficients of a (real-value) function.
Definition: valeur.h:287
Lorene::Sym_tensor
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:223
Lorene::Tensor::derive_cov
const Tensor & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
Definition: tensor.C:1002
Lorene::Metric
Metric for tensor calculation.
Definition: metric.h:90
Lorene::Scalar::set_spectral_base
void set_spectral_base(const Base_val &)
Sets the spectral bases of the Valeur va
Definition: scalar.C:797
Lorene::Scalar::dsdy
const Scalar & dsdy() const
Returns of *this , where .
Definition: scalar_deriv.C:297
Lorene::Mg3d
Multi-domain grid.
Definition: grilles.h:273
Lorene::Tensor::del_deriv
virtual void del_deriv() const
Deletes the derived quantities.
Definition: tensor.C:398
Lorene::Tensor
Tensor handling.
Definition: tensor.h:288
Lorene::Vector::operator()
const Scalar & operator()(int) const
Readonly access to a component.
Definition: vector.C:305
Lorene::Scalar::dsdr
const Scalar & dsdr() const
Returns of *this .
Definition: scalar_deriv.C:113
Lorene::Valeur::c_cf
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:302
Lorene::Base_vect_cart
Cartesian vectorial bases (triads).
Definition: base_vect.h:201
Lorene::Mg3d::get_nt
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Vector::p_potential
Scalar * p_potential[N_MET_MAX]
The potential giving the gradient part in the Helmholtz decomposition of any 3D vector .
Definition: vector.h:198
Lorene::Map::get_bvect_spher
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:783
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Map::flat_met_spher
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
Definition: map.C:321
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::Metric::cov
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric.C:280
Lorene::Map::r
Coord r
r coordinate centered on the grid
Definition: map.h:718
Lorene::Tensor::cmp
Scalar ** cmp
Array of size n_comp of pointers onto the components.
Definition: tensor.h:315
Lorene::Tensor::get_place_met
int get_place_met(const Metric &) const
Returns the position of the pointer on metre in the array met_depend .
Definition: tensor.C:443
Lorene::Vector::divergence
const Scalar & divergence(const Metric &) const
The divergence of this with respect to a Metric .
Definition: vector.C:381
Lorene::Mg3d::get_type_r
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition: grilles.h:474
Lorene::Scalar::srstdsdp
const Scalar & srstdsdp() const
Returns of *this .
Definition: scalar_deriv.C:177
Lorene::Tensor::met_depend
const Metric * met_depend[N_MET_MAX]
Array on the Metric 's which were used to compute derived quantities, like p_derive_cov ,...
Definition: tensor.h:327
Lorene::Scalar::check_dzpuis
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
Definition: scalar.C:873
Lorene::Map_af::integrale_surface
double integrale_surface(const Cmp &ci, double rayon) const
Performs the surface integration of ci on the sphere of radius rayon .
Definition: map_af_integ_surf.C:90
Lorene::Metric_flat
Flat metric for tensor calculation.
Definition: metric.h:261
Lorene::Scalar::exponential_filter_ylm
virtual void exponential_filter_ylm(int lzmin, int lzmax, int p, double alpha=-16.)
Applies an exponential filter to the spectral coefficients in the angular directions.
Definition: scalar_exp_filter.C:142
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
Lorene::Valeur::set_base_r
void set_base_r(int l, int base_r)
Sets the expansion basis for r ( ) functions in a given domain.
Definition: valeur.C:836
Lorene::Vector_divfree
Divergence-free vectors.
Definition: vector.h:724
Lorene::Tensor::annule_domain
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition: tensor.C:666
Lorene::Mtbl_cf::set
Tbl & set(int l)
Read/write of the Tbl containing the coefficients in a given domain.
Definition: mtbl_cf.h:294
R_CHEBP
#define R_CHEBP
base de Cheb. paire (rare) seulement
Definition: type_parite.h:168
Lorene::Vector::change_triad
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: vector_change_triad.C:75
Lorene::Scalar::div_tant
void div_tant()
Division by .
Definition: scalar_th_manip.C:111
R_CHEBI
#define R_CHEBI
base de Cheb. impaire (rare) seulement
Definition: type_parite.h:170
Lorene::Valeur::coef
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:148
Lorene::Vector::ope_killing_conf
Sym_tensor ope_killing_conf(const Metric &gam) const
Computes the conformal Killing operator associated with a given metric.
Definition: vector.C:462
Lorene::Scalar::get_spectral_va
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:601
Lorene::Matrice::inverse
Tbl inverse(const Tbl &sec_membre) const
Solves the linear system represented by the matrix.
Definition: matrice.C:424
Lorene::Base_vect::identify
virtual int identify() const =0
Returns a number to identify the sub-classe of Base_vect the object belongs to.
Lorene::Vector::potential
const Scalar & potential(const Metric &) const
Returns the potential in the Helmholtz decomposition.
Definition: vector.C:492
Lorene::Scalar::set_spectral_va
Valeur & set_spectral_va()
Returns va (read/write version)
Definition: scalar.h:604
Lorene::Mg3d::pseudo_base_vect_cart
Base_val ** pseudo_base_vect_cart() const
Returns the standard spectral bases for the Cartesian components of a pseudo-vector.
Definition: mg3d_std_base.C:794
Lorene::Tbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:361
Lorene::Scalar::mult_r
void mult_r()
Multiplication by r everywhere; dzpuis is not changed.
Definition: scalar_r_manip.C:208
Lorene::Tensor::set
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition: tensor.C:654
Lorene::Vector::set_vr_eta_mu
void set_vr_eta_mu(const Scalar &vr_i, const Scalar &eta_i, const Scalar &mu_i)
Defines the components through potentials and (see members p_eta and p_mu ), as well as the compon...
Definition: vector_etamu.C:189
Lorene::Tensor::mp
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:295
Lorene::Vector::curl
const Vector_divfree curl() const
The curl of this with respect to a (flat) Metric .
Definition: vector.C:402
Lorene::Vector::set
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:296
Lorene::Metric::con
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Definition: metric.C:290
Lorene::Matrice::set_band
void set_band(int up, int low) const
Calculate the band storage of *std.
Definition: matrice.C:364
Lorene::Base_vect_spher
Spherical orthonormal vectorial bases (triads).
Definition: base_vect.h:308
Lorene::Vector::mu
virtual const Scalar & mu() const
Gives the field such that the angular components of the vector are written:
Definition: vector_etamu.C:98
Lorene::Vector::p_div_free
Vector_divfree * p_div_free[N_MET_MAX]
The divergence-free vector of the Helmholtz decomposition of any 3D vector .
Definition: vector.h:205
Lorene::Tbl
Basic array class.
Definition: tbl.h:161
Lorene::Map::flat_met_cart
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
Definition: map.C:331
Lorene::Scalar::dsdz
const Scalar & dsdz() const
Returns of *this , where .
Definition: scalar_deriv.C:328
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Scalar::dsdx
const Scalar & dsdx() const
Returns of *this , where .
Definition: scalar_deriv.C:266
Lorene::Tensor::triad
const Base_vect * triad
Vectorial basis (triad) with respect to which the tensor components are defined.
Definition: tensor.h:303
Lorene::Mg3d::get_nr
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
Lorene::Vector::decompose_div
void decompose_div(const Metric &) const
Makes the Helmholtz decomposition (see documentation of p_potential ) of this with respect to a given...
Definition: vector.C:516
Lorene::Map_af
Affine radial mapping.
Definition: map.h:2027
Lorene::Map::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::Scalar::annule
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition: scalar.C:391
Lorene::Mg3d::pseudo_base_vect_spher
Base_val ** pseudo_base_vect_spher() const
Returns the standard spectral bases for the spherical components of a pseudo-vector.
Definition: mg3d_std_base.C:1114
Lorene::Mg3d::std_base_vect_cart
Base_val ** std_base_vect_cart() const
Returns the standard spectral bases for the Cartesian components of a vector.
Definition: mg3d_std_base.C:177
Lorene::Scalar::set_etat_zero
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:324
Lorene::Vector::exponential_filter_ylm
virtual void exponential_filter_ylm(int lzmin, int lzmax, int p, double alpha=-16.)
Applies exponential filters to all components (see Scalar::exponential_filter_ylm ).
Definition: vector.C:867
Lorene::Scalar::set_dzpuis
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:808
Lorene::Vector
Tensor field of valence 1.
Definition: vector.h:188
Lorene::Vector::~Vector
virtual ~Vector()
Destructor.
Definition: vector.C:208
Lorene::Tensor::derive_con
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
Lorene::Tensor::divergence
const Tensor & divergence(const Metric &gam) const
Computes the divergence of this with respect to some metric .
Definition: tensor.C:1055
Lorene::Vector::set_der_met_0x0
void set_der_met_0x0(int) const
Sets all the i-th components of met_depend in the class Vector (p_potential , etc....
Definition: vector.C:258
Lorene::Vector::derive_lie
Vector derive_lie(const Vector &v) const
Computes the Lie derivative of this with respect to some vector field v.
Definition: vector.C:392
Lorene::Scalar::dsdt
const Scalar & dsdt() const
Returns of *this .
Definition: scalar_deriv.C:208
Lorene::Tensor::type_indice
Itbl type_indice
1D array of integers (class Itbl ) of size valence containing the type of each index: COV for a cov...
Definition: tensor.h:310
Lorene::Vector::operator=
virtual void operator=(const Vector &a)
Assignment from a Vector.
Definition: vector.C:267
Lorene::Scalar::get_etat
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:554
Lorene::Tensor::trace
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.
Definition: tensor_calculus.C:94
Lorene::Vector::pseudo_spectral_base
virtual void pseudo_spectral_base()
Sets the standard spectal bases of decomposition for each component for a pseudo_vector.
Definition: vector.C:346
Lorene::Vector::eta
virtual const Scalar & eta() const
Gives the field such that the angular components of the vector are written:
Definition: vector_etamu.C:66
Lorene::Scalar::derive_con
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
Lorene::Tensor::valence
int valence
Valence of the tensor (0 = scalar, 1 = vector, etc...)
Definition: tensor.h:298
Lorene::Matrice
Matrix handling.
Definition: matrice.h:152
Lorene::Scalar::exponential_filter_r
virtual void exponential_filter_r(int lzmin, int lzmax, int p, double alpha=-16.)
Applies an exponential filter to the spectral coefficients in the radial direction.
Definition: scalar_exp_filter.C:60
Lorene::Map_af::get_beta
const double * get_beta() const
Returns the pointer on the array beta.
Definition: map_af.C:481
Lorene::Scalar::std_spectral_base
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
Lorene::Vector::del_derive_met
virtual void del_derive_met(int) const
Logical destructor of the derivatives depending on the i-th element of met_depend in the class Vector...
Definition: vector.C:242
Lorene::Vector::flux
double flux(double radius, const Metric &met) const
Computes the flux of the vector accross a sphere r = const.
Definition: vector.C:807
Lorene::Tensor::del_derive_met
virtual void del_derive_met(int) const
Logical destructor of the derivatives depending on the i-th element of met_depend .
Definition: tensor.C:414
Lorene::Vector::p_A
Scalar * p_A
Field defined by.
Definition: vector.h:241
Lorene::Scalar::poisson
Scalar poisson() const
Solves the scalar Poisson equation with *this as a source.
Definition: scalar_pde.C:136
Lorene::Vector::p_mu
Scalar * p_mu
Field such that the angular components of the vector are written:
Definition: vector.h:233
Lorene::Vector::del_deriv
virtual void del_deriv() const
Deletes the derived quantities.
Definition: vector.C:219
Lorene::Valeur::base
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:305
Lorene::Mg3d::std_base_vect_spher
Base_val ** std_base_vect_spher() const
Returns the standard spectral bases for the spherical components of a vector.
Definition: mg3d_std_base.C:499
Lorene::Scalar::get_dzpuis
int get_dzpuis() const
Returns dzpuis.
Definition: scalar.h:557
Lorene::Map_af::integrale_surface_infini
double integrale_surface_infini(const Cmp &ci) const
Performs the surface integration of ci at infinity.
Definition: map_af_integ_surf.C:196
Lorene::Valeur::ylm
void ylm()
Computes the coefficients of *this.
Definition: valeur_ylm.C:138
Lorene::Matrice::set
double & set(int j, int i)
Read/write of a particuliar element.
Definition: matrice.h:277
Lorene::Vector::set_der_0x0
void set_der_0x0() const
Sets the pointers on derived quantities to 0x0.
Definition: vector.C:232
Lorene::Tensor::n_comp
int n_comp
Number of stored components, depending on the symmetry.
Definition: tensor.h:312
Lorene::Map_af::get_alpha
const double * get_alpha() const
Returns the pointer on the array alpha.
Definition: map_af.C:477
Lorene::Map
Base class for coordinate mappings.
Definition: map.h:670
Lorene::Scalar::div_r
void div_r()
Division by r everywhere; dzpuis is not changed.
Definition: scalar_r_manip.C:133
Lorene::Mtbl_cf
Coefficients storage for the multi-domain spectral method.
Definition: mtbl_cf.h:186
Lorene::Base_vect
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
Lorene::Tensor::compute_derive_lie
void compute_derive_lie(const Vector &v, Tensor &resu) const
Computes the Lie derivative of this with respect to some vector field v (protected method; the public...
Definition: tensor_calculus.C:339
Lorene::Vector::p_eta
Scalar * p_eta
Field such that the angular components of the vector are written:
Definition: vector.h:219
Lorene::Valeur::ylm_i
void ylm_i()
Inverse of ylm()
Definition: valeur_ylm_i.C:131
Lorene::Vector::div_free
const Vector_divfree & div_free(const Metric &) const
Returns the div-free vector in the Helmholtz decomposition.
Definition: vector.C:504
Lorene::Matrice::set_lu
void set_lu() const
Calculate the LU-representation, assuming the band-storage has been done.
Definition: matrice.C:392
Lorene::Vector::std_spectral_base
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316