LORENE
map_af_radius.C
1 /*
2  * Methods of the class Map_af relative to the function
3  * r = R_l(xi, theta', phi')
4  */
5 
6 /*
7  * Copyright (c) 1999-2001 Eric Gourgoulhon
8  *
9  * This file is part of LORENE.
10  *
11  * LORENE is free software; you can redistribute it and/or modify
12  * it under the terms of the GNU General Public License as published by
13  * the Free Software Foundation; either version 2 of the License, or
14  * (at your option) any later version.
15  *
16  * LORENE is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19  * GNU General Public License for more details.
20  *
21  * You should have received a copy of the GNU General Public License
22  * along with LORENE; if not, write to the Free Software
23  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
24  *
25  */
26 
27 
28 char map_af_radius_C[] = "$Header: /cvsroot/Lorene/C++/Source/Map/map_af_radius.C,v 1.9 2014/10/13 08:53:03 j_novak Exp $" ;
29 
30 /*
31  * $Id: map_af_radius.C,v 1.9 2014/10/13 08:53:03 j_novak Exp $
32  * $Log: map_af_radius.C,v $
33  * Revision 1.9 2014/10/13 08:53:03 j_novak
34  * Lorene classes and functions now belong to the namespace Lorene.
35  *
36  * Revision 1.8 2014/10/06 15:13:12 j_novak
37  * Modified #include directives to use c++ syntax.
38  *
39  * Revision 1.7 2013/06/05 15:10:42 j_novak
40  * Suppression of FINJAC sampling in r. This Jacobi(0,2) base is now
41  * available by setting colloc_r to BASE_JAC02 in the Mg3d constructor.
42  *
43  * Revision 1.6 2008/09/01 08:12:03 j_novak
44  * Improved test on the [rmin, rmax] interval.
45  *
46  * Revision 1.5 2007/12/11 15:28:14 jl_cornou
47  * Jacobi(0,2) polynomials partially implemented
48  *
49  * Revision 1.4 2006/09/13 13:59:21 j_novak
50  * Higher tolerance thereshold for Map_af::val_lx
51  *
52  * Revision 1.3 2006/07/10 07:44:51 j_novak
53  * Correction of the comparison between rmin and rr (now has to be greater than
54  * some threshold).
55  *
56  * Revision 1.2 2006/07/05 12:36:51 n_vasset
57  * Added a test on rmin to see whether the point lies in the computational domains.
58  *
59  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
60  * LORENE
61  *
62  * Revision 1.5 1999/12/16 14:19:08 eric
63  * Introduction de l'argument const Param& par dans val_lx et val_lx_jk.
64  * (en remplacement de l'argument Tbl& param).
65  *
66  * Revision 1.4 1999/12/07 14:51:37 eric
67  * val_r_kj --> val_r_jk
68  * val_lx_kj -->val_lx_jk
69  * Changement ordre des arguments val_r, val_lx
70  *
71  * Revision 1.3 1999/12/06 16:47:21 eric
72  * Surcharge de val_lx avec la version sans param.
73  *
74  * Revision 1.2 1999/12/06 15:34:06 eric
75  * Ajout des fonctions val_r_kj et val_lx_kj.
76  *
77  * Revision 1.1 1999/12/06 13:12:16 eric
78  * Initial revision
79  *
80  *
81  * $Header: /cvsroot/Lorene/C++/Source/Map/map_af_radius.C,v 1.9 2014/10/13 08:53:03 j_novak Exp $
82  *
83  */
84 
85 #include <cmath>
86 
87 // Headers Lorene
88 #include "map.h"
89 
90  //------------------------------//
91  // val_r //
92  //------------------------------//
93 
94 
95 namespace Lorene {
96 double Map_af::val_r(int l, double xi, double, double) const {
97 
98  assert( l>=0 ) ;
99  assert( l<mg->get_nzone() ) ;
100 
101  double resu ;
102 
103  switch( mg->get_type_r(l) ) {
104 
105  case FIN: case RARE: {
106  resu = alpha[l] * xi + beta[l] ;
107  break ;
108  }
109 
110  case UNSURR: {
111  resu = double(1) / ( alpha[l] * xi + beta[l] ) ;
112  break ;
113  }
114 
115  default: {
116  cout << "Map_af::val_r: unknown type_r ! " << endl ;
117  abort () ;
118  }
119  }
120 
121  return resu ;
122 }
123 
124  //------------------------------//
125  // val_lx //
126  //------------------------------//
127 
128 void Map_af::val_lx(double rr, double, double, int& lz, double& xi) const {
129 
130  // In which domain is located r ?
131  // ----------------------------
132  int nz = mg->get_nzone() ;
133  lz = - 1 ;
134 
135  for (int l=0; l<nz; l++) {
136 
137  double rmax = alpha[l] + beta[l] ;
138  double rmin = beta[l] - alpha[l] ;
139  if (mg->get_type_r(l) == RARE) rmin = 0. ;
140  if (mg->get_type_r(l) == UNSURR) {
141  rmin = double(1)/rmin ;
142  rmax = double(1)/rmax ;
143  }
144  if ((rr - rmin >= -1.e-14*fabs(rmin)) && ( rr <= rmax )) {
145  lz = l ;
146  break ;
147  }
148  } // fin de la boucle sur les zones
149 
150  if (lz == -1) { // On n'a pas trouve la zone
151  cout.precision(16);
152  cout.setf(ios::showpoint);
153  cout << "Map_af::val_lx: the domain containing r = " << rr <<
154  " has not been found ! "
155  << endl ;
156  for (int l=0; l<nz; l++) {
157  double rmin = -alpha[l] + beta[l] ;
158  if (mg->get_type_r(l) == UNSURR) rmin = double(1)/rmin ;
159  if (mg->get_type_r(l) == RARE) rmin = 0. ;
160  cout << "domain " << l << " : r_min = " << rmin ;
161  double rmax = alpha[l] + beta[l] ;
162  if (mg->get_type_r(l) == UNSURR) rmax = double(1)/rmax ;
163  cout << " : r_max = " << rmax << endl ;
164  }
165  abort () ;
166  }
167 
168  // Computation of xi
169  // -----------------
170 
171  switch( mg->get_type_r(lz) ) {
172 
173  case FIN: case RARE: {
174  xi = ( rr - beta[lz] ) / alpha[lz] ;
175  break ;
176  }
177 
178  case UNSURR: {
179  xi = ( double(1)/rr - beta[lz] ) / alpha[lz] ;
180  break ;
181  }
182 
183  default: {
184  cout << "Map_af::val_lx: unknown type_r ! " << endl ;
185  abort () ;
186  }
187  }
188 
189 }
190 
191 
192 void Map_af::val_lx(double rr, double, double, const Param&,
193  int& lz, double& xi) const {
194 
195  val_lx(rr, 0., 0., lz, xi) ;
196 
197 }
198 
199 
200  //------------------------------//
201  // val_r_jk //
202  //------------------------------//
203 
204 
205 double Map_af::val_r_jk(int l, double xi, int, int) const {
206 
207  return val_r(l, xi, 0., 0.) ;
208 
209 }
210 
211  //------------------------------//
212  // val_lx_jk //
213  //------------------------------//
214 
215 void Map_af::val_lx_jk(double rr, int, int, const Param& par,
216  int& l, double& xi) const {
217 
218  val_lx(rr, 0., 0., par, l, xi) ;
219 
220 }
221 
222 
223 }
Lorene::Map_af::val_r_jk
virtual double val_r_jk(int l, double xi, int j, int k) const
Returns the value of the radial coordinate r for a given and a given collocation point in in a give...
Definition: map_af_radius.C:205
Lorene
Lorene prototypes.
Definition: app_hor.h:64
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::Map_af::val_r
virtual double val_r(int l, double xi, double theta, double pphi) const
Returns the value of the radial coordinate r for a given in a given domain.
Definition: map_af_radius.C:96
Lorene::Map_af::beta
double * beta
Array (size: mg->nzone ) of the values of in each domain.
Definition: map.h:2035
Lorene::Map_af::val_lx_jk
virtual void val_lx_jk(double rr, int j, int k, const Param &par, int &l, double &xi) const
Computes the domain index l and the value of corresponding to a point of arbitrary r but collocation...
Definition: map_af_radius.C:215
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Map_af::val_lx
virtual void val_lx(double rr, double theta, double pphi, int &l, double &xi) const
Computes the domain index l and the value of corresponding to a point given by its physical coordina...
Definition: map_af_radius.C:128
Lorene::Param
Parameter storage.
Definition: param.h:125
Lorene::Map_af::alpha
double * alpha
Array (size: mg->nzone ) of the values of in each domain.
Definition: map.h:2033
Lorene::Map::mg
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined
Definition: map.h:676