LORENE
map_af_poisson.C
1 /*
2  * Method of the class Map_af for the resolution of the scalar Poisson
3  * equation
4  */
5 
6 /*
7  * Copyright (c) 1999-2001 Eric Gourgoulhon
8  * Copyright (c) 1999-2001 Philippe Grandclement
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: map_af_poisson.C,v 1.7 2016/12/05 16:17:57 j_novak Exp $
33  * $Log: map_af_poisson.C,v $
34  * Revision 1.7 2016/12/05 16:17:57 j_novak
35  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
36  *
37  * Revision 1.6 2014/10/13 08:53:02 j_novak
38  * Lorene classes and functions now belong to the namespace Lorene.
39  *
40  * Revision 1.5 2005/08/25 12:14:09 p_grandclement
41  * Addition of a new method to solve the scalar Poisson equation, based on a multi-domain Tau-method
42  *
43  * Revision 1.4 2004/05/06 15:25:39 e_gourgoulhon
44  * The case dzpuis=5 with null value in CED is well treated now.
45  *
46  * Revision 1.3 2004/02/20 10:55:23 j_novak
47  * The versions dzpuis 5 -> 3 has been improved and polished. Should be
48  * operational now...
49  *
50  * Revision 1.2 2004/02/06 10:53:52 j_novak
51  * New dzpuis = 5 -> dzpuis = 3 case (not ready yet).
52  *
53  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
54  * LORENE
55  *
56  * Revision 1.9 2000/05/22 13:46:48 phil
57  * ajout du cas dzpuis = 3
58  *
59  * Revision 1.8 2000/02/09 14:44:24 eric
60  * Traitement de dzpuis ameliore.
61  *
62  * Revision 1.7 1999/12/22 16:37:10 eric
63  * Ajout de pot.set_dzpuis(0) a la fin.
64  *
65  * Revision 1.6 1999/12/22 15:11:03 eric
66  * Remplacement du test source.get_mp() == this par
67  * source.get_mp()->get_mg() == mg
68  * (idem pour pot),
69  * afin de permettre l'appel par Map_et::poisson.
70  *
71  * Revision 1.5 1999/12/21 13:02:37 eric
72  * Changement de prototype de la routine poisson : la solution est
73  * desormais passee en argument (et non plus en valeur de retour)
74  * pour s'adapter au prototype general de la fonction virtuelle
75  * Map::poisson.
76  *
77  * Revision 1.4 1999/12/21 10:06:29 eric
78  * Ajout de l'argument (muet) Param&.
79  *
80  * Revision 1.3 1999/12/07 16:48:50 phil
81  * On fait ylm_i avant de quitter
82  *
83  * Revision 1.2 1999/12/02 16:12:22 eric
84  * *** empty log message ***
85  *
86  * Revision 1.1 1999/12/02 14:30:07 eric
87  * Initial revision
88  *
89  *
90  * $Header: /cvsroot/Lorene/C++/Source/Map/map_af_poisson.C,v 1.7 2016/12/05 16:17:57 j_novak Exp $
91  *
92  */
93 
94 // Header Lorene:
95 #include "map.h"
96 #include "cmp.h"
97 
98 namespace Lorene {
99 Mtbl_cf sol_poisson(const Map_af&, const Mtbl_cf&, int, bool match = true) ;
100 Mtbl_cf sol_poisson_tau(const Map_af&, const Mtbl_cf&, int) ;
101 //*****************************************************************************
102 
103 void Map_af::poisson(const Cmp& source, Param& , Cmp& pot) const {
104 
105  assert(source.get_etat() != ETATNONDEF) ;
106  assert(source.get_mp()->get_mg() == mg) ;
107  assert(pot.get_mp()->get_mg() == mg) ;
108 
109  assert( source.check_dzpuis(2) || source.check_dzpuis(4)
110  || source.check_dzpuis(3) || source.check_dzpuis(5) ) ;
111 
112  bool match = true ;
113 
114  int dzpuis ;
115 
116  if ( (source.dz_nonzero()) || (source.get_dzpuis() > 3)) { //##awkward??
117  dzpuis = source.get_dzpuis() ;
118  }
119  else{
120  dzpuis = 4 ;
121  }
122 
123  match = !(dzpuis == 5) ;
124 
125  // Spherical harmonic expansion of the source
126  // ------------------------------------------
127 
128  const Valeur& sourva = source.va ;
129 
130  if (sourva.get_etat() == ETATZERO) {
131  pot.set_etat_zero() ;
132  return ;
133  }
134 
135  // Spectral coefficients of the source
136  assert(sourva.get_etat() == ETATQCQ) ;
137 
138  Valeur rho(sourva.get_mg()) ;
139  sourva.coef() ;
140  rho = *(sourva.c_cf) ; // copy of the coefficients of the source
141 
142  rho.ylm() ; // spherical harmonic transforms
143 
144  // Call to the Mtbl_cf version
145  // ---------------------------
146  Mtbl_cf resu = sol_poisson(*this, *(rho.c_cf), dzpuis, match) ;
147 
148  // Final result returned as a Cmp
149  // ------------------------------
150 
151  pot.set_etat_zero() ; // to call Cmp::del_t().
152 
153  pot.set_etat_qcq() ;
154 
155  pot.va = resu ;
156  (pot.va).ylm_i() ; // On repasse en base standard.
157 
158  (dzpuis == 5) ? pot.set_dzpuis(3) : pot.set_dzpuis(0) ;
159 
160 }
161 
162 
163  //----------------------
164  // Tau version method
165  //---------------------
166 
167 
168 void Map_af::poisson_tau(const Cmp& source, Param& , Cmp& pot) const {
169 
170  assert(source.get_etat() != ETATNONDEF) ;
171  assert(source.get_mp()->get_mg() == mg) ;
172  assert(pot.get_mp()->get_mg() == mg) ;
173 
174  assert( source.check_dzpuis(2) || source.check_dzpuis(4)
175  || source.check_dzpuis(3)) ;
176 
177  int dzpuis ;
178 
179  if ( (source.dz_nonzero()) || (source.get_dzpuis() > 3)) { //##awkward??
180  dzpuis = source.get_dzpuis() ;
181  }
182  else{
183  dzpuis = 4 ;
184  }
185 
186  // Spherical harmonic expansion of the source
187  // ------------------------------------------
188 
189  const Valeur& sourva = source.va ;
190 
191  if (sourva.get_etat() == ETATZERO) {
192  pot.set_etat_zero() ;
193  return ;
194  }
195 
196  // Spectral coefficients of the source
197  assert(sourva.get_etat() == ETATQCQ) ;
198 
199  Valeur rho(sourva.get_mg()) ;
200  sourva.coef() ;
201  rho = *(sourva.c_cf) ; // copy of the coefficients of the source
202 
203  rho.ylm() ; // spherical harmonic transforms
204 
205  // Call to the Mtbl_cf version
206  // ---------------------------
207 
208  Mtbl_cf resu = sol_poisson_tau(*this, *(rho.c_cf), dzpuis) ;
209 
210  // Final result returned as a Cmp
211  // ------------------------------
212 
213  pot.set_etat_zero() ; // to call Cmp::del_t().
214 
215  pot.set_etat_qcq() ;
216 
217  pot.va = resu ;
218  (pot.va).ylm_i() ; // On repasse en base standard.
219  pot.set_dzpuis(0) ;
220 }
221 
222 }
const Map * get_mp() const
Returns the mapping.
Definition: cmp.h:901
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:312
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:151
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
int get_etat() const
Returns the logical state.
Definition: cmp.h:899
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
const Mg3d * get_mg() const
Returns the Mg3d on which the this is defined.
Definition: valeur.h:763
int get_etat() const
Returns the logical state.
Definition: valeur.h:760
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: cmp.C:292
virtual void poisson(const Cmp &source, Param &par, Cmp &uu) const
Computes the solution of a scalar Poisson equation.
bool dz_nonzero() const
Returns true if the last domain is compactified and *this is not zero in this domain.
Definition: cmp.C:663
Parameter storage.
Definition: param.h:125
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: cmp.C:307
virtual void poisson_tau(const Cmp &source, Param &par, Cmp &uu) const
Computes the solution of a scalar Poisson equation using a Tau method.
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined.
Definition: map.h:694
int get_dzpuis() const
Returns dzpuis.
Definition: cmp.h:903
Coefficients storage for the multi-domain spectral method.
Definition: mtbl_cf.h:196
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: cmp.C:718
void set_dzpuis(int)
Set a value to dzpuis.
Definition: cmp.C:657
Valeur va
The numerical value of the Cmp.
Definition: cmp.h:464