LORENE
scalar_asymptot.C
1 /*
2  * Function Scalar::asymptot
3  *
4  */
5 
6 /*
7  * Copyright (c) 2003 Eric Gourgoulhon & Jerome Novak
8  *
9  * Copyright (c) 1999-2002 Eric Gourgoulhon (for previous class Cmp)
10  * Copyright (c) 1999-2001 Philippe Grandclement (for previous class Cmp)
11  *
12  * This file is part of LORENE.
13  *
14  * LORENE is free software; you can redistribute it and/or modify
15  * it under the terms of the GNU General Public License as published by
16  * the Free Software Foundation; either version 2 of the License, or
17  * (at your option) any later version.
18  *
19  * LORENE is distributed in the hope that it will be useful,
20  * but WITHOUT ANY WARRANTY; without even the implied warranty of
21  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22  * GNU General Public License for more details.
23  *
24  * You should have received a copy of the GNU General Public License
25  * along with LORENE; if not, write to the Free Software
26  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
27  *
28  */
29 
30 
31 char scalar_asymptot_C[] = "$Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_asymptot.C,v 1.5 2014/10/13 08:53:46 j_novak Exp $" ;
32 
33 /*
34  * $Id: scalar_asymptot.C,v 1.5 2014/10/13 08:53:46 j_novak Exp $
35  * $Log: scalar_asymptot.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 2004/02/21 17:05:14 e_gourgoulhon
43  * Method Scalar::point renamed Scalar::val_grid_point.
44  * Method Scalar::set_point renamed Scalar::set_grid_point.
45  *
46  * Revision 1.2 2003/10/08 14:24:09 j_novak
47  * replaced mult_r_zec with mult_r_ced
48  *
49  * Revision 1.1 2003/09/25 07:18:00 j_novak
50  * Method asymptot implemented.
51  *
52  * $Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_asymptot.C,v 1.5 2014/10/13 08:53:46 j_novak Exp $
53  *
54  */
55 
56 // Headers C
57 #include <cmath>
58 
59 // Headers Lorene
60 #include "tensor.h"
61 
62 namespace Lorene {
63 Valeur** Scalar::asymptot(int n0, const int flag) const {
64 
65  assert(n0 >= 0) ;
66  const Mg3d& mg = *(mp->get_mg()) ;
67  int nz = mg.get_nzone() ;
68  int nzm1 = nz-1 ;
69  assert(mg.get_type_r(nzm1) == UNSURR) ;
70  int np = mg.get_np(nzm1) ;
71  int nt = mg.get_nt(nzm1) ;
72  int nr = mg.get_nr(nzm1) ;
73  int nrm1 = nr-1 ;
74 
75  Valeur** vu = new Valeur*[n0+1] ;
76  for (int h=0; h<=n0; h++) {
77  vu[h] = new Valeur(mg.get_angu()) ;
78  }
79 
80  Scalar uu = *this ;
81 
82  int precis = 2 ;
83 
84  // The terms 1/r^h with h < dzpuis are null :
85  // -----------------------------------------
86  for (int h=0; h<dzpuis; h++) {
87 
88  vu[h]->set_etat_zero() ;
89 
90  }
91 
92  // Terms 1/r^h with h >= dzpuis :
93  // -----------------------------
94  for (int h=dzpuis; h<=n0; h++) {
95 
96  // Value on the sphere S^2 at infinity
97  // -----------------------------------
98  vu[h]->set_etat_c_qcq() ;
99  vu[h]->c->set_etat_qcq() ;
100  for (int l=0; l<nzm1; l++) {
101  vu[h]->c->t[l]->set_etat_zero() ;
102  }
103  vu[h]->c->t[nzm1]->set_etat_qcq() ;
104 
105  for (int k=0; k<np; k++) {
106  for (int j=0; j<nt; j++) {
107  vu[h]->set(nzm1, k, j, 0) = uu.val_grid_point(nzm1, k, j, nrm1) ;
108  }
109  }
110 
111  vu[h]->set_base( uu.va.base ) ;
112 
113  // Printing
114  // --------
115  if (flag != 0) {
116  cout << "Term in 1/r^" << h << endl ;
117  cout << "-------------" << endl ;
118 
119  double ndec = 0 ;
120  double vmin = (*vu[h])(nzm1, 0, 0, 0) ;
121  double vmax = vmin ;
122 
123  cout << " Values at the point (phi_k, theta_j) : " << endl ;
124  cout.precision(precis) ;
125  cout.setf(ios::showpoint) ;
126  for (int k=0; k<np; k++) {
127  cout << " k=" << k << " : " ;
128  for (int j=0; j<nt; j++) {
129  double xx = (*vu[h])(nzm1, k, j, 0) ;
130  cout << " " << setw(precis) << xx ;
131  ndec += fabs(xx) ;
132  vmin = ( xx < vmin ) ? xx : vmin ;
133  vmax = ( xx > vmax ) ? xx : vmax ;
134  }
135  cout << endl;
136  }
137  ndec /= np*nt ;
138  cout << "Minimum value on S^2 : " << vmin << endl ;
139  cout << "Maximum value on S^2 : " << vmax << endl ;
140  cout << "L^1 norm on S^2 : " << ndec << endl ;
141  }
142  // The value at infinity is substracted
143  // ------------------------------------
144  for (int k=0; k<np; k++) {
145  for (int j=0; j<nt; j++) {
146  double v_inf = (*vu[h])(nzm1, k, j, 0) ;
147  for (int i=0; i<nr; i++) {
148  uu.set_grid_point(nzm1, k, j, i) -= v_inf ;
149  }
150  }
151  }
152 
153  // Mutliplication by r
154  // -------------------
155 
156  uu.mult_r_ced() ;
157 
158  } // End of loop on h (development order)
159 
160  return vu ;
161 
162 }
163 }
Lorene::Scalar::mult_r_ced
void mult_r_ced()
Multiplication by r in the compactified external domain (CED), the dzpuis flag is not changed.
Definition: scalar_r_manip.C:269
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::Valeur
Values and coefficients of a (real-value) function.
Definition: valeur.h:287
Lorene::Valeur::set_etat_c_qcq
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
Lorene::Scalar::va
Valeur va
The numerical value of the Scalar
Definition: scalar.h:405
Lorene::Mg3d
Multi-domain grid.
Definition: grilles.h:273
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::Mg3d::get_angu
const Mg3d * get_angu() const
Returns the pointer on the associated angular grid.
Definition: mg3d.C:473
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Scalar::set_grid_point
double & set_grid_point(int l, int k, int j, int i)
Setting the value of the field at a given grid point.
Definition: scalar.h:684
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
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Lorene::Valeur::c
Mtbl * c
Values of the function at the points of the multi-grid
Definition: valeur.h:299
Lorene::Tbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:361
Lorene::Scalar::asymptot
Valeur ** asymptot(int n, const int flag=0) const
Asymptotic expansion at r = infinity.
Definition: scalar_asymptot.C:63
Lorene::Tensor::mp
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition: tensor.h:295
Lorene::Tbl::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: tbl.C:347
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
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::Valeur::set_base
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:810
Lorene::Scalar::dzpuis
int dzpuis
Power of r by which the quantity represented by this must be divided in the compactified external d...
Definition: scalar.h:403
Lorene::Mtbl::t
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's.
Definition: mtbl.h:132
Lorene::Valeur::base
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:305
Lorene::Valeur::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: valeur.C:689
Lorene::Valeur::set
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition: valeur.h:363
Lorene::Scalar::val_grid_point
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
Definition: scalar.h:637
Lorene::Mtbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:299