28char ope_pvect_r_mat_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Ope_elementary/Ope_pois_vect_r/ope_pvect_r_mat.C,v 1.3 2014/10/13 08:53:34 j_novak Exp $" ;
52#include "type_parite.h"
76Matrice _ope_pvect_r_mat_pas_prevu(
int,
int,
double,
int) ;
77Matrice _ope_pvect_r_mat_r_chebp(
int,
int,
double,
int) ;
78Matrice _ope_pvect_r_mat_r_chebi(
int,
int,
double,
int) ;
79Matrice _ope_pvect_r_mat_r_chebu(
int,
int,
double,
int) ;
80Matrice _ope_pvect_r_mat_r_cheb(
int,
int,
double,
int) ;
86Matrice _ope_pvect_r_mat_pas_prevu(
int n,
int l,
double echelle,
int puis) {
87 cout <<
"laplacien vect_r pas prevu..." << endl ;
88 cout <<
"n : " << n << endl ;
89 cout <<
"l : " << l << endl ;
90 cout <<
"puissance : " << puis << endl ;
91 cout <<
"echelle : " << echelle << endl ;
104Matrice _ope_pvect_r_mat_r_chebp (
int n,
int l,
double,
int) {
106 const int nmax = 100 ;
107 static Matrice* tab[nmax] ;
108 static int nb_dejafait = 0 ;
109 static int l_dejafait[nmax] ;
110 static int nr_dejafait[nmax] ;
115 for (
int conte=0 ; conte<nb_dejafait ; conte ++)
116 if ((l_dejafait[conte] == l) && (nr_dejafait[conte] == n))
121 if (nb_dejafait >= nmax) {
122 cout <<
"_ope_pvect_r_mat_r_chebp : trop de matrices" << endl ;
128 l_dejafait[nb_dejafait] = l ;
129 nr_dejafait[nb_dejafait] = n ;
138 double* vect =
new double[n] ;
140 for (
int i=0 ; i<n ; i++) {
141 for (
int j=0 ; j<n ; j++)
144 d2sdx2_1d (n, &vect,
R_CHEBP) ;
146 for (
int j=0 ; j<n ; j++)
147 dd.set(j, i) = vect[j] ;
150 for (
int i=0 ; i<n ; i++) {
151 for (
int j=0 ; j<n ; j++)
154 sxdsdx_1d (n, &vect,
R_CHEBP) ;
155 for (
int j=0 ; j<n ; j++)
156 xd.set(j, i) = vect[j] ;
160 for (
int i=0 ; i<n ; i++) {
161 for (
int j=0 ; j<n ; j++)
165 for (
int j=0 ; j<n ; j++)
166 xx.set(j, i) = vect[j] ;
173 res = dd+4*xd+(2-l*(l+1))*xx ;
175 res = dd + 2*xd - 2*xx ;
176 tab[nb_dejafait] =
new Matrice(res) ;
183 return *tab[indice] ;
193Matrice _ope_pvect_r_mat_r_chebi (
int n,
int l,
double,
int) {
195 const int nmax = 100 ;
196 static Matrice* tab[nmax] ;
197 static int nb_dejafait = 0 ;
198 static int l_dejafait[nmax] ;
199 static int nr_dejafait[nmax] ;
204 for (
int conte=0 ; conte<nb_dejafait ; conte ++)
205 if ((l_dejafait[conte] == l) && (nr_dejafait[conte] == n))
210 if (nb_dejafait >= nmax) {
211 cout <<
"_ope_pvect_r_mat_r_chebi : trop de matrices" << endl ;
217 l_dejafait[nb_dejafait] = l ;
218 nr_dejafait[nb_dejafait] = n ;
227 double* vect =
new double[n] ;
229 for (
int i=0 ; i<n ; i++) {
230 for (
int j=0 ; j<n ; j++)
233 d2sdx2_1d (n, &vect,
R_CHEBI) ;
234 for (
int j=0 ; j<n ; j++)
235 dd.set(j, i) = vect[j] ;
238 for (
int i=0 ; i<n ; i++) {
239 for (
int j=0 ; j<n ; j++)
242 sxdsdx_1d (n, &vect,
R_CHEBI) ;
243 for (
int j=0 ; j<n ; j++)
244 xd.set(j, i) = vect[j] ;
247 for (
int i=0 ; i<n ; i++) {
248 for (
int j=0 ; j<n ; j++)
252 for (
int j=0 ; j<n ; j++)
253 xx.set(j, i) = vect[j] ;
260 res = dd+4*xd+(2-l*(l+1))*xx ;
262 res = dd + 2*xd - 2*xx ;
263 tab[nb_dejafait] =
new Matrice(res) ;
270 return *tab[indice] ;
280Matrice _ope_pvect_r_mat_r_chebu(
int n,
int l,
double,
int puis) {
283 cout <<
"_ope_pvect_r_mat_r_chebu : only the case dzpuis = 4 "
284 <<
'\n' <<
"is implemented! \n"
285 <<
"dzpuis = " << puis << endl ;
289 const int nmax = 200 ;
290 static Matrice* tab[nmax] ;
291 static int nb_dejafait = 0 ;
292 static int l_dejafait[nmax] ;
293 static int nr_dejafait[nmax] ;
298 for (
int conte=0 ; conte<nb_dejafait ; conte ++)
299 if ((l_dejafait[conte] == l) && (nr_dejafait[conte] == n))
304 if (nb_dejafait >= nmax) {
305 cout <<
"_ope_pvect_r_mat_r_chebu : trop de matrices" << endl ;
310 l_dejafait[nb_dejafait] = l ;
311 nr_dejafait[nb_dejafait] = n ;
320 double* vect =
new double[n] ;
322 for (
int i=0 ; i<n ; i++) {
323 for (
int j=0 ; j<n ; j++)
326 d2sdx2_1d (n, &vect,
R_CHEBU) ;
327 for (
int j=0 ; j<n ; j++)
328 dd.set(j, i) = vect[j] ;
331 for (
int i=0 ; i<n ; i++) {
332 for (
int j=0 ; j<n ; j++)
336 sxm1_1d_cheb (n, vect) ;
337 for (
int j=0 ; j<n ; j++)
338 xd.set(j, i) = vect[j] ;
341 for (
int i=0 ; i<n ; i++) {
342 for (
int j=0 ; j<n ; j++)
346 for (
int j=0 ; j<n ; j++)
347 xx.set(j, i) = vect[j] ;
356 res = dd - 2*xd + (2 -l*(l+1))*xx ;
357 tab[nb_dejafait] =
new Matrice(res) ;
364 return *tab[indice] ;
373Matrice _ope_pvect_r_mat_r_cheb (
int n,
int l,
double echelle,
int) {
375 const int nmax = 200 ;
376 static Matrice* tab[nmax] ;
377 static int nb_dejafait = 0 ;
378 static int l_dejafait[nmax] ;
379 static int nr_dejafait[nmax] ;
380 static double vieux_echelle = 0;
383 if (vieux_echelle != echelle) {
384 for (
int i=0 ; i<nb_dejafait ; i++) {
386 nr_dejafait[i] = -1 ;
391 vieux_echelle = echelle ;
397 for (
int conte=0 ; conte<nb_dejafait ; conte ++)
398 if ((l_dejafait[conte] == l) && (nr_dejafait[conte] == n))
403 if (nb_dejafait >= nmax) {
404 cout <<
"_ope_pvect_r_mat_r_cheb : trop de matrices" << endl ;
410 l_dejafait[nb_dejafait] = l ;
411 nr_dejafait[nb_dejafait] = n ;
420 double* vect =
new double[n] ;
422 for (
int i=0 ; i<n ; i++) {
423 for (
int j=0 ; j<n ; j++)
426 d2sdx2_1d (n, &vect,
R_CHEB) ;
427 for (
int j=0 ; j<n ; j++)
428 dd.set(j, i) = vect[j]*echelle*echelle ;
431 for (
int i=0 ; i<n ; i++) {
432 for (
int j=0 ; j<n ; j++)
435 d2sdx2_1d (n, &vect,
R_CHEB) ;
436 multx_1d (n, &vect,
R_CHEB) ;
437 for (
int j=0 ; j<(n>i+1 ? i+1 : n) ; j++)
438 dd.set(j, i) += 2*echelle*vect[j] ;
441 for (
int i=0 ; i<n ; i++) {
442 for (
int j=0 ; j<n ; j++)
445 d2sdx2_1d (n, &vect,
R_CHEB) ;
446 multx_1d (n, &vect,
R_CHEB) ;
447 multx_1d (n, &vect,
R_CHEB) ;
448 for (
int j=0 ; j<(n>i+1 ? i+1 : n) ; j++)
449 dd.set(j, i) += vect[j] ;
452 for (
int i=0 ; i<n ; i++) {
453 for (
int j=0 ; j<n ; j++)
456 sxdsdx_1d (n, &vect,
R_CHEB) ;
457 for (
int j=0 ; j<n ; j++)
458 xd.set(j, i) = vect[j]*echelle ;
461 for (
int i=0 ; i<n ; i++) {
462 for (
int j=0 ; j<n ; j++)
465 sxdsdx_1d (n, &vect,
R_CHEB) ;
466 multx_1d (n, &vect,
R_CHEB) ;
467 for (
int j=0 ; j<(n>i+1 ? i+1 : n) ; j++)
468 xd.set(j, i) += vect[j] ;
471 for (
int i=0 ; i<n ; i++) {
472 for (
int j=0 ; j<n ; j++)
475 sx2_1d (n, &vect,
R_CHEB) ;
476 for (
int j=0 ; j<n ; j++)
477 xx.set(j, i) = vect[j] ;
486 res = dd + 4*xd + (2 - l*(l+1))*xx ;
487 tab[nb_dejafait] =
new Matrice(res) ;
494 return *tab[indice] ;
502Matrice ope_pvect_r_mat(
int n,
int l,
double echelle,
int puis,
int base_r)
506 static Matrice (*ope_pvect_r_mat[
MAX_BASE])(int, int, double, int) ;
513 ope_pvect_r_mat[i] = _ope_pvect_r_mat_pas_prevu ;
516 ope_pvect_r_mat[
R_CHEB >>
TRA_R] = _ope_pvect_r_mat_r_cheb ;
517 ope_pvect_r_mat[
R_CHEBU >>
TRA_R] = _ope_pvect_r_mat_r_chebu ;
518 ope_pvect_r_mat[
R_CHEBP >>
TRA_R] = _ope_pvect_r_mat_r_chebp ;
519 ope_pvect_r_mat[
R_CHEBI >>
TRA_R] = _ope_pvect_r_mat_r_chebi ;
522 Matrice res(ope_pvect_r_mat[base_r](n, l, echelle, puis)) ;
#define MAX_BASE
Nombre max. de bases differentes.
#define R_CHEBU
base de Chebychev ordinaire (fin), dev. en 1/r
#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