LORENE
comb_lin_cpt.C
1 /*
2  * Copyright (c) 2000-2001 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 as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * LORENE is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with LORENE; if not, write to the Free Software
18  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19  *
20  */
21 
22 
23 
24 
25 /*
26  * $Id: comb_lin_cpt.C,v 1.5 2016/12/05 16:18:09 j_novak Exp $
27  * $Log: comb_lin_cpt.C,v $
28  * Revision 1.5 2016/12/05 16:18:09 j_novak
29  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
30  *
31  * Revision 1.4 2014/10/13 08:53:28 j_novak
32  * Lorene classes and functions now belong to the namespace Lorene.
33  *
34  * Revision 1.3 2014/10/06 15:16:08 j_novak
35  * Modified #include directives to use c++ syntax.
36  *
37  * Revision 1.2 2002/10/16 14:37:11 j_novak
38  * Reorganization of #include instructions of standard C++, in order to
39  * use experimental version 3 of gcc.
40  *
41  * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
42  * LORENE
43  *
44  * Revision 2.1 2000/11/22 19:29:50 eric
45  * Changement de nom_C en comb_lin_cpt_C
46  * Nettoyage des includes
47  *
48  * Revision 2.0 2000/03/16 16:25:18 phil
49  * *** empty log message ***
50  *
51  *
52  * $Header: /cvsroot/Lorene/C++/Source/Non_class_members/PDE/comb_lin_cpt.C,v 1.5 2016/12/05 16:18:09 j_novak Exp $
53  *
54  */
55 
56 // Headers C
57 #include <cstdlib>
58 
59 // Headers Lorene
60 #include "matrice.h"
61 
62 /*
63  * Gestion des CL permettant de mettre a bande les operateurs associes a poisson
64  * compact. Version pour les sources et les matrices des operateurs.
65  */
66 
67 
68 // Version Matrice --> Matrice
69 namespace Lorene {
70 Matrice _cl_cpt_pas_prevu (const Matrice &source, int) {
71  cout << "Combinaison lineaire pas prevu..." << endl ;
72  cout << "Source : " << source << endl ;
73  abort() ;
74  return source;
75 }
76 
77 
78  //-------------------
79  //-- R_CHEBP -----
80  //-------------------
81 
82 
83 Matrice _cl_cpt_r_chebp (const Matrice &source, int) {
84 
85  int n = source.get_dim(0) ;
86  assert (n == source.get_dim(1)) ;
87 
88  Matrice barre(source) ;
89  int dirac = 1 ;
90  for (int i=0 ; i<n-2 ; i++) {
91  for (int j=0 ; j<n ; j++)
92  barre.set(i, j) = ((1+dirac)*source(i, j)-source(i+2, j))/(i+1) ;
93  if (i==0) dirac = 0 ;
94  }
95 
96  Matrice res(barre) ;
97  for (int i=0 ; i<n-4 ; i++)
98  for (int j=0 ; j<n ; j++)
99  res.set(i, j) = barre(i, j)-barre(i+2, j) ;
100 
101  res.set_band(4, 1) ;
102  res.set_lu() ;
103  return res ;
104 }
105 
106  //-------------------
107  //-- R_CHEBI -----
108  //-------------------
109 
110 
111 Matrice _cl_cpt_r_chebi (const Matrice &source, int l) {
112  int n = source.get_dim(0) ;
113  assert (n == source.get_dim(1)) ;
114 
115  Matrice barre(source) ;
116  for (int i=0 ; i<n-2 ; i++)
117  for (int j=0 ; j<n ; j++)
118  barre.set(i, j) = (source(i, j)-source(i+1, j))/(i+1) ;
119 
120  Matrice res(barre) ;
121  for (int i=0 ; i<n-4 ; i++)
122  for (int j=0 ; j<n ; j++)
123  res.set(i, j) = barre(i, j)-barre(i+2, j) ;
124 
125  if (l==1)
126  res.set_band(3, 0) ;
127  else
128  res.set_band(3, 1) ;
129  res.set_lu() ;
130  return res ;
131 
132 }
133 
134 
135  //-------------------------
136  //- La routine a appeler ---
137  //---------------------------
138 
139 Matrice combinaison_cpt (const Matrice &source, int l, int base_r) {
140 
141  // Routines de derivation
142  static Matrice (*combinaison_cpt[MAX_BASE])
143  (const Matrice &, int) ;
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  combinaison_cpt[i] = _cl_cpt_pas_prevu ;
151  }
152  // Les routines existantes
153  combinaison_cpt[R_CHEBP >> TRA_R] = _cl_cpt_r_chebp ;
154  combinaison_cpt[R_CHEBI >> TRA_R] = _cl_cpt_r_chebi ;
155  }
156 
157  Matrice res(combinaison_cpt[base_r](source, l)) ;
158  return res ;
159 }
160 
161  //--------------------------------------------------------------
162  // Version Tbl
163  //--------------------------------------------------------------
164 
165 
166 
167 Tbl _cl_cpt_pas_prevu(const Tbl& tb) {
168  cout << "combinaison_nul_pas_prevu " << endl ;
169  cout << "tb : " << tb << endl ;
170  abort() ;
171  return tb ;
172 }
173 
174 
175  //-------------------
176  //-- R_CHEBP -----
177  //-------------------
178 Tbl _cl_cpt_r_chebp(const Tbl& tb) {
179 
180  assert (tb.get_etat() != ETATNONDEF) ;
181  int n=tb.get_dim(0) ;
182 
183  Tbl barre(tb) ;
184  int dirac = 1 ;
185  for (int i=0 ; i<n-2 ; i++) {
186  barre.set(i) = ((1+dirac)*tb(i)-tb(i+2))/(i+1) ;
187  if (i==0) dirac = 0 ;
188  }
189 
190  Tbl res(barre) ;
191  for (int i=0 ; i<n-4 ; i++)
192  res.set(i) = barre(i)-barre(i+2) ;
193 
194  return res ;
195 }
196 
197 
198 
199  //-------------------
200  //-- R_CHEBI -----
201  //-------------------
202 Tbl _cl_cpt_r_chebi(const Tbl& tb) {
203 
204  assert (tb.get_etat() != ETATNONDEF) ;
205  int n=tb.get_dim(0) ;
206 
207  Tbl barre(tb) ;
208  for (int i=0 ; i<n-2 ; i++)
209  barre.set(i) = (tb(i)-tb(i+1))/(i+1) ;
210 
211  Tbl res(barre) ;
212  for (int i=0 ; i<n-4 ; i++)
213  res.set(i) = barre(i)-barre(i+2) ;
214 
215  return res ;
216 }
217 
218 
219  //----------------------------
220  //- Routine a appeler ---
221  //------------------------------
222 
223 Tbl combinaison_cpt (const Tbl &source, int base_r) {
224 
225  // Routines de derivation
226  static Tbl (*combinaison_cpt[MAX_BASE])(const Tbl&) ;
227  static int nap = 0 ;
228 
229  // Premier appel
230  if (nap==0) {
231  nap = 1 ;
232  for (int i=0 ; i<MAX_BASE ; i++) {
233  combinaison_cpt[i] = _cl_cpt_pas_prevu ;
234  }
235  // Les routines existantes
236  combinaison_cpt[R_CHEBP >> TRA_R] = _cl_cpt_r_chebp ;
237  combinaison_cpt[R_CHEBI >> TRA_R] = _cl_cpt_r_chebi ;
238  }
239 
240  Tbl res(combinaison_cpt[base_r](source)) ;
241  return res ;
242 }
243 
244 }
Lorene prototypes.
Definition: app_hor.h:67
#define TRA_R
Translation en R, used for a bitwise shift (in hex)
Definition: type_parite.h:158
#define R_CHEBI
base de Cheb. impaire (rare) seulement
Definition: type_parite.h:170
#define R_CHEBP
base de Cheb. paire (rare) seulement
Definition: type_parite.h:168
int get_dim(int i) const
Returns the dimension of the matrix.
Definition: matrice.C:263
#define MAX_BASE
Nombre max. de bases differentes.
Definition: type_parite.h:144