LORENE
ope_pois_vect_r.C
1 /*
2  * Methods of class Ope_pois_vect_r
3  *
4  * (see file ope_elementary.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 2004 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 
28 char ope_pois_vect_r_C[] = "$Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_pois_vect_r/ope_pois_vect_r.C,v 1.2 2014/10/13 08:53:34 j_novak Exp $" ;
29 
30 /*
31  * $Id: ope_pois_vect_r.C,v 1.2 2014/10/13 08:53:34 j_novak Exp $
32  * $Log: ope_pois_vect_r.C,v $
33  * Revision 1.2 2014/10/13 08:53:34 j_novak
34  * Lorene classes and functions now belong to the namespace Lorene.
35  *
36  * Revision 1.1 2004/05/10 15:28:22 j_novak
37  * First version of functions for the solution of the r-component of the
38  * vector Poisson equation.
39  *
40  *
41  * $Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_pois_vect_r/ope_pois_vect_r.C,v 1.2 2014/10/13 08:53:34 j_novak Exp $
42  *
43  */
44 
45 #include"type_parite.h"
46 #include "ope_elementary.h"
47 
48 namespace Lorene {
49 Matrice ope_pvect_r_mat(int , int , double , int, int ) ;
50 Tbl sh_pvect_r(int, int, double, int) ;
51 Matrice cl_pvect_r (const Matrice&, int, double, int, int) ;
52 Matrice nondeg_pvect_r (const Matrice&, int, double, int, int) ;
53 Tbl val_solh (int, double, double, int) ;
54 
55 // Standard constructor :
56 Ope_pois_vect_r::Ope_pois_vect_r (int nbr, int baser, double alf, double bet, int lq, int dz):
57  Ope_poisson(nbr, baser, alf, bet, lq, dz)
58 {
59  assert (dzpuis == 4) ;
60 }
61 
62 // Constructor by copy :
64 {
65  assert (dzpuis == 4) ;
66 
67 }
68 
69 // Destructor :
71 
72 // True functions :
74  if (ope_mat != 0x0)
75  delete ope_mat ;
76 
77  ope_mat = new Matrice
78  (ope_pvect_r_mat(nr, l_quant, beta/alpha, dzpuis, base_r)) ;
79 }
80 
82  if (ope_mat == 0x0)
83  do_ope_mat() ;
84 
85  if (ope_cl != 0x0)
86  delete ope_cl ;
87 
88  ope_cl = new Matrice
89  (cl_pvect_r(*ope_mat, l_quant, beta/alpha, dzpuis, base_r)) ;
90 }
91 
93  if (ope_cl == 0x0)
94  do_ope_cl() ;
95 
96  if (non_dege != 0x0)
97  delete non_dege ;
98 
99  non_dege = new Matrice
100  (nondeg_pvect_r(*ope_cl, l_quant, beta/alpha, dzpuis, base_r)) ;
101 }
102 
104 
105  int l1 = ( (l_quant == 0) ? 1 : l_quant - 1 ) ;
106  int l2 = l_quant + 1 ;
107 
108  Tbl valeurs1 (val_solh (l1, alpha, beta, base_r)) ;
109  Tbl valeurs2 (val_solh (l2, alpha, beta, base_r)) ;
110 
111  assert (valeurs1.get_ndim() == valeurs2.get_ndim()) ;
112 
113  if (valeurs1.get_ndim() == 2) {
114  // cas 2 sh
115  s_one_plus = valeurs1(0,0) ;
116  s_one_minus = valeurs1(0,1) ;
117  ds_one_plus = valeurs1(0,2) ;
118  ds_one_minus = valeurs1(0,3) ;
119 
120  s_two_plus = valeurs2(1,0) ;
121  s_two_minus = valeurs2(1,1) ;
122  ds_two_plus = valeurs2(1,2) ;
123  ds_two_minus = valeurs2(1,3) ;
124  }
125  else {
126  // cas 1 sh :
127  Tbl& valeurs = (base_r == R_CHEBU) ? valeurs2 : valeurs1 ;
128  s_one_plus = valeurs(0) ;
129  s_one_minus = valeurs(1) ;
130  ds_one_plus = valeurs(2) ;
131  ds_one_minus = valeurs(3) ;
132  }
133 
134  return sh_pvect_r(nr, l_quant, beta/alpha, base_r) ;
135 }
136 
137 }
Matrix handling.
Definition: matrice.h:152
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 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 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.
int nr
Number of radial points.
Class for the operator of the r component of the vector Poisson equation.
virtual Tbl get_solh() const
Computes the homogeneous solutions(s).
virtual void do_ope_cl() const
Computes the banded-matrix of the operator.
Ope_pois_vect_r(int nbr, int baser, double alf, double bet, int lq, int dz)
Standard constructor.
virtual void do_non_dege() const
Computes the non-degenerated matrix of the operator.
virtual void do_ope_mat() const
Computes the matrix of the operator.
virtual ~Ope_pois_vect_r()
Destructor.
Class for the operator of the Poisson equation (i.e.
int dzpuis
the associated dzpuis, if in the compactified domain.
int l_quant
quantum number
Basic array class.
Definition: tbl.h:161
int get_ndim() const
Gives the number of dimensions (ie dim.ndim)
Definition: tbl.h:400
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
Definition: type_parite.h:180
Lorene prototypes.
Definition: app_hor.h:64