LORENE
ope_poisson.C
1/*
2 * Copyright (c) 2003 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
21char ope_poisson_C[] = "$Header: /cvsroot/Lorene/C++/Source/Ope_elementary/ope_poisson.C,v 1.3 2014/10/13 08:53:33 j_novak Exp $" ;
22
23/*
24 * $Id: ope_poisson.C,v 1.3 2014/10/13 08:53:33 j_novak Exp $
25 * $Log: ope_poisson.C,v $
26 * Revision 1.3 2014/10/13 08:53:33 j_novak
27 * Lorene classes and functions now belong to the namespace Lorene.
28 *
29 * Revision 1.2 2004/06/14 15:07:11 j_novak
30 * New methods for the construction of the elliptic operator appearing in
31 * the vector Poisson equation (acting on eta).
32 *
33 * Revision 1.1 2003/12/11 14:48:50 p_grandclement
34 * Addition of ALL (and that is a lot !) the files needed for the general elliptic solver ... UNDER DEVELOPEMENT...
35 *
36 *
37 * $Header: /cvsroot/Lorene/C++/Source/Ope_elementary/ope_poisson.C,v 1.3 2014/10/13 08:53:33 j_novak Exp $
38 *
39 */
40
41#include "proto.h"
42#include "ope_elementary.h"
43
44// Standard constructor :
45namespace Lorene {
46Ope_poisson::Ope_poisson (int nbr, int baser, double alf, double bet, int lq, int dz):
47 Ope_elementary(nbr, baser, alf, bet), l_quant (lq),
48 dzpuis (dz) {
49
50 assert ((dzpuis==2) || (dzpuis==3) || (dzpuis==4)) ;
51}
52
53// Constructor by copy :
55 l_quant (so.l_quant), dzpuis (so.dzpuis) {
56
57 assert ((dzpuis==2) || (dzpuis==3) || (dzpuis==4)) ;
58}
59
60// Destructor :
62
63// True functions :
65 if (ope_mat != 0x0)
66 delete ope_mat ;
67
68 ope_mat = new Matrice
69 (laplacien_mat(nr, l_quant, beta/alpha, dzpuis, base_r)) ;
70}
71
73 if (ope_mat == 0x0)
74 do_ope_mat() ;
75
76 if (ope_cl != 0x0)
77 delete ope_cl ;
78
79 ope_cl = new Matrice
80 (combinaison(*ope_mat, l_quant, beta/alpha, dzpuis, base_r)) ;
81}
82
84 if (ope_cl == 0x0)
85 do_ope_cl() ;
86
87 if (non_dege != 0x0)
88 delete non_dege ;
89
90 non_dege = new Matrice
91 (prepa_nondege(*ope_cl, l_quant, beta/alpha, dzpuis, base_r)) ;
92}
93
95
96 if (non_dege == 0x0)
97 do_non_dege() ;
98
100
101 Tbl valeurs (val_solp (res, alpha, base_r)) ;
102 sp_plus = valeurs(0) ;
103 sp_minus = valeurs(1) ;
104 dsp_plus = valeurs(2) ;
105 dsp_minus = valeurs(3) ;
106
107 return res ;
108}
109
111
112 Tbl valeurs (val_solh (l_quant, alpha, beta, base_r)) ;
113 if (valeurs.get_ndim() == 2) {
114 // cas 2 sh
115 s_one_plus = valeurs(0,0) ;
116 s_one_minus = valeurs(0,1) ;
117 ds_one_plus = valeurs(0,2) ;
118 ds_one_minus = valeurs(0,3) ;
119
120 s_two_plus = valeurs(1,0) ;
121 s_two_minus = valeurs(1,1) ;
122 ds_two_plus = valeurs(1,2) ;
123 ds_two_minus = valeurs(1,3) ;
124 }
125 else {
126 // cas 1 sh :
127 s_one_plus = valeurs(0) ;
128 s_one_minus = valeurs(1) ;
129 ds_one_plus = valeurs(2) ;
130 ds_one_minus = valeurs(3) ;
131 }
132
133 return solh(nr, l_quant, beta/alpha, base_r) ;
134}
135
137
138 if (ope_mat != 0x0) {
139 delete ope_mat ;
140 ope_mat = 0x0 ;
141 }
142
143 if (ope_cl != 0x0) {
144 delete ope_cl ;
145 ope_cl = 0x0 ;
146 }
147
148 if (non_dege != 0x0) {
149 delete non_dege ;
150 non_dege = 0x0 ;
151 }
152 l_quant ++ ;
153}
154
156
157 assert(l_quant > 0) ;
158
159 if (ope_mat != 0x0) {
160 delete ope_mat ;
161 ope_mat = 0x0 ;
162 }
163
164 if (ope_cl != 0x0) {
165 delete ope_cl ;
166 ope_cl = 0x0 ;
167 }
168
169 if (non_dege != 0x0) {
170 delete non_dege ;
171 non_dege = 0x0 ;
172 }
173 l_quant -- ;
174}
175}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Matrix handling.
Definition matrice.h:152
Basic class for elementary elliptic operators.
Matrice * ope_mat
Pointer on the matrix representation of the operator.
double ds_two_minus
Value of the derivative of the second homogeneous solution at the inner boundary.
double s_two_plus
Value of the second homogeneous solution at the outer boundary.
double s_one_minus
Value of the first homogeneous solution at the inner boundary.
double beta
Parameter of the associated mapping.
double dsp_plus
Value of the derivative of the particular solution at the outer boundary.
double ds_one_plus
Value of the derivative of the first homogeneous solution at the outer boundary.
double ds_one_minus
Value of the derivative of the first homogeneous solution at the inner boundary.
double alpha
Parameter of the associated mapping.
double s_two_minus
Value of the second homogeneous solution at the inner boundary.
int base_r
Radial basis of decomposition.
double s_one_plus
Value of the first homogeneous solution at the outer boundary.
Matrice * ope_cl
Pointer on the banded-matrix of the operator.
double sp_minus
Value of the particular solution at the inner boundary.
double ds_two_plus
Value of the derivative of the second homogeneous solution at the outer boundary.
Matrice * non_dege
Pointer on the non-degenerated matrix of the operator.
double sp_plus
Value of the particular solution at the outer boundary.
int nr
Number of radial points.
double dsp_minus
Value of the derivative of the particular solution at the inner boundary.
Class for the operator of the Poisson equation (i.e.
virtual Tbl get_solh() const
Computes the homogeneous solutions(s).
virtual ~Ope_poisson()
Destructor.
Definition ope_poisson.C:61
int dzpuis
the associated dzpuis, if in the compactified domain.
virtual void dec_l_quant()
Decreases the quatum number l by one unit.
virtual void do_non_dege() const
Computes the non-degenerated matrix of the operator.
Definition ope_poisson.C:83
virtual void do_ope_cl() const
Computes the banded-matrix of the operator.
Definition ope_poisson.C:72
virtual void do_ope_mat() const
Computes the matrix of the operator.
Definition ope_poisson.C:64
virtual void inc_l_quant()
Increases the quatum number l by one unit.
int l_quant
quantum number
virtual Tbl get_solp(const Tbl &so) const
Computes the particular solution, given the source so .
Definition ope_poisson.C:94
Ope_poisson(int nbr, int baser, double alf, double bet, int lq, int dz)
Standard constructor.
Definition ope_poisson.C:46
Basic array class.
Definition tbl.h:161
Lorene prototypes.
Definition app_hor.h:64