LORENE
ope_pois_tens_rr.C
1 /*
2  * Methods of class Ope_pois_tens_rr
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_tens_rr_C[] = "$Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_pois_vect_r/ope_pois_tens_rr.C,v 1.2 2014/10/13 08:53:34 j_novak Exp $" ;
29 
30 /*
31  * $Id: ope_pois_tens_rr.C,v 1.2 2014/10/13 08:53:34 j_novak Exp $
32  * $Log: ope_pois_tens_rr.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/12/23 16:30:15 j_novak
37  * New files and class for the solution of the rr component of the tensor Poisson
38  * equation.
39  *
40  *
41  *
42  * $Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_pois_vect_r/ope_pois_tens_rr.C,v 1.2 2014/10/13 08:53:34 j_novak Exp $
43  *
44  */
45 
46 #include"type_parite.h"
47 #include "ope_elementary.h"
48 
49 namespace Lorene {
50 Matrice ope_ptens_rr_mat(int , int , double , int, int ) ;
51 Tbl sh_ptens_rr(int, int, double, int) ;
52 Matrice cl_ptens_rr (const Matrice&, int, double, int, int) ;
53 Matrice nondeg_ptens_rr (const Matrice&, int, double, int, int) ;
54 Tbl val_solh (int, double, double, int) ;
55 
56 // Standard constructor :
57 Ope_pois_tens_rr::Ope_pois_tens_rr (int nbr, int baser, double alf, double bet, int lq, int dz):
58  Ope_poisson(nbr, baser, alf, bet, lq, dz)
59 {
60  assert (dzpuis == 4) ;
61  assert (l_quant > 1) ;
62 }
63 
64 // Constructor by copy :
66 {
67  assert (dzpuis == 4) ;
68  assert (l_quant > 1) ;
69 }
70 
71 // Destructor :
73 
74 // True functions :
76  if (ope_mat != 0x0)
77  delete ope_mat ;
78 
79  ope_mat = new Matrice
80  (ope_ptens_rr_mat(nr, l_quant, beta/alpha, dzpuis, base_r)) ;
81 }
82 
84  if (ope_mat == 0x0)
85  do_ope_mat() ;
86 
87  if (ope_cl != 0x0)
88  delete ope_cl ;
89 
90  ope_cl = new Matrice
91  (cl_ptens_rr(*ope_mat, l_quant, beta/alpha, dzpuis, base_r)) ;
92 }
93 
95  if (ope_cl == 0x0)
96  do_ope_cl() ;
97 
98  if (non_dege != 0x0)
99  delete non_dege ;
100 
101  non_dege = new Matrice
102  (nondeg_ptens_rr(*ope_cl, l_quant, beta/alpha, dzpuis, base_r)) ;
103 }
104 
106 
107  assert (l_quant > 1) ;
108  int l1 = l_quant - 2 ;
109  int l2 = l_quant + 2 ;
110 
111  Tbl valeurs1 (val_solh (l1, alpha, beta, base_r)) ;
112  Tbl valeurs2 (val_solh (l2, alpha, beta, base_r)) ;
113 
114  assert (valeurs1.get_ndim() == valeurs2.get_ndim()) ;
115 
116  if (valeurs1.get_ndim() == 2) {
117  // cas 2 sh
118  s_one_plus = valeurs1(0,0) ;
119  s_one_minus = valeurs1(0,1) ;
120  ds_one_plus = valeurs1(0,2) ;
121  ds_one_minus = valeurs1(0,3) ;
122 
123  s_two_plus = valeurs2(1,0) ;
124  s_two_minus = valeurs2(1,1) ;
125  ds_two_plus = valeurs2(1,2) ;
126  ds_two_minus = valeurs2(1,3) ;
127  }
128  else {
129  // cas 1 sh :
130  Tbl& valeurs = (base_r == R_CHEBU) ? valeurs2 : valeurs1 ;
131  s_one_plus = valeurs(0) ;
132  s_one_minus = valeurs(1) ;
133  ds_one_plus = valeurs(2) ;
134  ds_one_minus = valeurs(3) ;
135  }
136 
137  return sh_ptens_rr(nr, l_quant, beta/alpha, base_r) ;
138 }
139 
140 }
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 rr component of the divergence-free tensor Poisson equation.
virtual void do_ope_cl() const
Computes the banded-matrix of the operator.
virtual Tbl get_solh() const
Computes the homogeneous solutions(s).
virtual void do_ope_mat() const
Computes the matrix of the operator.
virtual void do_non_dege() const
Computes the non-degenerated matrix of the operator.
virtual ~Ope_pois_tens_rr()
Destructor.
Ope_pois_tens_rr(int nbr, int baser, double alf, double bet, int lq, int dz)
Standard constructor.
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