LORENE
ope_poisson_pseudo_1d_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_poisson_pseudo_1d_solh_C[] = "$Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_poisson_pseudo_1d/ope_poisson_pseudo_1d_solh.C,v 1.3 2014/10/13 08:53:35 j_novak Exp $" ;
22
23/*
24 * $Id: ope_poisson_pseudo_1d_solh.C,v 1.3 2014/10/13 08:53:35 j_novak Exp $
25 * $Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_poisson_pseudo_1d/ope_poisson_pseudo_1d_solh.C,v 1.3 2014/10/13 08:53:35 j_novak Exp $
26 *
27 */
28#include <cmath>
29#include <cstdlib>
30
31#include "proto.h"
32#include "ope_elementary.h"
33
34
35 //------------------------------------
36 // Routine pour les cas non prevus --
37 //------------------------------------
38namespace Lorene {
39Tbl _solh_poisson_pseudo_1d_pas_prevu (int, int,double, double, Tbl&) {
40
41 cout << " Solution homogene pas prevue ..... : "<< endl ;
42 exit(-1) ;
43 Tbl res(1) ;
44 return res;
45}
46
47
48 //-------------------
49 //-- R_CHEB ------
50 //-------------------
51
52Tbl _solh_poisson_pseudo_1d_r_cheb (int n, int l, double alpha, double beta,
53 Tbl& val_lim) {
54
55
56 double echelle = beta / alpha ;
57 int expo_un = l ;
58 int expo_deux = -l+1 ;
59
60 val_lim.set(0,0) = pow(echelle-1, double(expo_un)) ;
61 val_lim.set(0,1) = double(expo_un) * pow(echelle-1, double(expo_un-1))/alpha ;
62 val_lim.set(0,2) = pow(echelle+1, double(expo_un)) ;
63 val_lim.set(0,3) = double(expo_un) * pow(echelle+1, double(expo_un-1))/alpha ;
64
65
66 val_lim.set(1,0) = pow(echelle-1, double(expo_deux)) ;
67 val_lim.set(1,1) = double(expo_deux) * pow(echelle-1, double(expo_deux-1))/alpha ;
68 val_lim.set(1,2) = pow(echelle+1, double(expo_deux)) ;
69 val_lim.set(1,3) = double(expo_deux) * pow(echelle+1, double(expo_deux-1))/alpha ;
70
71
72 Tbl res(2, n) ;
73 res.set_etat_qcq() ;
74 double* coloc = new double[n] ;
75
76 int * deg = new int[3] ;
77 deg[0] = 1 ;
78 deg[1] = 1 ;
79 deg[2] = n ;
80
81 //Construction de la premiere solution homogene :
82
83 for (int i=0 ; i<n ; i++)
84 coloc[i] = pow(echelle-cos(M_PI*i/(n-1)), double(expo_un)) ;
85
86 cfrcheb(deg, deg, coloc, deg, coloc) ;
87 for (int i=0 ; i<n ;i++)
88 res.set(0, i) = coloc[i] ;
89
90 // construction de la seconde solution homogene :
91 for (int i=0 ; i<n ; i++)
92 coloc[i] = pow(echelle-cos(M_PI*i/(n-1)), double(expo_deux)) ;
93
94 cfrcheb(deg, deg, coloc, deg, coloc) ;
95 for (int i=0 ; i<n ;i++)
96 res.set(1, i) = coloc[i] ;
97
98
99 delete [] coloc ;
100 delete [] deg ;
101
102 return res ;
103}
104
105
106 //-------------------
107 //-- R_CHEBP ------
108 //-------------------
109
110Tbl _solh_poisson_pseudo_1d_r_chebp (int n, int l, double alpha,
111 double, Tbl& val_lim) {
112
113
114 val_lim.set(0,0) = (l!=0) ? 1 : 0 ;
115 val_lim.set(0,1) = (l!=1) ? 0 : 1 ;
116 val_lim.set(0,2) = 1. ;
117 val_lim.set(0,3) = double(l)/alpha ;
118
119
120 assert (div(l, 2).rem ==0) ;
121
122 Tbl res(n) ;
123 res.set_etat_qcq() ;
124 double* coloc = new double[n] ;
125
126 int * deg = new int[3] ;
127 deg[0] = 1 ;
128 deg[1] = 1 ;
129 deg[2] = n ;
130
131 for (int i=0 ; i<n ; i++)
132 coloc[i] = pow(sin(M_PI*i/2/(n-1)), double(l)) ;
133
134 cfrchebp(deg, deg, coloc, deg, coloc) ;
135 for (int i=0 ; i<n ;i++)
136 res.set(i) = coloc[i] ;
137
138 delete [] coloc ;
139 delete [] deg ;
140
141 return res ;
142}
143
144 //-------------------
145 //-- R_CHEBI -----
146 //-------------------
147
148Tbl _solh_poisson_pseudo_1d_r_chebi (int n, int l,
149 double alpha, double, Tbl& val_lim) {
150
151 val_lim.set(0,0) = 0 ;
152 val_lim.set(0,1) = (l!=1) ? 0 : 1 ;
153 val_lim.set(0,2) = 1. ;
154 val_lim.set(0,3) = double(l)/alpha ;
155
156
157
158 assert (div(l, 2).rem == 1) ;
159
160 Tbl res(n) ;
161 res.set_etat_qcq() ;
162 double* coloc = new double[n] ;
163
164 int * deg = new int[3] ;
165 deg[0] = 1 ;
166 deg[1] = 1 ;
167 deg[2] = n ;
168
169 for (int i=0 ; i<n ; i++)
170 coloc[i] = pow(sin(M_PI*i/2/(n-1)), double(l)) ;
171
172 cfrchebi(deg, deg, coloc, deg, coloc) ;
173 for (int i=0 ; i<n ;i++)
174 res.set(i) = coloc[i] ;
175
176 delete [] coloc ;
177 delete [] deg ;
178
179 return res ;
180}
181
182
184
185 // Routines de derivation
187 static int nap = 0 ;
188
189 // Premier appel
190 if (nap==0) {
191 nap = 1 ;
192 for (int i=0 ; i<MAX_BASE ; i++) {
193 solh_poisson_pseudo_1d[i] = _solh_poisson_pseudo_1d_pas_prevu ;
194 }
195 // Les routines existantes
196 solh_poisson_pseudo_1d[R_CHEB >> TRA_R] = _solh_poisson_pseudo_1d_r_cheb ;
197 solh_poisson_pseudo_1d[R_CHEBP >> TRA_R] = _solh_poisson_pseudo_1d_r_chebp ;
198 solh_poisson_pseudo_1d[R_CHEBI >> TRA_R] = _solh_poisson_pseudo_1d_r_chebi ;
199 }
200
201 Tbl val_lim (2,4) ;
202 val_lim.set_etat_qcq() ;
204
205 s_one_minus = val_lim(0,0) ;
206 ds_one_minus = val_lim(0,1) ;
207 s_one_plus = val_lim(0,2) ;
208 ds_one_plus = val_lim(0,3) ;
209
210 s_two_minus = val_lim(1,0) ;
211 ds_two_minus = val_lim(1,1) ;
212 s_two_plus = val_lim(1,2) ;
213 ds_two_plus = val_lim(1,3) ;
214
215 return res ;
216}
217}
Tbl & set(int l)
Read/write of the value in a given domain.
Definition cmp.h:724
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.
virtual Tbl get_solh() const
Computes the homogeneous solutions(s).
Basic array class.
Definition tbl.h:161
double & set(int i)
Read/write of a particular element (index i) (1D case)
Definition tbl.h:281
Cmp sin(const Cmp &)
Sine.
Definition cmp_math.C:69
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_CHEBI
base de Cheb. impaire (rare) seulement
#define TRA_R
Translation en R, used for a bitwise shift (in hex)
#define R_CHEB
base de Chebychev ordinaire (fin)
#define R_CHEBP
base de Cheb. paire (rare) seulement
Lorene prototypes.
Definition app_hor.h:64