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