LORENE
map_af_poisson_angu.C
1/*
2 * Resolution of the angular Poisson equation.
3 *
4 * (see file map.h for documentation)
5 *
6 */
7
8/*
9 * Copyright (c) 2003-2005 Eric Gourgoulhon & Jerome Novak
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 version 2
15 * as published by the Free Software Foundation.
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
28char map_af_poisson_angu_C[] = "$Header: /cvsroot/Lorene/C++/Source/Map/map_af_poisson_angu.C,v 1.4 2014/10/13 08:53:03 j_novak Exp $" ;
29
30/*
31 * $Id: map_af_poisson_angu.C,v 1.4 2014/10/13 08:53:03 j_novak Exp $
32 * $Log: map_af_poisson_angu.C,v $
33 * Revision 1.4 2014/10/13 08:53:03 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.3 2005/04/04 21:31:31 e_gourgoulhon
37 * Added argument lambda to method poisson_angu
38 * to deal with the generalized angular Poisson equation:
39 * Lap_ang u + lambda u = source.
40 *
41 * Revision 1.2 2003/10/16 08:49:23 j_novak
42 * Added a flag to decide wether the output is in the Ylm or in the standard base.
43 *
44 * Revision 1.1 2003/10/15 21:11:26 e_gourgoulhon
45 * Added method poisson_angu.
46 *
47 *
48 *
49 * $Header: /cvsroot/Lorene/C++/Source/Map/map_af_poisson_angu.C,v 1.4 2014/10/13 08:53:03 j_novak Exp $
50 *
51 */
52
53// Lorene headers
54#include "tensor.h"
55#include "param.h"
56
57namespace Lorene {
58void Map_af::poisson_angu(const Scalar& source, Param& par, Scalar& uu,
59 double lambda) const {
60
61 assert(source.get_etat() != ETATNONDEF) ;
62
63 assert(&(source.get_mp()) == this ) ;
64 assert(&(uu.get_mp()) == this ) ;
65
66 // Spherical harmonic expansion of the source
67 // ------------------------------------------
68
69 const Valeur& sourva = source.get_spectral_va() ;
70
71 if (sourva.get_etat() == ETATZERO) {
72 uu.set_etat_zero() ;
73 return ;
74 }
75
76 // Spectral coefficients of the source
77 assert(sourva.get_etat() == ETATQCQ) ;
78 sourva.coef() ;
79
80 Valeur resu(mg) ;
81 resu = *(sourva.c_cf) ; // copy of the coefficients of the source
82
83 resu.ylm() ; // spherical harmonic transform
84
85 // Call to the Mtbl_cf version
86 // ---------------------------
87 (resu.c_cf)->poisson_angu(lambda) ;
88
89 if (par.get_n_int() == 0) resu.ylm_i() ; // Back to standard bases
90 //in the case of no flag present
91 // in the Param
92
93 // Final result returned as a Scalar
94 // ---------------------------------
95
96 uu = resu ;
97
98 uu.set_dzpuis( source.get_dzpuis() ) ; // dzpuis unchanged
99}
100
101}
virtual void poisson_angu(const Scalar &source, Param &par, Scalar &uu, double lambda=0) const
Computes the solution of the generalized angular Poisson equation.
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined
Definition map.h:676
Parameter storage.
Definition param.h:125
int get_n_int() const
Returns the number of stored int 's addresses.
Definition param.C:239
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
int get_dzpuis() const
Returns dzpuis.
Definition scalar.h:557
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition scalar.h:601
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition scalar.h:554
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition scalar.C:808
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.
void ylm_i()
Inverse of ylm()
const Map & get_mp() const
Returns the mapping.
Definition tensor.h:861
Lorene prototypes.
Definition app_hor.h:64