LORENE
ope_vorton_solh.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
21char ope_vorton_solh_C[] = "$Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_vorton/ope_vorton_solh.C,v 1.5 2014/10/13 08:53:37 j_novak Exp $" ;
22
23/*
24 * $Id: ope_vorton_solh.C,v 1.5 2014/10/13 08:53:37 j_novak Exp $
25 * $Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_vorton/ope_vorton_solh.C,v 1.5 2014/10/13 08:53:37 j_novak Exp $
26 *
27 */
28#include <cmath>
29#include <cstdlib>
30
31#include "proto.h"
32#include "ope_elementary.h"
33
34 //------------------------------------
35 // Routine pour les cas non prevus --
36 //------------------------------------
37namespace Lorene {
38Tbl _solh_vorton_pas_prevu (int, int, double, double, Tbl&) {
39
40 cout << " Solution homogene pas prevue ..... : "<< endl ;
41 exit(-1) ;
42 Tbl res(1) ;
43 return res;
44}
45
46 //-------------------
47 //-- R_CHEBU -----
48 //-------------------
49
50Tbl _solh_vorton_r_chebu (int n, int l, double alpha, double, Tbl& val_lim) {
51
52 double l_one = -double(l) ;
53 double rminus = -0.5/alpha ;
54
55 Tbl res(n) ;
56 res.set_etat_qcq() ;
57 double* coloc = new double[n] ;
58
59 int * deg = new int[3] ;
60 deg[0] = 1 ;
61 deg[1] = 1 ;
62 deg[2] = n ;
63
64 //Construction de la premiere solution homogene :
65 for (int i=0 ; i<n ; i++)
66 coloc[i] = pow(1./alpha/(-cos(M_PI*i/(n-1))-1) , l_one) ;
67
68 cfrcheb(deg, deg, coloc, deg, coloc) ;
69 for (int i=0 ; i<n ;i++)
70 res.set(i) = coloc[i] ;
71
72 delete [] coloc ;
73 delete [] deg ;
74
75 val_lim.set(0,0) = pow(rminus, l_one) ;
76 val_lim.set(0,1) = l_one*pow(rminus, l_one-1.) ;
77 val_lim.set(0,2) = 0. ;
78 val_lim.set(0,3) = 0. ;
79 val_lim /= sqrt(double(2)) ;
80
81 return res ;
82}
83
84
85 //-------------------
86 //-- R_CHEB ------
87 //-------------------
88
89Tbl _solh_vorton_r_cheb (int n, int l, double alpha, double beta, Tbl& val_lim) {
90
91
92 double l_one = double(l+1) ;
93 double l_two = double(-l) ;
94 double rminus = beta - alpha ;
95 double rplus = beta + alpha ;
96
97 Tbl res(2, n) ;
98 res.set_etat_qcq() ;
99 double* coloc = new double[n] ;
100
101 int * deg = new int[3] ;
102 deg[0] = 1 ;
103 deg[1] = 1 ;
104 deg[2] = n ;
105
106 //Construction de la premiere solution homogene :
107 for (int i=0 ; i<n ; i++)
108 coloc[i] = pow(alpha*(-cos(M_PI*i/(n-1))) + beta, l_one) ;
109
110 cfrcheb(deg, deg, coloc, deg, coloc) ;
111 for (int i=0 ; i<n ;i++)
112 res.set(0, i) = coloc[i] ;
113
114 // construction de la seconde solution homogene :
115 for (int i=0 ; i<n ; i++)
116 coloc[i] = pow(alpha*(-cos(M_PI*i/(n-1))) + beta, l_two) ;
117
118 cfrcheb(deg, deg, coloc, deg, coloc) ;
119 for (int i=0 ; i<n ;i++)
120 res.set(1, i) = coloc[i] ;
121
122 delete [] coloc ;
123 delete [] deg ;
124
125 val_lim.set(0,0) = pow(rminus, l_one) ;
126 val_lim.set(0,1) = l_one*pow(rminus, l_one-1) ;
127 val_lim.set(0,2) = pow(rplus, l_one) ;
128 val_lim.set(0,3) = l_one*pow(rplus, l_one-1) ;
129
130 val_lim.set(1,0) = pow(rminus, l_two) ;
131 val_lim.set(1,1) = l_two*pow(rminus, l_two-1) ;
132 val_lim.set(1,2) = pow(rplus, l_two) ;
133 val_lim.set(1,3) = l_two*pow(rplus, l_two-1) ;
134 val_lim /= sqrt(double(2)) ;
135
136 return res ;
137}
138
139
141
142 // Routines de derivation
143 static Tbl (*solh_vorton[MAX_BASE]) (int, int, double, double, Tbl&) ;
144 static int nap = 0 ;
145
146 // Premier appel
147 if (nap==0) {
148 nap = 1 ;
149 for (int i=0 ; i<MAX_BASE ; i++) {
150 solh_vorton[i] = _solh_vorton_pas_prevu ;
151 }
152 // Les routines existantes
153 solh_vorton[R_CHEB >> TRA_R] = _solh_vorton_r_cheb ;
154 solh_vorton[R_CHEBU >> TRA_R] = _solh_vorton_r_chebu ;
155 }
156
157 Tbl val_lim (2 ,4) ;
158 val_lim.set_etat_qcq() ;
160
161
162 s_one_minus = val_lim(0,0) ;
163 ds_one_minus = val_lim(0,1) ;
164 s_one_plus = val_lim(0,2) ;
165 ds_one_plus = val_lim(0,3) ;
166
167 if (res.get_ndim()>1) {
168 s_two_minus = val_lim(1,0) ;
169 ds_two_minus = val_lim(1,1) ;
170 s_two_plus = val_lim(1,2) ;
171 ds_two_plus = val_lim(1,3) ;
172 }
173
174 return res ;
175}
176}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
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.
double ds_two_plus
Value of the derivative of the second homogeneous solution at the outer boundary.
int nr
Number of radial points.
int l_quant
quantum number
virtual Tbl get_solh() const
Computes the homogeneous solutions(s).
Basic array class.
Definition tbl.h:161
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Cmp cos(const Cmp &)
Cosine.
Definition cmp_math.C:94
#define MAX_BASE
Nombre max. de bases differentes.
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
#define TRA_R
Translation en R, used for a bitwise shift (in hex)
#define R_CHEB
base de Chebychev ordinaire (fin)
Lorene prototypes.
Definition app_hor.h:64