LORENE
param_elliptic_pseudo_1d.C
1 /*
2  * Copyright (c) 2004 Philippe Grandclement
3  *
4  * This file is part of LORENE.
5  *
6  * LORENE is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License version 2
8  * as published by the Free Software Foundation.
9  *
10  * LORENE is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with LORENE; if not, write to the Free Software
17  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18  *
19  */
20 
21 
22 
23 /*
24  * $Id: param_elliptic_pseudo_1d.C,v 1.5 2018/11/16 14:34:37 j_novak Exp $
25  * $Log: param_elliptic_pseudo_1d.C,v $
26  * Revision 1.5 2018/11/16 14:34:37 j_novak
27  * Changed minor points to avoid some compilation warnings.
28  *
29  * Revision 1.4 2016/12/05 16:18:14 j_novak
30  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
31  *
32  * Revision 1.3 2014/10/13 08:53:37 j_novak
33  * Lorene classes and functions now belong to the namespace Lorene.
34  *
35  * Revision 1.2 2014/10/06 15:13:15 j_novak
36  * Modified #include directives to use c++ syntax.
37  *
38  * Revision 1.1 2004/08/24 09:14:49 p_grandclement
39  * Addition of some new operators, like Poisson in 2d... It now requieres the
40  * GSL library to work.
41  *
42  * Also, the way a variable change is stored by a Param_elliptic is changed and
43  * no longer uses Change_var but rather 2 Scalars. The codes using that feature
44  * will requiere some modification. (It should concern only the ones about monopoles)
45  *
46  *
47  * $Header: /cvsroot/Lorene/C++/Source/Param_elliptic/param_elliptic_pseudo_1d.C,v 1.5 2018/11/16 14:34:37 j_novak Exp $
48  *
49  */
50 
51 #include "headcpp.h"
52 
53 #include <cmath>
54 #include <cstdlib>
55 
56 #include "param_elliptic.h"
57 #include "base_val.h"
58 #include "map.h"
59 #include "ope_elementary.h"
60 #include "change_var.h"
61 #include "scalar.h"
62 
63 
64 namespace Lorene {
66 
67  if (type_map != MAP_AFF) {
68  cout << "set_poisson_pseudo_1d only defined for an affine mapping..." << endl ;
69  abort() ;
70  }
71  else {
72 
73  int nz = get_mp().get_mg()->get_nzone() ;
74 
75  int nr ;
76  double alpha, beta ;
77  int m_quant, l_quant, base_r_1d ;
78 
79  int conte = 0 ;
80  for (int l=0 ; l<nz ; l++) {
81 
82  nr = get_mp().get_mg()->get_nr(l) ;
83  alpha = get_alpha (l) ;
84  beta = get_beta (l) ;
85 
86  for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
87  for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
88  if (operateurs[conte] != 0x0)
89  delete operateurs[conte] ;
90  source.get_spectral_va().base.give_quant_numbers(l, k, j, m_quant, l_quant, base_r_1d) ;
91  if ((k!=1) && (l!=nz-1))
92  operateurs[conte] = new Ope_poisson_pseudo_1d (nr, base_r_1d, alpha, beta, l_quant) ;
93  else
94  operateurs[conte] = 0x0 ;
95  conte ++ ;
96  }
97  }
98  }
99 }
100 
101 void Param_elliptic::set_helmholtz_minus_pseudo_1d(int zone, double masse, Scalar& source) {
102 
103  int dzpuis = source.get_dzpuis() ;
104  assert (masse > 0) ;
105 
106  if (type_map != MAP_AFF) {
107  cout << "set_helmholtz_minus_pseudo_1d only defined for an affine mapping..." << endl ;
108  abort() ;
109  }
110  else {
111 
112  int nz = get_mp().get_mg()->get_nzone() ;
113  if (zone == nz-1)
114  source.check_dzpuis(2) ;
115  int nr ;
116  double alpha, beta ;
117  int m_quant, l_quant, base_r_1d ;
118 
119  int conte = 0 ;
120  for (int l=0 ; l<nz ; l++) {
121 
122  nr = get_mp().get_mg()->get_nr(l) ;
123  alpha = get_alpha (l) ;
124  beta = get_beta (l) ;
125 
126  for (int k=0 ; k<get_mp().get_mg()->get_np(l)+1 ; k++)
127  for (int j=0 ; j<get_mp().get_mg()->get_nt(l) ; j++) {
128  if (l==zone) {
129  if (operateurs[conte] != 0x0)
130  delete operateurs[conte] ;
132  (l, k, j, m_quant, l_quant, base_r_1d) ;
133  operateurs[conte] = new Ope_helmholtz_minus_pseudo_1d (nr, base_r_1d,
134  alpha, beta, l_quant, masse, dzpuis) ;
135  }
136  conte ++ ;
137  }
138  }
139  }
140 }
141 
142 }
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
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
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
Class for the operator of the modified Helmholtz equation in pseudo-1d.
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:315
int get_dzpuis() const
Returns dzpuis.
Definition: scalar.h:563
int type_map
Type of mapping either MAP_AFF or MAP_LOG.
const Map_radial & get_mp() const
Returns the mapping.
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Class for the operator of the Poisson equation in pseudo 1d.
Ope_elementary ** operateurs
Array on the elementary operators.
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
void set_poisson_pseudo_1d(Scalar &so)
Set the operator to everywhere but in the compactified domain.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
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: scalar.C:879
void set_helmholtz_minus_pseudo_1d(int zone, double mas, Scalar &so)
Set the operator to in one domain.
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:607