23char solp_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Non_class_members/PDE/solp.C,v 1.11 2014/10/13 08:53:31 j_novak Exp $" ;
120#include "type_parite.h"
141Tbl _solp_pas_prevu (
const Matrice &lap,
const Matrice &nondege,
double alpha,
142 double beta,
const Tbl &source,
int puis) {
143 cout <<
" Solution particuliere pas prevue ..... : "<< endl ;
144 cout <<
" Laplacien : " << endl << lap << endl ;
145 cout <<
" Non degenere : " << endl << nondege << endl ;
146 cout <<
" Source : " << &source << endl ;
147 cout <<
" Alpha : " << alpha << endl ;
148 cout <<
" Beta : " << beta << endl ;
149 cout <<
" Puiss : " << puis << endl ;
161Tbl _solp_r_cheb (
const Matrice &lap,
const Matrice &nondege,
double alpha,
162 double beta,
const Tbl &source,
int) {
165 int dege = n-nondege.get_dim(0) ;
168 Tbl source_aux(source*alpha*alpha) ;
169 Tbl xso(source_aux) ;
170 Tbl xxso(source_aux) ;
171 multx_1d(n, &xso.t,
R_CHEB) ;
172 multx_1d(n, &xxso.t,
R_CHEB) ;
173 multx_1d(n, &xxso.t,
R_CHEB) ;
174 source_aux = beta*beta/alpha/alpha*source_aux+2*beta/alpha*xso+xxso ;
175 source_aux = combinaison(source_aux, 0,
R_CHEB) ;
179 for (
int i=0 ; i<n-dege ; i++)
180 so.set(i) = source_aux(i) ;
182 Tbl auxi(nondege.inverse(so)) ;
186 for (
int i=dege ; i<n ; i++)
187 res.set(i) = auxi(i-dege) ;
189 for (
int i=0 ; i<dege ; i++)
198Tbl _solp_r_jaco02 (
const Matrice &lap,
const Matrice &nondege,
double alpha,
199 double,
const Tbl &source,
int) {
202 int dege = n-nondege.get_dim(0) ;
205 Tbl source_aux(source*alpha*alpha) ;
206 source_aux = combinaison(source_aux, 0,
R_JACO02) ;
210 for (
int i=0 ; i<n-dege ; i++)
211 so.set(i) = source_aux(i) ;
213 Tbl auxi(nondege.inverse(so)) ;
217 for (
int i=dege ; i<n ; i++)
218 res.set(i) = auxi(i-dege) ;
220 for (
int i=0 ; i<dege ; i++)
230Tbl _solp_r_chebp (
const Matrice &lap,
const Matrice &nondege,
double alpha,
231 double ,
const Tbl &source,
int) {
234 int dege = n-nondege.get_dim(0) ;
235 assert ((dege==2) || (dege == 1)) ;
236 Tbl source_aux(alpha*alpha*source) ;
237 source_aux = combinaison(source_aux, 0,
R_CHEBP) ;
241 for (
int i=0 ; i<n-dege ; i++)
242 so.set(i) = source_aux(i);
244 Tbl auxi(nondege.inverse(so)) ;
248 for (
int i=dege ; i<n ; i++)
249 res.set(i) = auxi(i-dege) ;
251 for (
int i=0 ; i<dege ; i++)
256 for (
int i=0 ; i<n ; i++)
259 else somme += res(i) ;
271Tbl _solp_r_chebi (
const Matrice &lap,
const Matrice &nondege,
double alpha,
272 double,
const Tbl &source,
int) {
276 int dege = n-nondege.get_dim(0) ;
277 assert ((dege == 2) || (dege ==1)) ;
278 Tbl source_aux(source*alpha*alpha) ;
279 source_aux = combinaison(source_aux, 0,
R_CHEBI) ;
283 for (
int i=0 ; i<n-dege ; i++)
284 so.set(i) = source_aux(i);
286 Tbl auxi(nondege.inverse(so)) ;
290 for (
int i=dege ; i<n ; i++)
291 res.set(i) = auxi(i-dege) ;
293 for (
int i=0 ; i<dege ; i++)
298 for (
int i=0 ; i<n ; i++)
300 somme -= (2*i+1)*res(i) ;
301 else somme += (2*i+1)*res(i) ;
314Tbl _solp_r_chebu (
const Matrice &lap,
const Matrice &nondege,
double alpha,
315 double,
const Tbl &source,
int puis) {
322 res = _solp_r_chebu_cinq (lap, nondege, source) ;
325 res = _solp_r_chebu_quatre (lap, nondege, alpha, source) ;
328 res = _solp_r_chebu_trois (lap, nondege, alpha, source) ;
331 res = _solp_r_chebu_deux (lap, nondege, source) ;
342Tbl _solp_r_chebu_quatre (
const Matrice &lap,
const Matrice &nondege,
double alpha,
346 int dege = n-nondege.get_dim(0) ;
347 assert ((dege==3) || (dege ==2)) ;
348 Tbl source_aux(source*alpha*alpha) ;
349 source_aux = combinaison(source_aux, 4,
R_CHEBU) ;
353 for (
int i=0 ; i<n-dege ; i++)
354 so.set(i) = source_aux(i);
356 Tbl auxi(nondege.inverse(so)) ;
360 for (
int i=dege ; i<n ; i++)
361 res.set(i) = auxi(i-dege) ;
363 for (
int i=0 ; i<dege ; i++)
368 for (
int i=0 ; i<n ; i++)
369 somme += i*i*res(i) ;
370 double somme_deux = somme ;
371 for (
int i=0 ; i<n ; i++)
372 somme_deux -= res(i) ;
373 res.set(1) = -somme ;
374 res.set(0) = somme_deux ;
379 for (
int i=0 ; i<n ; i++)
381 res.set(0) = -somme ;
387Tbl _solp_r_chebu_trois (
const Matrice &lap,
const Matrice &nondege,
double alpha,
391 int dege = n-nondege.get_dim(0) ;
394 Tbl source_aux(source*alpha) ;
395 source_aux = combinaison(source_aux, 3,
R_CHEBU) ;
399 for (
int i=0 ; i<n-dege ; i++)
400 so.set(i) = source_aux(i);
402 Tbl auxi(nondege.inverse(so)) ;
406 for (
int i=dege ; i<n ; i++)
407 res.set(i) = auxi(i-dege) ;
409 for (
int i=0 ; i<dege ; i++)
413 for (
int i=0 ; i<n ; i++)
415 res.set(0) = -somme ;
421Tbl _solp_r_chebu_deux (
const Matrice &lap,
const Matrice &nondege,
425 int dege = n-nondege.get_dim(0) ;
426 assert ((dege==1) || (dege ==2)) ;
427 Tbl source_aux(combinaison(source, 2,
R_CHEBU)) ;
431 for (
int i=0 ; i<n-dege ; i++)
432 so.set(i) = source_aux(i);
434 Tbl auxi(nondege.inverse(so)) ;
438 for (
int i=dege ; i<n ; i++)
439 res.set(i) = auxi(i-dege) ;
441 for (
int i=0 ; i<dege ; i++)
446 for (
int i=0 ; i<n ; i++)
449 res.set(0) = -somme ;
456Tbl _solp_r_chebu_cinq (
const Matrice &lap,
const Matrice &nondege,
460 int dege = n-nondege.get_dim(0) ;
462 Tbl source_aux(combinaison(source, 5,
R_CHEBU)) ;
466 for (
int i=0 ; i<n-dege ; i++)
467 so.set(i) = source_aux(i);
469 Tbl auxi(nondege.inverse(so)) ;
473 for (
int i=dege ; i<n ; i++)
474 res.set(i) = auxi(i-dege) ;
476 for (
int i=0 ; i<dege ; i++)
481 for (
int i=0 ; i<n ; i++)
484 res.set(0) = -somme ;
496Tbl solp(
const Matrice &lap,
const Matrice &nondege,
double alpha,
497 double beta,
const Tbl &source,
int puis,
int base_r) {
500 static Tbl (*solp[
MAX_BASE])(
const Matrice&,
const Matrice&, double, double,
508 solp[i] = _solp_pas_prevu ;
518 return solp[base_r](lap, nondege, alpha, beta, source, puis) ;
int get_dim(int i) const
Gives the i-th dimension (ie dim.dim[i])
#define MAX_BASE
Nombre max. de bases differentes.
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
#define R_JACO02
base de Jacobi(0,2) ordinaire (finjac)
#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