LORENE
scalar_import.C
1 /*
2  * Member function of the Scalar class for initiating a Scalar from
3  * a Scalar defined on another mapping.
4  */
5 
6 /*
7  * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
8  * Copyright (c) 1999-2001 Eric Gourgoulhon (Cmp version)
9  *
10  * This file is part of LORENE.
11  *
12  * LORENE is free software; you can redistribute it and/or modify
13  * it under the terms of the GNU General Public License as published by
14  * the Free Software Foundation; either version 2 of the License, or
15  * (at your option) any later version.
16  *
17  * LORENE is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU General Public License for more details.
21  *
22  * You should have received a copy of the GNU General Public License
23  * along with LORENE; if not, write to the Free Software
24  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25  *
26  */
27 
28 
29 
30 
31 /*
32  * $Id: scalar_import.C,v 1.6 2016/12/05 16:18:18 j_novak Exp $
33  * $Log: scalar_import.C,v $
34  * Revision 1.6 2016/12/05 16:18:18 j_novak
35  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
36  *
37  * Revision 1.5 2014/10/13 08:53:46 j_novak
38  * Lorene classes and functions now belong to the namespace Lorene.
39  *
40  * Revision 1.4 2014/10/06 15:16:15 j_novak
41  * Modified #include directives to use c++ syntax.
42  *
43  * Revision 1.3 2003/10/10 15:57:29 j_novak
44  * Added the state one (ETATUN) to the class Scalar
45  *
46  * Revision 1.2 2003/10/01 13:04:44 e_gourgoulhon
47  * The method Tensor::get_mp() returns now a reference (and not
48  * a pointer) onto a mapping.
49  *
50  * Revision 1.1 2003/09/25 09:07:05 j_novak
51  * Added the functions for importing from another mapping (to be tested).
52  *
53  *
54  * $Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_import.C,v 1.6 2016/12/05 16:18:18 j_novak Exp $
55  *
56  */
57 
58 // Headers C
59 #include <cmath>
60 
61 // Headers Lorene
62 #include "tensor.h"
63 #include "param.h"
64 #include "nbr_spx.h"
65 
66  //-------------------------------//
67  // Importation in all domains //
68  //-------------------------------//
69 
70 namespace Lorene {
71 void Scalar::import(const Scalar& ci) {
72 
73  int nz = mp->get_mg()->get_nzone() ;
74 
75  import(nz, ci) ;
76 
77 }
78 
79  //--------------------------------------//
80  // Importation in inner domains only //
81  //--------------------------------------//
82 
83 void Scalar::import(int nzet, const Scalar& cm_d) {
84 
85  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
86 
87  // Trivial case : mappings identical !
88  // -----------------------------------
89 
90  if (mp_d == mp) {
91  *this = cm_d ;
92  return ;
93  }
94 
95  // Relative orientation of the two mappings
96  // ----------------------------------------
97 
98  int align_rel = (mp->get_bvect_cart()).get_align()
99  * (mp_d->get_bvect_cart()).get_align() ;
100 
101  switch (align_rel) {
102 
103  case 1 : { // the two mappings have aligned Cartesian axis
104  import_align(nzet, cm_d) ;
105  break ;
106  }
107 
108  case -1 : { // the two mappings have anti-aligned Cartesian axis
109  import_anti(nzet, cm_d) ;
110  break ;
111  }
112 
113  case 0 : { // general case
114  import_gal(nzet, cm_d) ;
115  break ;
116  }
117 
118  default : {
119  cout << "Scalar::import : unexpected value of align_rel : "
120  << align_rel << endl ;
121  abort() ;
122  break ;
123  }
124 
125  }
126 
127 }
128 
129  //--------------------------------------//
130  // General case (axis not aligned) //
131  //--------------------------------------//
132 
133 
134 void Scalar::import_gal(int nzet, const Scalar& cm_d) {
135 
136  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
137 
138  // Trivial case : mappings identical !
139  // -----------------------------------
140 
141  if (mp_d == mp) {
142  *this = cm_d ;
143  return ;
144  }
145 
146  // Another trivial case : null Scalar
147  // -------------------------------
148 
149  if (cm_d.get_etat() == ETATZERO) {
150  set_etat_zero() ;
151  return ;
152  }
153 
154  if (cm_d.get_etat() == ETATUN) {
155  set_etat_one() ;
156  return ;
157  }
158 
159  // Protections
160  // -----------
161 
162  assert(cm_d.get_etat() != ETATNONDEF) ;
163 
164  if (cm_d.get_dzpuis() != 0) {
165  cout <<
166  "Scalar::import : the dzpuis of the Scalar to be imported must be zero !"
167  << endl ;
168  abort() ;
169  }
170 
171 
172  const Mg3d* mg_a = mp->get_mg() ;
173  int nz_a = mg_a->get_nzone() ;
174  assert(nzet <= nz_a) ;
175 
176 
177  // General case :
178  // -------------
179  assert(cm_d.get_etat() == ETATQCQ) ;
180  const Valeur& va_d = cm_d.get_spectral_va() ;
181  va_d.coef() ; // The coefficients are required
182 
183 
184  // Preparations for storing the result in *this
185  // --------------------------------------------
186  del_t() ; // delete all previously computed derived quantities
187 
188  set_etat_qcq() ; // Set the state to ETATQCQ
189 
190  va.set_etat_c_qcq() ; // Allocates the memory for the Mtbl va.c
191  // if it does not exist already
192  va.c->set_etat_qcq() ; // Allocates the memory for the Tbl's in each
193  // domain if they do not exist already
194 
195 
196  // Absolute coordinates of the origin of the Departure mapping
197  double xo_d = mp_d->get_ori_x() ;
198  double yo_d = mp_d->get_ori_y() ;
199  double zo_d = mp_d->get_ori_z() ;
200 
201  // Orientation relative to the Absolute frame of the Departure mapping
202  double rot_phi_d = mp_d->get_rot_phi() ;
203 
204  // Orientation relative to the Absolute frame of the Arrival mapping
205  double rot_phi_a = mp->get_rot_phi() ;
206 
207  // r, theta, phi, X, Y and Z on the Arrival mapping
208  // update of the corresponding Coord's if necessary
209 
210  if ( (mp->r).c == 0x0 ) (mp->r).fait() ;
211  if ( (mp->tet).c == 0x0 ) (mp->tet).fait() ;
212  if ( (mp->phi).c == 0x0 ) (mp->phi).fait() ;
213  if ( (mp->xa).c == 0x0 ) (mp->xa).fait() ;
214  if ( (mp->ya).c == 0x0 ) (mp->ya).fait() ;
215  if ( (mp->za).c == 0x0 ) (mp->za).fait() ;
216 
217  const Mtbl* mr_a = (mp->r).c ;
218  const Mtbl* mtet_a = (mp->tet).c ;
219  const Mtbl* mphi_a = (mp->phi).c ;
220  const Mtbl* mxa_a = (mp->xa).c ;
221  const Mtbl* mya_a = (mp->ya).c ;
222  const Mtbl* mza_a = (mp->za).c ;
223 
224  Param par_precis ; // Required precision in the method Map::val_lx
225  int nitermax = 100 ; // Maximum number of iteration in the secant method
226  int niter ;
227  double precis = 1e-15 ; // Absolute precision in the secant method
228  par_precis.add_int(nitermax) ;
229  par_precis.add_int_mod(niter) ;
230  par_precis.add_double(precis) ;
231 
232 
233  // Loop of the Arrival domains where the computation is to be performed
234  // --------------------------------------------------------------------
235 
236  for (int l=0; l < nzet; l++) {
237 
238  int nr = mg_a->get_nr(l) ;
239  int nt = mg_a->get_nt(l) ;
240  int np = mg_a->get_np(l) ;
241 
242 
243  const double* pr_a = mr_a->t[l]->t ; // Pointer on the values of r
244  const double* ptet_a = mtet_a->t[l]->t ; // Pointer on the values of theta
245  const double* pphi_a = mphi_a->t[l]->t ; // Pointer on the values of phi
246  const double* pxa_a = mxa_a->t[l]->t ; // Pointer on the values of X
247  const double* pya_a = mya_a->t[l]->t ; // Pointer on the values of Y
248  const double* pza_a = mza_a->t[l]->t ; // Pointer on the values of Z
249 
250  (va.c->t[l])->set_etat_qcq() ; // Allocates the array of double to
251  // store the result
252  double* ptx = (va.c->t[l])->t ; // Pointer on the allocated array
253 
254 
255  // Loop on all the grid points in the considered arrival domain:
256 
257  for (int k=0 ; k<np ; k++) {
258  for (int j=0 ; j<nt ; j++) {
259  for (int i=0 ; i<nr ; i++) {
260 
261  double r = *pr_a ;
262  double rd, tetd, phid ;
263  if (r == __infinity) {
264  rd = r ;
265  tetd = *ptet_a ;
266  phid = *pphi_a + rot_phi_a - rot_phi_d ;
267  if (phid < 0) phid += 2*M_PI ;
268  }
269  else {
270  // Coordinates in a Cartesian frame centered on
271  // the Departure mapping and whose axes are
272  // parallel to those of the Absolue Frame
273  double xd = *pxa_a - xo_d ;
274  double yd = *pya_a - yo_d ;
275  double zd = *pza_a - zo_d ;
276 
277  // Spherical coordinates on the Departure mapping
278  double rhod2 = xd*xd + yd*yd ;
279  double rhod = sqrt( rhod2 ) ;
280  rd = sqrt(rhod2 + zd*zd) ;
281  tetd = atan2(rhod, zd) ;
282  phid = atan2(yd, xd) - rot_phi_d ; // (rotation)
283  if (phid < 0) phid += 2*M_PI ;
284  }
285 
286 
287  // NB: to increase the efficiency, the method Scalar::val_point
288  // is not invoked; the method Mtbl_cf::val_point is
289  // called directly instead.
290 
291  // Value of the grid coordinates (l,xi) corresponding to
292  // (rd,tetd,phid) :
293 
294  int ld ; // domain index
295  double xxd ; // radial coordinate xi in [0,1] or [-1,1]
296  mp_d->val_lx(rd, tetd, phid, par_precis, ld, xxd) ;
297 
298  // Value of the Departure Scalar at the obtained point:
299  *ptx = va_d.c_cf->val_point(ld, xxd, tetd, phid) ;
300 
301  // Next point :
302  ptx++ ;
303  pr_a++ ;
304  ptet_a++ ;
305  pphi_a++ ;
306  pxa_a++ ;
307  pya_a++ ;
308  pza_a++ ;
309 
310  }
311  }
312  }
313 
314 
315  } // End of the loop on the Arrival domains
316 
317  // In the remaining domains, *this is set to zero:
318  // ----------------------------------------------
319 
320  if (nzet < nz_a) {
321  annule(nzet, nz_a - 1) ;
322  }
323 
324  // Treatment of dzpuis
325  // -------------------
326 
327 
328  set_dzpuis(0) ;
329 
330 }
331 
332 
333  //-----------------------------------------//
334  // Case of Cartesian axis anti-aligned //
335  //-----------------------------------------//
336 
337 
338 void Scalar::import_anti(int nzet, const Scalar& cm_d) {
339 
340  // Trivial case : null Scalar
341  // ------------------------
342 
343  if (cm_d.get_etat() == ETATZERO) {
344  set_etat_zero() ;
345  return ;
346  }
347 
348  if (cm_d.get_etat() == ETATUN) {
349  set_etat_one() ;
350  return ;
351  }
352 
353  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
354 
355  // Protections
356  // -----------
357  int align = (mp->get_bvect_cart()).get_align() ;
358 
359  assert( align * (mp_d->get_bvect_cart()).get_align() == -1 ) ;
360 
361  assert(cm_d.get_etat() == ETATQCQ) ;
362 
363  if (cm_d.get_dzpuis() != 0) {
364  cout <<
365  "Scalar::import : the dzpuis of the Scalar to be imported must be zero !"
366  << endl ;
367  abort() ;
368  }
369 
370 
371  const Mg3d* mg_a = mp->get_mg() ;
372  int nz_a = mg_a->get_nzone() ;
373  assert(nzet <= nz_a) ;
374 
375  const Valeur& va_d = cm_d.get_spectral_va() ;
376  va_d.coef() ; // The coefficients are required
377 
378 
379  // Preparations for storing the result in *this
380  // --------------------------------------------
381  del_t() ; // delete all previously computed derived quantities
382 
383  set_etat_qcq() ; // Set the state to ETATQCQ
384 
385  va.set_etat_c_qcq() ; // Allocates the memory for the Mtbl va.c
386  // if it does not exist already
387  va.c->set_etat_qcq() ; // Allocates the memory for the Tbl's in each
388  // domain if they do not exist already
389 
390 
391  // Departure (x,y,z) coordinates of the origin of the Arrival mapping :
392 
393  double xx_a, yy_a, zz_a ;
394  if (align == 1) {
395  xx_a = mp_d->get_ori_x() - mp->get_ori_x() ;
396  yy_a = mp_d->get_ori_y() - mp->get_ori_y() ;
397  }
398  else {
399  xx_a = mp->get_ori_x() - mp_d->get_ori_x() ;
400  yy_a = mp->get_ori_y() - mp_d->get_ori_y() ;
401  }
402  zz_a = mp->get_ori_z() - mp_d->get_ori_z() ;
403 
404 
405  // r, theta, phi, x, y and z on the Arrival mapping
406  // update of the corresponding Coord's if necessary
407 
408  if ( (mp->r).c == 0x0 ) (mp->r).fait() ;
409  if ( (mp->tet).c == 0x0 ) (mp->tet).fait() ;
410  if ( (mp->phi).c == 0x0 ) (mp->phi).fait() ;
411  if ( (mp->x).c == 0x0 ) (mp->x).fait() ;
412  if ( (mp->y).c == 0x0 ) (mp->y).fait() ;
413  if ( (mp->z).c == 0x0 ) (mp->z).fait() ;
414 
415  const Mtbl* mr_a = (mp->r).c ;
416  const Mtbl* mtet_a = (mp->tet).c ;
417  const Mtbl* mphi_a = (mp->phi).c ;
418  const Mtbl* mx_a = (mp->x).c ;
419  const Mtbl* my_a = (mp->y).c ;
420  const Mtbl* mz_a = (mp->z).c ;
421 
422  Param par_precis ; // Required precision in the method Map::val_lx
423  int nitermax = 100 ; // Maximum number of iteration in the secant method
424  int niter ;
425  double precis = 1e-15 ; // Absolute precision in the secant method
426  par_precis.add_int(nitermax) ;
427  par_precis.add_int_mod(niter) ;
428  par_precis.add_double(precis) ;
429 
430 
431  // Loop of the Arrival domains where the computation is to be performed
432  // --------------------------------------------------------------------
433 
434  for (int l=0; l < nzet; l++) {
435 
436  int nr = mg_a->get_nr(l) ;
437  int nt = mg_a->get_nt(l) ;
438  int np = mg_a->get_np(l) ;
439 
440 
441  const double* pr_a = mr_a->t[l]->t ; // Pointer on the values of r
442  const double* ptet_a = mtet_a->t[l]->t ; // Pointer on the values of theta
443  const double* pphi_a = mphi_a->t[l]->t ; // Pointer on the values of phi
444  const double* px_a = mx_a->t[l]->t ; // Pointer on the values of X
445  const double* py_a = my_a->t[l]->t ; // Pointer on the values of Y
446  const double* pz_a = mz_a->t[l]->t ; // Pointer on the values of Z
447 
448  (va.c->t[l])->set_etat_qcq() ; // Allocates the array of double to
449  // store the result
450  double* ptx = (va.c->t[l])->t ; // Pointer on the allocated array
451 
452 
453  // Loop on all the grid points in the considered arrival domain:
454 
455  for (int k=0 ; k<np ; k++) {
456  for (int j=0 ; j<nt ; j++) {
457  for (int i=0 ; i<nr ; i++) {
458 
459  double r = *pr_a ;
460  double rd, tetd, phid ;
461  if (r == __infinity) {
462  rd = r ;
463  tetd = *ptet_a ;
464  phid = *pphi_a + M_PI ;
465  if (phid < 0) phid += 2*M_PI ;
466  }
467  else {
468 
469  // Cartesian coordinates on the Departure mapping
470  double xd = - *px_a + xx_a ;
471  double yd = - *py_a + yy_a ;
472  double zd = *pz_a + zz_a ;
473 
474  // Spherical coordinates on the Departure mapping
475  double rhod2 = xd*xd + yd*yd ;
476  double rhod = sqrt( rhod2 ) ;
477  rd = sqrt(rhod2 + zd*zd) ;
478  tetd = atan2(rhod, zd) ;
479  phid = atan2(yd, xd) ;
480  if (phid < 0) phid += 2*M_PI ;
481  }
482 
483 
484  // NB: to increase the efficiency, the method Scalar::val_point
485  // is not invoked; the method Mtbl_cf::val_point is
486  // called directly instead.
487 
488  // Value of the grid coordinates (l,xi) corresponding to
489  // (rd,tetd,phid) :
490 
491  int ld ; // domain index
492  double xxd ; // radial coordinate xi in [0,1] or [-1,1]
493  mp_d->val_lx(rd, tetd, phid, par_precis, ld, xxd) ;
494 
495  // Value of the Departure Scalar at the obtained point:
496  *ptx = va_d.c_cf->val_point(ld, xxd, tetd, phid) ;
497 
498  // Next point :
499  ptx++ ;
500  pr_a++ ;
501  ptet_a++ ;
502  pphi_a++ ;
503  px_a++ ;
504  py_a++ ;
505  pz_a++ ;
506 
507  }
508  }
509  }
510 
511 
512  } // End of the loop on the Arrival domains
513 
514  // In the remaining domains, *this is set to zero:
515  // ----------------------------------------------
516 
517  if (nzet < nz_a) {
518  annule(nzet, nz_a - 1) ;
519  }
520 
521  // Treatment of dzpuis
522  // -------------------
523 
524 
525  set_dzpuis(0) ;
526 
527 }
528 
529 
530  //-------------------------------------//
531  // Case of aligned Cartesian axis //
532  //-------------------------------------//
533 
534 
535 void Scalar::import_align(int nzet, const Scalar& cm_d) {
536 
537  // Trivial case : null Scalar
538  // ------------------------
539 
540  if (cm_d.get_etat() == ETATZERO) {
541  set_etat_zero() ;
542  return ;
543  }
544  if (cm_d.get_etat() == ETATUN) {
545  set_etat_one() ;
546  return ;
547  }
548 
549  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
550 
551  // Protections
552  // -----------
553  int align = (mp->get_bvect_cart()).get_align() ;
554 
555  assert( align * (mp_d->get_bvect_cart()).get_align() == 1 ) ;
556 
557  assert(cm_d.get_etat() == ETATQCQ) ;
558 
559  if (cm_d.get_dzpuis() != 0) {
560  cout <<
561  "Scalar::import : the dzpuis of the Scalar to be imported must be zero !"
562  << endl ;
563  abort() ;
564  }
565 
566 
567  const Mg3d* mg_a = mp->get_mg() ;
568  int nz_a = mg_a->get_nzone() ;
569  assert(nzet <= nz_a) ;
570 
571  const Valeur& va_d = cm_d.get_spectral_va() ;
572  va_d.coef() ; // The coefficients are required
573 
574 
575  // Preparations for storing the result in *this
576  // --------------------------------------------
577  del_t() ; // delete all previously computed derived quantities
578 
579  set_etat_qcq() ; // Set the state to ETATQCQ
580 
581  va.set_etat_c_qcq() ; // Allocates the memory for the Mtbl va.c
582  // if it does not exist already
583  va.c->set_etat_qcq() ; // Allocates the memory for the Tbl's in each
584  // domain if they do not exist already
585 
586 
587  // Departure (x,y,z) coordinates of the origin of the Arrival mapping :
588 
589  double xx_a, yy_a, zz_a ;
590  if (align == 1) {
591  xx_a = mp->get_ori_x() - mp_d->get_ori_x() ;
592  yy_a = mp->get_ori_y() - mp_d->get_ori_y() ;
593  }
594  else {
595  xx_a = mp_d->get_ori_x() - mp->get_ori_x() ;
596  yy_a = mp_d->get_ori_y() - mp->get_ori_y() ;
597  }
598  zz_a = mp->get_ori_z() - mp_d->get_ori_z() ;
599 
600 
601  // r, theta, phi, x, y and z on the Arrival mapping
602  // update of the corresponding Coord's if necessary
603 
604  if ( (mp->r).c == 0x0 ) (mp->r).fait() ;
605  if ( (mp->tet).c == 0x0 ) (mp->tet).fait() ;
606  if ( (mp->phi).c == 0x0 ) (mp->phi).fait() ;
607  if ( (mp->x).c == 0x0 ) (mp->x).fait() ;
608  if ( (mp->y).c == 0x0 ) (mp->y).fait() ;
609  if ( (mp->z).c == 0x0 ) (mp->z).fait() ;
610 
611  const Mtbl* mr_a = (mp->r).c ;
612  const Mtbl* mtet_a = (mp->tet).c ;
613  const Mtbl* mphi_a = (mp->phi).c ;
614  const Mtbl* mx_a = (mp->x).c ;
615  const Mtbl* my_a = (mp->y).c ;
616  const Mtbl* mz_a = (mp->z).c ;
617 
618  Param par_precis ; // Required precision in the method Map::val_lx
619  int nitermax = 100 ; // Maximum number of iteration in the secant method
620  int niter ;
621  double precis = 1e-15 ; // Absolute precision in the secant method
622  par_precis.add_int(nitermax) ;
623  par_precis.add_int_mod(niter) ;
624  par_precis.add_double(precis) ;
625 
626 
627  // Loop of the Arrival domains where the computation is to be performed
628  // --------------------------------------------------------------------
629 
630  for (int l=0; l < nzet; l++) {
631 
632  int nr = mg_a->get_nr(l) ;
633  int nt = mg_a->get_nt(l) ;
634  int np = mg_a->get_np(l) ;
635 
636 
637  const double* pr_a = mr_a->t[l]->t ; // Pointer on the values of r
638  const double* ptet_a = mtet_a->t[l]->t ; // Pointer on the values of theta
639  const double* pphi_a = mphi_a->t[l]->t ; // Pointer on the values of phi
640  const double* px_a = mx_a->t[l]->t ; // Pointer on the values of X
641  const double* py_a = my_a->t[l]->t ; // Pointer on the values of Y
642  const double* pz_a = mz_a->t[l]->t ; // Pointer on the values of Z
643 
644  (va.c->t[l])->set_etat_qcq() ; // Allocates the array of double to
645  // store the result
646  double* ptx = (va.c->t[l])->t ; // Pointer on the allocated array
647 
648 
649  // Loop on all the grid points in the considered arrival domain:
650 
651  for (int k=0 ; k<np ; k++) {
652  for (int j=0 ; j<nt ; j++) {
653  for (int i=0 ; i<nr ; i++) {
654 
655  double r = *pr_a ;
656  double rd, tetd, phid ;
657  if (r == __infinity) {
658  rd = r ;
659  tetd = *ptet_a ;
660  phid = *pphi_a ;
661  }
662  else {
663 
664  // Cartesian coordinates on the Departure mapping
665  double xd = *px_a + xx_a ;
666  double yd = *py_a + yy_a ;
667  double zd = *pz_a + zz_a ;
668 
669  // Spherical coordinates on the Departure mapping
670  double rhod2 = xd*xd + yd*yd ;
671  double rhod = sqrt( rhod2 ) ;
672  rd = sqrt(rhod2 + zd*zd) ;
673  tetd = atan2(rhod, zd) ;
674  phid = atan2(yd, xd) ;
675  if (phid < 0) phid += 2*M_PI ;
676  }
677 
678 
679  // NB: to increase the efficiency, the method Scalar::val_point
680  // is not invoked; the method Mtbl_cf::val_point is
681  // called directly instead.
682 
683  // Value of the grid coordinates (l,xi) corresponding to
684  // (rd,tetd,phid) :
685 
686  int ld ; // domain index
687  double xxd ; // radial coordinate xi in [0,1] or [-1,1]
688  mp_d->val_lx(rd, tetd, phid, par_precis, ld, xxd) ;
689 
690  // Value of the Departure Scalar at the obtained point:
691  *ptx = va_d.c_cf->val_point(ld, xxd, tetd, phid) ;
692 
693  // Next point :
694  ptx++ ;
695  pr_a++ ;
696  ptet_a++ ;
697  pphi_a++ ;
698  px_a++ ;
699  py_a++ ;
700  pz_a++ ;
701 
702  }
703  }
704  }
705 
706 
707  } // End of the loop on the Arrival domains
708 
709  // In the remaining domains, *this is set to zero:
710  // ----------------------------------------------
711 
712  if (nzet < nz_a) {
713  annule(nzet, nz_a - 1) ;
714  }
715 
716  // Treatment of dzpuis
717  // -------------------
718 
719 
720  set_dzpuis(0) ;
721 
722 }
723 }
Coord xa
Absolute x coordinate.
Definition: map.h:748
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:312
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Definition: param.C:249
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:151
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:330
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition: scalar.C:397
Multi-domain array.
Definition: mtbl.h:118
double get_ori_y() const
Returns the y coordinate of the origin.
Definition: map.h:788
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
Base class for coordinate mappings.
Definition: map.h:688
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:786
void import_gal(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings do not have a part...
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:560
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition: map.h:793
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: scalar.C:359
friend Scalar sqrt(const Scalar &)
Square root.
Definition: scalar_math.C:266
Coord tet
coordinate centered on the grid
Definition: map.h:737
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:814
Coord phi
coordinate centered on the grid
Definition: map.h:738
void import_anti(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned ...
int get_dzpuis() const
Returns dzpuis.
Definition: scalar.h:563
double * t
The array of double.
Definition: tbl.h:176
Mtbl * c
Values of the function at the points of the multi-grid.
Definition: valeur.h:309
Parameter storage.
Definition: param.h:125
void set_etat_one()
Sets the logical state to ETATUN (one).
Definition: scalar.C:340
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Valeur va
The numerical value of the Scalar.
Definition: scalar.h:411
void del_t()
Logical destructor.
Definition: scalar.C:285
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:302
void import_align(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Carte...
double val_point(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
void import(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
Definition: scalar_import.C:71
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
Coord ya
Absolute y coordinate.
Definition: map.h:749
Multi-domain grid.
Definition: grilles.h:279
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
Coord y
y coordinate centered on the grid
Definition: map.h:745
Coord za
Absolute z coordinate.
Definition: map.h:750
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
Definition: param.C:318
Coord x
x coordinate centered on the grid
Definition: map.h:744
void set_etat_c_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl c )...
Definition: valeur.C:704
double get_ori_z() const
Returns the z coordinate of the origin.
Definition: map.h:790
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:301
Tbl ** t
Array (size nzone ) of pointers on the Tbl &#39;s.
Definition: mtbl.h:132
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:874
Coord z
z coordinate centered on the grid
Definition: map.h:746
void add_int_mod(int &n, int position=0)
Adds the address of a new modifiable int to the list.
Definition: param.C:388
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:607
Coord r
r coordinate centered on the grid
Definition: map.h:736
virtual void val_lx(double rr, double theta, double pphi, int &l, double &xi) const =0
Computes the domain index l and the value of corresponding to a point given by its physical coordina...