LORENE
scalar_import_symy.C
1 /*
2  * Member function of the Scalar class for initiating a Scalar from
3  * a Scalar defined on another mapping.
4  * Case where both Scalar's are symmetric with respect to their y=0 plane.
5  */
6 
7 /*
8  * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
9  * Copyright (c) 1999-2001 Eric Gourgoulhon (Cmp version)
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 scalar_import_symy_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_import_symy.C,v 1.5 2014/10/13 08:53:46 j_novak Exp $" ;
31 
32 
33 /*
34  * $Id: scalar_import_symy.C,v 1.5 2014/10/13 08:53:46 j_novak Exp $
35  * $Log: scalar_import_symy.C,v $
36  * Revision 1.5 2014/10/13 08:53:46 j_novak
37  * Lorene classes and functions now belong to the namespace Lorene.
38  *
39  * Revision 1.4 2014/10/06 15:16:15 j_novak
40  * Modified #include directives to use c++ syntax.
41  *
42  * Revision 1.3 2003/10/10 15:57:29 j_novak
43  * Added the state one (ETATUN) to the class Scalar
44  *
45  * Revision 1.2 2003/10/01 13:04:44 e_gourgoulhon
46  * The method Tensor::get_mp() returns now a reference (and not
47  * a pointer) onto a mapping.
48  *
49  * Revision 1.1 2003/09/25 09:07:05 j_novak
50  * Added the functions for importing from another mapping (to be tested).
51  *
52  *
53  * $Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_import_symy.C,v 1.5 2014/10/13 08:53:46 j_novak Exp $
54  *
55  */
56 
57 
58 
59 // Headers C
60 #include <cmath>
61 
62 // Headers Lorene
63 #include "tensor.h"
64 #include "param.h"
65 #include "nbr_spx.h"
66 
67  //-------------------------------//
68  // Importation in all domains //
69  //-------------------------------//
70 
71 namespace Lorene {
72 void Scalar::import_symy(const Scalar& ci) {
73 
74  int nz = mp->get_mg()->get_nzone() ;
75 
76  import_symy(nz, ci) ;
77 
78 }
79 
80  //--------------------------------------//
81  // Importation in inner domains only //
82  //--------------------------------------//
83 
84 void Scalar::import_symy(int nzet, const Scalar& cm_d) {
85 
86  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
87 
88  // Trivial case : mappings identical !
89  // -----------------------------------
90 
91  if (mp_d == mp) {
92  *this = cm_d ;
93  return ;
94  }
95 
96  // Relative orientation of the two mappings
97  // ----------------------------------------
98 
99  int align_rel = (mp->get_bvect_cart()).get_align()
100  * (mp_d->get_bvect_cart()).get_align() ;
101 
102  switch (align_rel) {
103 
104  case 1 : { // the two mappings have aligned Cartesian axis
105  import_align_symy(nzet, cm_d) ;
106  break ;
107  }
108 
109  case -1 : { // the two mappings have anti-aligned Cartesian axis
110  import_anti_symy(nzet, cm_d) ;
111  break ;
112  }
113 
114  default : {
115  cout << "Scalar::import_symy : unexpected value of align_rel : "
116  << align_rel << endl ;
117  abort() ;
118  break ;
119  }
120 
121  }
122 
123 }
124 
125 
126  //-----------------------------------------//
127  // Case of Cartesian axis anti-aligned //
128  //-----------------------------------------//
129 
130 
131 void Scalar::import_anti_symy(int nzet, const Scalar& cm_d) {
132 
133  // Trivial case : null Scalar
134  // ------------------------
135 
136  if (cm_d.get_etat() == ETATZERO) {
137  set_etat_zero() ;
138  return ;
139  }
140  if (cm_d.get_etat() == ETATUN) {
141  set_etat_one() ;
142  return ;
143  }
144 
145  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
146 
147  // Protections
148  // -----------
149  int align = (mp->get_bvect_cart()).get_align() ;
150 
151  assert( align * (mp_d->get_bvect_cart()).get_align() == -1 ) ;
152 
153  assert(cm_d.get_etat() == ETATQCQ) ;
154 
155  if (cm_d.get_dzpuis() != 0) {
156  cout <<
157  "Scalar::import_anti_symy : the dzpuis of the Scalar to be imported"
158  << " must be zero !" << endl ;
159  abort() ;
160  }
161 
162 
163  const Mg3d* mg_a = mp->get_mg() ;
164  assert(mg_a->get_type_p() == NONSYM) ;
165 
166  int nz_a = mg_a->get_nzone() ;
167  assert(nzet <= nz_a) ;
168 
169  const Valeur& va_d = cm_d.get_spectral_va() ;
170  va_d.coef() ; // The coefficients are required
171 
172 
173  // Preparations for storing the result in *this
174  // --------------------------------------------
175  del_t() ; // delete all previously computed derived quantities
176 
177  set_etat_qcq() ; // Set the state to ETATQCQ
178 
179  va.set_etat_c_qcq() ; // Allocates the memory for the Mtbl va.c
180  // if it does not exist already
181  va.c->set_etat_qcq() ; // Allocates the memory for the Tbl's in each
182  // domain if they do not exist already
183 
184 
185  // Departure (x,y,z) coordinates of the origin of the Arrival mapping :
186 
187  double xx_a, yy_a, zz_a ;
188  if (align == 1) {
189  xx_a = mp_d->get_ori_x() - mp->get_ori_x() ;
190  yy_a = mp_d->get_ori_y() - mp->get_ori_y() ;
191  }
192  else {
193  xx_a = mp->get_ori_x() - mp_d->get_ori_x() ;
194  yy_a = mp->get_ori_y() - mp_d->get_ori_y() ;
195  }
196  zz_a = mp->get_ori_z() - mp_d->get_ori_z() ;
197 
198 
199  // r, theta, phi, x, y and z on the Arrival mapping
200  // update of the corresponding Coord's if necessary
201 
202  if ( (mp->r).c == 0x0 ) (mp->r).fait() ;
203  if ( (mp->tet).c == 0x0 ) (mp->tet).fait() ;
204  if ( (mp->phi).c == 0x0 ) (mp->phi).fait() ;
205  if ( (mp->x).c == 0x0 ) (mp->x).fait() ;
206  if ( (mp->y).c == 0x0 ) (mp->y).fait() ;
207  if ( (mp->z).c == 0x0 ) (mp->z).fait() ;
208 
209  const Mtbl* mr_a = (mp->r).c ;
210  const Mtbl* mtet_a = (mp->tet).c ;
211  const Mtbl* mphi_a = (mp->phi).c ;
212  const Mtbl* mx_a = (mp->x).c ;
213  const Mtbl* my_a = (mp->y).c ;
214  const Mtbl* mz_a = (mp->z).c ;
215 
216  Param par_precis ; // Required precision in the method Map::val_lx
217  int nitermax = 100 ; // Maximum number of iteration in the secant method
218  int niter ;
219  double precis = 1e-15 ; // Absolute precision in the secant method
220  par_precis.add_int(nitermax) ;
221  par_precis.add_int_mod(niter) ;
222  par_precis.add_double(precis) ;
223 
224 
225  // Loop of the Arrival domains where the computation is to be performed
226  // --------------------------------------------------------------------
227 
228  for (int l=0; l < nzet; l++) {
229 
230  int nr = mg_a->get_nr(l) ;
231  int nt = mg_a->get_nt(l) ;
232  int np = mg_a->get_np(l) ;
233 
234 
235  const double* pr_a = mr_a->t[l]->t ; // Pointer on the values of r
236  const double* ptet_a = mtet_a->t[l]->t ; // Pointer on the values of theta
237  const double* pphi_a = mphi_a->t[l]->t ; // Pointer on the values of phi
238  const double* px_a = mx_a->t[l]->t ; // Pointer on the values of X
239  const double* py_a = my_a->t[l]->t ; // Pointer on the values of Y
240  const double* pz_a = mz_a->t[l]->t ; // Pointer on the values of Z
241 
242  (va.c->t[l])->set_etat_qcq() ; // Allocates the array of double to
243  // store the result
244  double* ptx = (va.c->t[l])->t ; // Pointer on the allocated array
245 
246 
247  // Loop on half the grid points in the considered arrival domain
248  // (the other half will be obtained by symmetry with respect to
249  // the y=0 plane).
250 
251  for (int k=0 ; k<np/2+1 ; k++) { // np/2+1 : half the grid
252  for (int j=0 ; j<nt ; j++) {
253  for (int i=0 ; i<nr ; i++) {
254 
255  double r = *pr_a ;
256  double rd, tetd, phid ;
257  if (r == __infinity) {
258  rd = r ;
259  tetd = *ptet_a ;
260  phid = *pphi_a + M_PI ;
261  if (phid < 0) phid += 2*M_PI ;
262  }
263  else {
264 
265  // Cartesian coordinates on the Departure mapping
266  double xd = - *px_a + xx_a ;
267  double yd = - *py_a + yy_a ;
268  double zd = *pz_a + zz_a ;
269 
270  // Spherical coordinates on the Departure mapping
271  double rhod2 = xd*xd + yd*yd ;
272  double rhod = sqrt( rhod2 ) ;
273  rd = sqrt(rhod2 + zd*zd) ;
274  tetd = atan2(rhod, zd) ;
275  phid = atan2(yd, xd) ;
276  if (phid < 0) phid += 2*M_PI ;
277  }
278 
279 
280  // NB: to increase the efficiency, the method Scalar::val_point
281  // is not invoked; the method Mtbl_cf::val_point is
282  // called directly instead.
283 
284  // Value of the grid coordinates (l,xi) corresponding to
285  // (rd,tetd,phid) :
286 
287  int ld ; // domain index
288  double xxd ; // radial coordinate xi in [0,1] or [-1,1]
289  mp_d->val_lx(rd, tetd, phid, par_precis, ld, xxd) ;
290 
291  // Value of the Departure Scalar at the obtained point:
292  *ptx = va_d.c_cf->val_point_symy(ld, xxd, tetd, phid) ;
293 
294  // Next point :
295  ptx++ ;
296  pr_a++ ;
297  ptet_a++ ;
298  pphi_a++ ;
299  px_a++ ;
300  py_a++ ;
301  pz_a++ ;
302 
303  }
304  }
305  }
306 
307  // The remaining points are obtained by symmetry with rspect to the
308  // y=0 plane
309 
310  for (int k=np/2+1 ; k<np ; k++) {
311 
312  // pointer on the value (already computed) at the point symmetric
313  // with respect to the plane y=0
314  double* ptx_symy = (va.c->t[l])->t + (np-k)*nt*nr ;
315 
316  // copy :
317  for (int j=0 ; j<nt ; j++) {
318  for (int i=0 ; i<nr ; i++) {
319  *ptx = *ptx_symy ;
320  ptx++ ;
321  ptx_symy++ ;
322  }
323  }
324  }
325 
326 
327  } // End of the loop on the Arrival domains
328 
329  // In the remaining domains, *this is set to zero:
330  // ----------------------------------------------
331 
332  if (nzet < nz_a) {
333  annule(nzet, nz_a - 1) ;
334  }
335 
336  // Treatment of dzpuis
337  // -------------------
338 
339  set_dzpuis(0) ;
340 
341 }
342 
343 
344  //-------------------------------------//
345  // Case of aligned Cartesian axis //
346  //-------------------------------------//
347 
348 
349 void Scalar::import_align_symy(int nzet, const Scalar& cm_d) {
350 
351  // Trivial case : null Scalar
352  // ------------------------
353 
354  if (cm_d.get_etat() == ETATZERO) {
355  set_etat_zero() ;
356  return ;
357  }
358  if (cm_d.get_etat() == ETATUN) {
359  set_etat_one() ;
360  return ;
361  }
362 
363  const Map* mp_d = &(cm_d.get_mp()) ; // Departure mapping
364 
365  // Protections
366  // -----------
367  int align = (mp->get_bvect_cart()).get_align() ;
368 
369  assert( align * (mp_d->get_bvect_cart()).get_align() == 1 ) ;
370 
371  assert(cm_d.get_etat() == ETATQCQ) ;
372 
373  if (cm_d.get_dzpuis() != 0) {
374  cout <<
375  "Scalar::import_align_symy : the dzpuis of the Scalar to be imported"
376  << " must be zero !" << endl ;
377  abort() ;
378  }
379 
380 
381  const Mg3d* mg_a = mp->get_mg() ;
382  assert(mg_a->get_type_p() == NONSYM) ;
383 
384  int nz_a = mg_a->get_nzone() ;
385  assert(nzet <= nz_a) ;
386 
387  const Valeur& va_d = cm_d.get_spectral_va() ;
388  va_d.coef() ; // The coefficients are required
389 
390 
391  // Preparations for storing the result in *this
392  // --------------------------------------------
393  del_t() ; // delete all previously computed derived quantities
394 
395  set_etat_qcq() ; // Set the state to ETATQCQ
396 
397  va.set_etat_c_qcq() ; // Allocates the memory for the Mtbl va.c
398  // if it does not exist already
399  va.c->set_etat_qcq() ; // Allocates the memory for the Tbl's in each
400  // domain if they do not exist already
401 
402 
403  // Departure (x,y,z) coordinates of the origin of the Arrival mapping :
404 
405  double xx_a, yy_a, zz_a ;
406  if (align == 1) {
407  xx_a = mp->get_ori_x() - mp_d->get_ori_x() ;
408  yy_a = mp->get_ori_y() - mp_d->get_ori_y() ;
409  }
410  else {
411  xx_a = mp_d->get_ori_x() - mp->get_ori_x() ;
412  yy_a = mp_d->get_ori_y() - mp->get_ori_y() ;
413  }
414  zz_a = mp->get_ori_z() - mp_d->get_ori_z() ;
415 
416 
417  // r, theta, phi, x, y and z on the Arrival mapping
418  // update of the corresponding Coord's if necessary
419 
420  if ( (mp->r).c == 0x0 ) (mp->r).fait() ;
421  if ( (mp->tet).c == 0x0 ) (mp->tet).fait() ;
422  if ( (mp->phi).c == 0x0 ) (mp->phi).fait() ;
423  if ( (mp->x).c == 0x0 ) (mp->x).fait() ;
424  if ( (mp->y).c == 0x0 ) (mp->y).fait() ;
425  if ( (mp->z).c == 0x0 ) (mp->z).fait() ;
426 
427  const Mtbl* mr_a = (mp->r).c ;
428  const Mtbl* mtet_a = (mp->tet).c ;
429  const Mtbl* mphi_a = (mp->phi).c ;
430  const Mtbl* mx_a = (mp->x).c ;
431  const Mtbl* my_a = (mp->y).c ;
432  const Mtbl* mz_a = (mp->z).c ;
433 
434  Param par_precis ; // Required precision in the method Map::val_lx
435  int nitermax = 100 ; // Maximum number of iteration in the secant method
436  int niter ;
437  double precis = 1e-15 ; // Absolute precision in the secant method
438  par_precis.add_int(nitermax) ;
439  par_precis.add_int_mod(niter) ;
440  par_precis.add_double(precis) ;
441 
442 
443  // Loop of the Arrival domains where the computation is to be performed
444  // --------------------------------------------------------------------
445 
446  for (int l=0; l < nzet; l++) {
447 
448  int nr = mg_a->get_nr(l) ;
449  int nt = mg_a->get_nt(l) ;
450  int np = mg_a->get_np(l) ;
451 
452 
453  const double* pr_a = mr_a->t[l]->t ; // Pointer on the values of r
454  const double* ptet_a = mtet_a->t[l]->t ; // Pointer on the values of theta
455  const double* pphi_a = mphi_a->t[l]->t ; // Pointer on the values of phi
456  const double* px_a = mx_a->t[l]->t ; // Pointer on the values of X
457  const double* py_a = my_a->t[l]->t ; // Pointer on the values of Y
458  const double* pz_a = mz_a->t[l]->t ; // Pointer on the values of Z
459 
460  (va.c->t[l])->set_etat_qcq() ; // Allocates the array of double to
461  // store the result
462  double* ptx = (va.c->t[l])->t ; // Pointer on the allocated array
463 
464 
465 
466  // Loop on half the grid points in the considered arrival domain
467  // (the other half will be obtained by symmetry with respect to
468  // the y=0 plane).
469 
470  for (int k=0 ; k<np/2+1 ; k++) { // np/2+1 : half the grid
471  for (int j=0 ; j<nt ; j++) {
472  for (int i=0 ; i<nr ; i++) {
473 
474  double r = *pr_a ;
475  double rd, tetd, phid ;
476  if (r == __infinity) {
477  rd = r ;
478  tetd = *ptet_a ;
479  phid = *pphi_a ;
480  }
481  else {
482 
483  // Cartesian coordinates on the Departure mapping
484  double xd = *px_a + xx_a ;
485  double yd = *py_a + yy_a ;
486  double zd = *pz_a + zz_a ;
487 
488  // Spherical coordinates on the Departure mapping
489  double rhod2 = xd*xd + yd*yd ;
490  double rhod = sqrt( rhod2 ) ;
491  rd = sqrt(rhod2 + zd*zd) ;
492  tetd = atan2(rhod, zd) ;
493  phid = atan2(yd, xd) ;
494  if (phid < 0) phid += 2*M_PI ;
495  }
496 
497 
498  // NB: to increase the efficiency, the method Scalar::val_point
499  // is not invoked; the method Mtbl_cf::val_point is
500  // called directly instead.
501 
502  // Value of the grid coordinates (l,xi) corresponding to
503  // (rd,tetd,phid) :
504 
505  int ld ; // domain index
506  double xxd ; // radial coordinate xi in [0,1] or [-1,1]
507  mp_d->val_lx(rd, tetd, phid, par_precis, ld, xxd) ;
508 
509  // Value of the Departure Scalar at the obtained point:
510  *ptx = va_d.c_cf->val_point_symy(ld, xxd, tetd, phid) ;
511 
512  // Next point :
513  ptx++ ;
514  pr_a++ ;
515  ptet_a++ ;
516  pphi_a++ ;
517  px_a++ ;
518  py_a++ ;
519  pz_a++ ;
520 
521  }
522  }
523  }
524 
525 
526  // The remaining points are obtained by symmetry with rspect to the
527  // y=0 plane
528 
529  for (int k=np/2+1 ; k<np ; k++) {
530 
531  // pointer on the value (already computed) at the point symmetric
532  // with respect to the plane y=0
533  double* ptx_symy = (va.c->t[l])->t + (np-k)*nt*nr ;
534 
535  // copy :
536  for (int j=0 ; j<nt ; j++) {
537  for (int i=0 ; i<nr ; i++) {
538  *ptx = *ptx_symy ;
539  ptx++ ;
540  ptx_symy++ ;
541  }
542  }
543  }
544 
545  } // End of the loop on the Arrival domains
546 
547  // In the remaining domains, *this is set to zero:
548  // ----------------------------------------------
549 
550  if (nzet < nz_a) {
551  annule(nzet, nz_a - 1) ;
552  }
553 
554  // Treatment of dzpuis
555  // -------------------
556 
557  set_dzpuis(0) ;
558 
559 }
560 }
int get_type_p() const
Returns the type of sampling in the direction: SYM : : symmetry with respect to the transformatio...
Definition: grilles.h:495
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:302
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Definition: param.C:246
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:148
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:324
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition: scalar.C:391
Multi-domain array.
Definition: mtbl.h:118
double get_ori_y() const
Returns the y coordinate of the origin.
Definition: map.h:770
Lorene prototypes.
Definition: app_hor.h:64
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Base class for coordinate mappings.
Definition: map.h:670
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:768
void import_anti_symy(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have anti-aligned ...
Values and coefficients of a (real-value) function.
Definition: valeur.h:287
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:554
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: scalar.C:353
friend Scalar sqrt(const Scalar &)
Square root.
Definition: scalar_math.C:263
Coord tet
coordinate centered on the grid
Definition: map.h:719
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:808
Coord phi
coordinate centered on the grid
Definition: map.h:720
int get_dzpuis() const
Returns dzpuis.
Definition: scalar.h:557
double * t
The array of double.
Definition: tbl.h:173
Mtbl * c
Values of the function at the points of the multi-grid.
Definition: valeur.h:299
Parameter storage.
Definition: param.h:125
void set_etat_one()
Sets the logical state to ETATUN (one).
Definition: scalar.C:334
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Valeur va
The numerical value of the Scalar.
Definition: scalar.h:405
void del_t()
Logical destructor.
Definition: scalar.C:279
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:299
void import_align_symy(int nzet, const Scalar &ci)
Assignment to another Scalar defined on a different mapping, when the two mappings have aligned Carte...
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
Multi-domain grid.
Definition: grilles.h:273
double val_point_symy(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 ...
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:791
Coord y
y coordinate centered on the grid
Definition: map.h:727
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
Definition: param.C:315
Coord x
x coordinate centered on the grid
Definition: map.h:726
void import_symy(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
void set_etat_c_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl c )...
Definition: valeur.C:701
double get_ori_z() const
Returns the z coordinate of the origin.
Definition: map.h:772
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:295
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:861
Coord z
z coordinate centered on the grid
Definition: map.h:728
void add_int_mod(int &n, int position=0)
Adds the address of a new modifiable int to the list.
Definition: param.C:385
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:601
Coord r
r coordinate centered on the grid
Definition: map.h:718
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...