1 #@ MODIF t_fonction Utilitai DATE 04/09/2007 AUTEUR DURAND C.DURAND
2 # -*- coding: iso-8859-1 -*-
3 # CONFIGURATION MANAGEMENT OF EDF VERSION
4 # ======================================================================
5 # COPYRIGHT (C) 1991 - 2005 EDF R&D WWW.CODE-ASTER.ORG
6 # THIS PROGRAM 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.
11 # THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
12 # WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
13 # MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
14 # GENERAL PUBLIC LICENSE FOR MORE DETAILS.
16 # YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
17 # ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
18 # 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
19 # ======================================================================
25 # -----------------------------------------------------------------------------
26 class FonctionError(Exception): pass
28 class ParametreError(FonctionError): pass # probleme de NOM_PARA
29 class InterpolationError(FonctionError): pass
30 class ProlongementError(FonctionError): pass
32 # -----------------------------------------------------------------------------
33 def interp(typ_i,val,x1,x2,y1,y2) :
34 """Interpolation linéaire/logarithmique entre un couple de valeurs
36 if typ_i==['LIN','LIN']: return y1+(y2-y1)*(val-x1)/(x2-x1)
37 if typ_i==['LIN','LOG']: return exp(log(y1)+(val-x1)*(log(y2)-log(y1))/(x2-x1))
38 if typ_i==['LOG','LOG']: return exp(log(y1)+(log(val)-log(x1))*(log(y2)-log(y1))/(log(x2)-log(x1)))
39 if typ_i==['LOG','LIN']: return y1+(log(val)-log(x1))*(y2-y1)/(log(x2)-log(x1))
41 if val==x1 : return y1
42 elif val==x2 : return y2
44 raise InterpolationError, "abscisse = %g, intervalle = [%g, %g]" % (val, x1, x2)
47 listb=list(Set(liste))
52 # -----------------------------------------------------------------------------
54 """Classe pour fonctions réelles, équivalent au type aster = fonction_sdaster
56 def __init__(self,vale_x,vale_y,para,nom='') :
57 """Création d'un objet fonction
58 - vale_x et vale_y sont des listes de réels de meme longueur
59 - para est un dictionnaire contenant les entrées PROL_DROITE, PROL_GAUCHE et INTERPOL (cf sd ASTER)
64 if pk!=['INTERPOL','NOM_PARA','NOM_RESU','PROL_DROITE','PROL_GAUCHE'] :
65 raise FonctionError, 'fonction : parametres incorrects'
66 if para['INTERPOL'] not in [['NON','NON'],['LIN','LIN'],['LIN','LOG'],['LOG','LOG'],['LOG','LIN'],] :
67 raise FonctionError, 'fonction : parametre INTERPOL incorrect'
68 if para['PROL_DROITE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
69 raise FonctionError, 'fonction : parametre PROL_DROITE incorrect'
70 if para['PROL_GAUCHE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
71 raise FonctionError, 'fonction : parametre PROL_GAUCHE incorrect'
72 self.vale_x = array(vale_x)
73 self.vale_y = array(vale_y)
75 if len(self.vale_x)!=len(self.vale_y) :
76 raise FonctionError, 'fonction : longueur abscisse <> longueur ordonnées'
77 if not is_ordo(self.vale_x) :
78 raise FonctionError, 'fonction : abscisses non strictement croissantes'
80 def __add__(self,other) :
81 """addition avec une autre fonction ou un nombre, par surcharge de l'opérateur +
83 if isinstance(other,t_fonction):
84 para=copy.copy(self.para)
85 vale_x,para['PROL_GAUCHE'],para['PROL_DROITE']=self.homo_support(other)
86 fff=self.evalfonc(vale_x)
87 ggg=other.evalfonc(vale_x)
88 if isinstance(self,t_fonction_c): return t_fonction_c(vale_x,fff.vale_y+ggg.vale_y,para)
89 else : return t_fonction(vale_x,fff.vale_y+ggg.vale_y,para)
90 elif type(other) in [types.FloatType,types.IntType,types.ComplexType] :
91 if isinstance(self,t_fonction_c): return t_fonction_c(self.vale_x,self.vale_y+other,self.para)
92 else : return t_fonction(self.vale_x,self.vale_y+other,self.para)
93 else: raise FonctionError, 'fonctions : erreur de type dans __add__'
95 def __mul__(self,other) :
96 """multiplication avec une autre fonction ou un nombre, par surcharge de l'opérateur *
98 if isinstance(other,t_fonction):
99 para=copy.copy(self.para)
100 vale_x,para['PROL_GAUCHE'],para['PROL_DROITE']=self.homo_support(other)
101 fff=self.evalfonc(vale_x)
102 ggg=other.evalfonc(vale_x)
103 if isinstance(self,t_fonction_c): return t_fonction_c(vale_x,fff.vale_y*ggg.vale_y,para)
104 else : return t_fonction(vale_x,fff.vale_y*ggg.vale_y,para)
105 elif type(other) in [types.FloatType,types.IntType] :
106 return t_fonction(self.vale_x,self.vale_y*other,self.para)
107 elif type(other) ==types.ComplexType :
108 return t_fonction_c(self.vale_x,self.vale_y*other,self.para)
109 else: raise FonctionError, 'fonctions : erreur de type dans __mul__'
112 """affichage de la fonction en double colonne
115 for i in range(len(self.vale_x)) :
116 texte.append('%f %f' % (self.vale_x[i],self.vale_y[i]))
117 return '\n'.join(texte)
119 def __getitem__(self,other) :
120 """composition de deux fonction F[G]=FoG=F(G(x))
122 para=copy.copy(self.para)
123 if other.para['NOM_RESU']!=self.para['NOM_PARA'] :
124 raise ParametreError,'''composition de fonctions : NOM_RESU1 et NOM_PARA2 incohérents '''
125 para['NOM_PARA']==other.para['NOM_PARA']
126 return t_fonction(other.vale_x,map(self,other.vale_y),para)
128 def __call__(self,val,tol=1.e-6):
129 """méthode pour évaluer f(x)
130 - tolérance, par défaut 1.e-6 en relatif sur la longueur de l'intervalle
131 - adjacent, pour capter les erreurs d'arrondi en cas de prolongement exclu
133 i=searchsorted(self.vale_x,val)
136 if self.para['PROL_GAUCHE']=='EXCLU' :
137 eps_g=(val-self.vale_x[0] )/(self.vale_x[1] -self.vale_x[0])
138 if abs(eps_g)<=tol : return self.vale_y[0]
139 else : raise ProlongementError, 'fonction évaluée hors du domaine de définition'
141 if self.para['PROL_GAUCHE']=='CONSTANT' : return self.vale_y[0]
142 if self.para['PROL_GAUCHE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val,self.vale_x[0],
147 if self.para['PROL_DROITE']=='EXCLU' :
148 eps_d=(val-self.vale_x[-1])/(self.vale_x[-1]-self.vale_x[-2])
149 if abs(eps_d)<=tol : return self.vale_y[-1]
150 else : raise ProlongementError, 'fonction évaluée hors du domaine de définition'
152 if self.para['PROL_DROITE']=='CONSTANT' : return self.vale_y[-1]
153 if self.para['PROL_DROITE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val,self.vale_x[-1],
158 return interp(self.para['INTERPOL'],val,self.vale_x[i-1],
163 def homo_support(self,other) :
164 """Renvoie le support d'abscisses homogénéisé entre self et other
165 i.e. si prolongement exclu, on retient plus grand min ou plus petit max, selon
166 si prolongement autorisé, on conserve les abscisses d'une fonction, extrapolantes
168 Pour les points intermédiaires : union et tri des valeurs des vale_x réunis.
170 if other.vale_x[0]>self.vale_x[0]:
171 if other.para['PROL_GAUCHE']!='EXCLU' : f_g=self
174 if self.para['PROL_GAUCHE'] !='EXCLU' : f_g=other
176 val_min =f_g.vale_x[0]
177 prol_gauche=f_g.para['PROL_GAUCHE']
178 if self.vale_x[-1]>other.vale_x[-1]:
179 if other.para['PROL_DROITE']!='EXCLU' : f_d=self
182 if self.para['PROL_DROITE'] !='EXCLU' : f_d=other
184 val_max =f_d.vale_x[-1]
185 prol_droite=f_d.para['PROL_DROITE']
186 vale_x=concatenate((self.vale_x,other.vale_x))
187 vale_x=clip(vale_x,val_min,val_max)
188 vale_x=sort(list(Set(vale_x)))
189 return vale_x, prol_gauche, prol_droite
191 def cut(self,rinf,rsup,prec,crit='RELATIF') :
192 """Renvoie la fonction self dont on a 'coupé' les extrémités en x=rinf et x=rsup
193 pour la recherche de rinf et rsup dans la liste d'abscisses :
194 prec=precision crit='absolu' ou 'relatif'
196 para=copy.copy(self.para)
197 para['PROL_GAUCHE']='EXCLU'
198 para['PROL_DROITE']='EXCLU'
199 if crit=='ABSOLU' : rinf_tab=greater(abs(self.vale_x-rinf),prec)
200 elif crit=='RELATIF': rinf_tab=greater(abs(self.vale_x-rinf),prec*rinf)
201 else : raise FonctionError, 'fonction : cut : critère absolu ou relatif'
202 if crit=='ABSOLU' : rsup_tab=greater(abs(self.vale_x-rsup),prec)
203 elif crit=='RELATIF': rsup_tab=greater(abs(self.vale_x-rsup),prec*rsup)
204 else : raise FonctionError, 'fonction : cut : critère absolu ou relatif'
205 if alltrue(rinf_tab) : i=searchsorted(self.vale_x,rinf)
206 else : i=rinf_tab.tolist().index(0)+1
207 if alltrue(rsup_tab) : j=searchsorted(self.vale_x,rsup)
208 else : j=rsup_tab.tolist().index(0)
209 vale_x=array([rinf,]+self.vale_x.tolist()[i:j]+[rsup,])
210 vale_y=array([self(rinf),]+self.vale_y.tolist()[i:j]+[self(rsup),])
211 return t_fonction(vale_x,vale_y,para)
213 def cat(self,other,surcharge) :
214 """renvoie une fonction concaténée avec une autre, avec règles de surcharge
216 para=copy.copy(self.para)
217 if self.para['INTERPOL']!=other.para['INTERPOL'] : raise FonctionError, 'concaténation de fonctions à interpolations différentes'
218 if min(self.vale_x)<min(other.vale_x) :
224 para['PROL_GAUCHE']=f1.para['PROL_GAUCHE']
225 if surcharge=='GAUCHE' :
226 i=searchsorted(f2.vale_x,f1.vale_x[-1])
227 if i!=len(f2.vale_x) : para['PROL_DROITE']=f2.para['PROL_DROITE']
228 else : para['PROL_DROITE']=f1.para['PROL_DROITE']
229 vale_x=array(f1.vale_x.tolist()+f2.vale_x[i:].tolist())
230 vale_y=array(f1.vale_y.tolist()+f2.vale_y[i:].tolist())
231 elif surcharge=='DROITE' :
232 i=searchsorted(f1.vale_x,f2.vale_x[0])
233 if i!=len(f1.vale_x) : para['PROL_DROITE']=f2.para['PROL_DROITE']
234 else : para['PROL_DROITE']=f1.para['PROL_DROITE']
235 vale_x=array(f1.vale_x[:i].tolist()+f2.vale_x.tolist())
236 vale_y=array(f1.vale_y[:i].tolist()+f2.vale_y.tolist())
237 return t_fonction(vale_x,vale_y,para)
240 """mise en forme de la fonction selon un vecteur unique (x1,y1,x2,y2,...)
242 __tab=array([self.vale_x,self.vale_y])
243 return ravel(transpose(__tab)).tolist()
246 """renvoie un dictionnaire des valeurs d'ordonnées min et max
248 val_min=min(self.vale_y)
249 val_max=max(self.vale_y)
251 vm['min']=[[self.vale_y[i],self.vale_x[i]] for i in range(len(self.vale_x))\
252 if self.vale_y[i]==val_min]
253 vm['max']=[[self.vale_y[i],self.vale_x[i]] for i in range(len(self.vale_x))\
254 if self.vale_y[i]==val_max]
257 for item in vm['min'] : item.reverse()
258 for item in vm['max'] : item.reverse()
261 def trapeze(self,coef) :
262 """renvoie la primitive de la fonction, calculée avec la constante d'intégration 'coef'
264 trapz = zeros(len(self.vale_y),Float)
266 trapz[1:] = (self.vale_y[1:]+self.vale_y[:-1])/2*(self.vale_x[1:]-self.vale_x[:-1])
268 para=copy.copy(self.para)
269 para['PROL_GAUCHE']='EXCLU'
270 para['PROL_DROITE']='EXCLU'
271 if para['NOM_RESU'][:4]=='VITE' : para['NOM_RESU']='DEPL'
272 elif para['NOM_RESU'][:4]=='ACCE' : para['NOM_RESU']='VITE'
273 return t_fonction(self.vale_x,prim_y,para)
275 def simpson(self,coef) :
276 """renvoie la primitive de la fonction, calculée avec la constante d'intégration 'coef'
278 para=copy.copy(self.para)
279 para['PROL_GAUCHE']='EXCLU'
280 para['PROL_DROITE']='EXCLU'
281 if para['NOM_RESU'][:4]=='VITE' : para['NOM_RESU']='DEPL'
282 elif para['NOM_RESU'][:4]=='ACCE' : para['NOM_RESU']='VITE'
285 h2 = self.vale_x[1] - self.vale_x[0]
286 tabl = [coef,coef +(fb+fm)*h2/2.]
287 prim_y = copy.copy(tabl)
291 for i in range(2,len(self.vale_x)) :
293 h2 = self.vale_x[i] - self.vale_x[i-1]
299 if h1==0. or h2==0. or abs( deltah / h1 ) <= eps :
302 epsi = deltah / (h1*h2)
303 ct = (1.-epsi*h2,2.+epsi*deltah,1.+epsi*h1)
304 tabl[iperm] = tabl[iperm] + (bma/6.)*(ct[0]*fa+ct[1]*fm+ct[2]*fb)
305 prim_y.append(tabl[iperm])
307 return t_fonction(self.vale_x,prim_y,para)
310 """renvoie la dérivée de la fonction
312 pas=self.vale_x[1:]-self.vale_x[:-1]
313 pentes=(self.vale_y[1:]-self.vale_y[:-1])/(self.vale_x[1:]-self.vale_x[:-1])
314 derive=(pentes[1:]*pas[1:]+pentes[:-1]*pas[:-1])/(pas[1:]+pas[:-1])
315 derv_y=[pentes[0]]+derive.tolist()+[pentes[-1]]
316 para=copy.copy(self.para)
317 para['PROL_GAUCHE']='EXCLU'
318 para['PROL_DROITE']='EXCLU'
319 if para['NOM_RESU'][:4]=='DEPL' : para['NOM_RESU']='VITE'
320 elif para['NOM_RESU'][:4]=='VITE' : para['NOM_RESU']='ACCE'
321 return t_fonction(self.vale_x,derv_y,para)
324 """renvoie l'inverse de la fonction
325 on intervertit vale_x et vale_y, on swape interpolation
327 para=copy.copy(self.para)
328 para['NOM_RESU']='TOUTRESU'
329 para['NOM_PARA']=self.para['NOM_PARA']
330 para['INTERPOL'].reverse()
331 if para['PROL_GAUCHE']=='CONSTANT' : para['PROL_GAUCHE']='EXCLU'
332 if para['PROL_DROITE']=='CONSTANT' : para['PROL_DROITE']='EXCLU'
335 if not is_ordo(vale_x) :
338 return t_fonction(vale_x,vale_y,para)
341 """renvoie la mm fonction avec valeur absolue des ordonnées
343 para=copy.copy(self.para)
344 if para['PROL_GAUCHE']=='LINEAIRE' : para['PROL_GAUCHE']='EXCLU'
345 if para['PROL_DROITE']=='LINEAIRE' : para['PROL_DROITE']='EXCLU'
346 return t_fonction(self.vale_x,absolute(self.vale_y),para)
348 def evalfonc(self,liste_val) :
349 """renvoie la mm fonction interpolée aux points définis par la liste 'liste_val'
351 return self.__class__(liste_val,map(self,liste_val),self.para)
353 def suppr_tend(self) :
354 """pour les corrections d'accélérogrammes
355 suppression de la tendance moyenne d'une fonction
357 para=copy.copy(self.para)
358 xy=sum(self.vale_x*self.vale_y)
361 xx=sum(self.vale_x*self.vale_x)
363 a1 = ( n*xy - x0*y0) / (n*xx - x0*x0)
364 a0 = (xx*x0 - x0*xy) / (n*xx - x0*x0)
365 return t_fonction(self.vale_x,self.vale_y-a1*self.vale_x-a0,self.para)
368 """norme de la fonction
371 __ex=__ex.trapeze(0.)
372 return sqrt(__ex.vale_y[-1])
374 def fft(self,methode) :
375 """renvoie la transformée de Fourier rapide FFT
378 para=copy.copy(self.para)
379 para['NOM_PARA']='FREQ'
380 if self.para['NOM_PARA']!='INST' :
381 raise ParametreError, 'fonction réelle : FFT : NOM_PARA=INST pour une transformée directe'
382 pas = self.vale_x[1]-self.vale_x[0]
383 for i in range(1,len(self.vale_x)) :
384 ecart = abs(((self.vale_x[i]-self.vale_x[i-1])-pas)/pas)
386 raise FonctionError, 'fonction réelle : FFT : la fonction doit etre à pas constant'
387 n=int(log(len(self.vale_x))/log(2))
388 if methode=='TRONCATURE' :
389 vale_y=self.vale_y[:2**n]
390 elif methode=='PROL_ZERO' :
391 vale_y=self.vale_y.tolist()
392 vale_y=vale_y+[0.]*(2**(n+1)-len(self.vale_x))
394 elif methode=='COMPLET' :
397 pasfreq =1./(pas*(len(vect)-1))
398 vale_x =[pasfreq*i for i in range(len(vect))]
400 return t_fonction_c(vale_x,vale_y,para)
403 # -----------------------------------------------------------------------------
404 class t_fonction_c(t_fonction) :
405 """Classe pour fonctions complexes, équivalent au type aster = fonction_c
408 """mise en forme de la fonction selon un vecteur unique (x1,yr1,yi1,x2,yr2,yr2,...)
410 __tab=array([self.vale_x,self.vale_y.real,self.vale_y.imag])
411 return ravel(transpose(__tab)).tolist()
414 """affichage de la fonction en double colonne
417 for i in range(len(self.vale_x)) :
418 texte.append('%f %f + %f .j' % (self.vale_x[i],self.vale_y[i].real,self.vale_y[i].imag))
419 return '\n'.join(texte)
421 def fft(self,methode,syme) :
422 """renvoie la transformée de Fourier rapide FFT (sens inverse)
425 para=copy.copy(self.para)
426 para['NOM_PARA']='INST'
427 if self.para['NOM_PARA']!='FREQ' :
428 raise ParametreError, 'fonction complexe : FFT : NOM_PARA=FREQ pour une transformée directe'
429 pas = self.vale_x[1]-self.vale_x[0]
430 for i in range(1,len(self.vale_x)) :
431 ecart = abs(((self.vale_x[i]-self.vale_x[i-1])-pas)/pas)
433 raise FonctionError, 'fonction complexe : FFT : la fonction doit etre à pas constant'
434 n=int(log(len(self.vale_x))/log(2))
435 if syme=='OUI' and len(self.vale_x)==2**n :
436 vale_fonc=self.vale_y
437 elif syme=='NON' and len(self.vale_x)==2**n :
438 vale_fonc=self.vale_y.tolist()
439 vale_fon1=vale_fonc[:]
441 vale_fonc=vale_fonc+vale_fon1
442 vale_fonc=array(vale_fonc)
443 elif syme=='NON' and len(self.vale_x)!=2**n and methode=='PROL_ZERO' :
444 vale_fonc=self.vale_y.tolist()+[complex(0.)]*(2**(n+1)-len(self.vale_x))
445 vale_fon1=vale_fonc[:]
447 vale_fonc=vale_fonc+vale_fon1
448 vale_fonc=array(vale_fonc)
449 elif syme=='NON' and len(self.vale_x)!=2**n and methode=='TRONCATURE' :
450 vale_fonc=self.vale_y[:2**n]
451 vale_fonc=vale_fonc.tolist()
452 vale_fon1=vale_fonc[:]
454 vale_fonc=vale_fonc+vale_fon1
455 vale_fonc=array(vale_fonc)
456 if syme=='OUI' and methode!='COMPLET' and len(self.vale_x)!=2**n :
457 raise FonctionError, 'fonction complexe : FFT : syme=OUI et nombre de points<>2**n'
458 if methode!='COMPLET' :
459 part1=vale_fonc[:len(vale_fonc)/2+1]
460 part2=vale_fonc[1:len(vale_fonc)/2]
461 if methode=='COMPLET' :
462 mil=int(len(self.vale_x)/2)
463 mil2=int((len(self.vale_x)+1)/2)
464 vale_fonc=self.vale_y
466 mil=int(len(self.vale_x)/2)
467 mil2=int((len(self.vale_x)+1)/2)
469 part2=vale_fonc[1:-1]
472 vale_fonc=array(vale_fonc.tolist()+part2)
473 mil=int(len(self.vale_x)/2)*2
474 mil2=int((len(self.vale_x)+1)/2)*2-1
475 part1=vale_fonc[:mil+1]
476 part2=vale_fonc[1:mil2]
477 part2=conjugate(part2)
480 vale_fonc=array(part1.tolist()+part2)
481 vect=FFT.inverse_fft(vale_fonc)
483 pasfreq =1./(pas*(len(vect)-1))
484 vale_x =[pasfreq*i for i in range(len(vect))]
485 pas2 =(1./self.vale_x[-1])*((len(self.vale_x))/float(len(vect)))
487 return t_fonction(vale_x,vale_y,para)
490 # -----------------------------------------------------------------------------
492 """Classe pour nappes, équivalent au type aster = nappe_sdaster
494 def __init__(self,vale_para,l_fonc,para,nom='') :
495 """Création d'un objet nappe
496 - vale_para est la liste de valeur des parametres (mot clé PARA dans DEFI_NAPPE)
497 - para est un dictionnaire contenant les entrées PROL_DROITE, PROL_GAUCHE et INTERPOL (cf sd ASTER)
498 - l_fonc est la liste des fonctions, de cardinal égal à celui de vale_para
503 if pk!=['INTERPOL','NOM_PARA','NOM_PARA_FONC','NOM_RESU','PROL_DROITE','PROL_GAUCHE'] :
504 raise FonctionError, 'nappe : parametres incorrects'
505 if para['INTERPOL'] not in [['NON','NON'],['LIN','LIN'],
506 ['LIN','LOG'],['LOG','LOG'],['LOG','LIN'],] :
507 raise FonctionError, 'nappe : parametre INTERPOL incorrect'
508 if para['PROL_DROITE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
509 raise FonctionError, 'nappe : parametre PROL_DROITE incorrect'
510 if para['PROL_GAUCHE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
511 raise FonctionError, 'nappe : parametre PROL_GAUCHE incorrect'
512 self.vale_para = array(vale_para)
513 if type(l_fonc) not in (types.ListType,types.TupleType) :
514 raise FonctionError, 'nappe : la liste de fonctions fournie n est pas une liste'
515 if len(l_fonc)!=len(vale_para) :
516 raise FonctionError, 'nappe : nombre de fonctions différent du nombre de valeurs du paramètre'
518 if not isinstance(f,t_fonction) and not isinstance(f,t_fonction_c) :
519 raise FonctionError, 'nappe : les fonctions fournies ne sont pas du bon type'
523 def __call__(self,val1,val2):
524 """méthode pour évaluer nappe(val1,val2)
526 i=searchsorted(self.vale_para,val1)
527 n=len(self.vale_para)
529 if val1==self.vale_para[0] : return self.l_fonc[0](val2)
530 if val1 <self.vale_para[0] :
531 if self.para['PROL_GAUCHE']=='EXCLU' : raise ParametreError, 'nappe évaluée hors du domaine de définition'
532 if self.para['PROL_GAUCHE']=='CONSTANT' : return self.l_fonc[0](val2)
533 if self.para['PROL_GAUCHE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val1,
536 self.l_fonc[0](val2),
537 self.l_fonc[1](val2))
539 if val1==self.vale_para[-1] : return self.l_fonc[-1](val2)
540 if val1 >self.vale_para[-1] :
541 if self.para['PROL_DROITE']=='EXCLU' : raise ParametreError, 'nappe évaluée hors du domaine de définition'
542 if self.para['PROL_DROITE']=='CONSTANT' : return self.l_fonc[-1](val2)
543 if self.para['PROL_DROITE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val1,
546 self.l_fonc[-1](val2),
547 self.l_fonc[-2](val2))
549 return interp(self.para['INTERPOL'],val1,self.vale_para[i-1],
551 self.l_fonc[i-1](val2),
552 self.l_fonc[i](val2))
554 def evalfonc(self, liste_val) :
555 """Renvoie la mm nappe dont les fonctions sont interpolées aux points définis
556 par la liste 'liste_val'.
559 for f in self.l_fonc:
560 f2 = f.evalfonc(liste_val)
562 return t_nappe(self.vale_para, l_fonc, self.para)
564 def __add__(self,other) :
565 """addition avec une autre nappe ou un nombre, par surcharge de l'opérateur +
568 if isinstance(other,t_nappe):
569 if self.vale_para!=other.vale_para : raise ParametreError, 'nappes à valeurs de paramètres différentes'
570 for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]+other.l_fonc[i])
571 elif type(other) in [types.FloatType,types.IntType] :
572 for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]+other)
573 else: raise FonctionError, 't_nappe : erreur de type dans __add__'
574 return t_nappe(self.vale_para,l_fonc,self.para)
576 def __mul__(self,other) :
577 """multiplication avec une autre fonction ou un nombre, par surcharge de l'opérateur *
580 if isinstance(other,t_nappe):
581 if self.vale_para!=other.vale_para : raise ParametreError, 'nappes à valeurs de paramètres différentes'
582 for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]*other.l_fonc[i])
583 elif type(other) in [types.FloatType,types.IntType] :
584 for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]*other)
585 else: raise FonctionError, 't_nappe : erreur de type dans __mul__'
586 return t_nappe(self.vale_para,l_fonc,self.para)
589 """affichage de la nappe en double colonne
592 for i in range(len(self.vale_para)) :
593 texte.append('paramètre : %f' % self.vale_para[i])
594 texte.append(repr(self.l_fonc[i]))
595 return '\n'.join(texte)
597 def homo_support(self,other) :
598 """Renvoie la nappe self avec un support union de celui de self et de other
599 le support est la discrétisation vale_para et les discrétisations des fonctions
603 vale_para=self.vale_para.tolist()+other.vale_para.tolist()
604 vale_para=list(Set(vale_para))
606 vale_para=array(vale_para)
608 for val in vale_para :
609 if val in self.vale_para:
610 l_fonc.append(self.l_fonc[searchsorted(self.vale_para, val)])
611 elif val in other.vale_para:
612 other_fonc=other.l_fonc[searchsorted(other.vale_para, val)]
613 new_vale_x=other_fonc.vale_x
614 new_para =other_fonc.para
615 new_vale_y=[self(val,x) for x in new_vale_x]
616 if isinstance(other_fonc, t_fonction):
617 l_fonc.append(t_fonction(new_vale_x, new_vale_y, new_para))
618 if isinstance(other_fonc, t_fonction_c):
619 l_fonc.append(t_fonction_c(new_vale_x, new_vale_y, new_para))
621 raise FonctionError, 'combinaison de nappes : incohérence'
622 return t_nappe(vale_para,l_fonc,self.para)
625 """renvoie un dictionnaire des valeurs d'ordonnées min et max
627 val_min=min([min(fonc.vale_y) for fonc in self.l_fonc])
628 val_max=max([max(fonc.vale_y) for fonc in self.l_fonc])
629 vm={'min':[],'max':[]}
630 for i in range(len(self.vale_para)) :
631 for j in range(len(self.l_fonc[i].vale_y)) :
632 y = self.l_fonc[i].vale_y[j]
633 if y==val_min : vm['min'].append([y,self.l_fonc[i].vale_x[j],self.vale_para[i]])
634 if y==val_max : vm['max'].append([y,self.l_fonc[i].vale_x[j],self.vale_para[i]])
637 for item in vm['min'] : item.reverse()
638 for item in vm['max'] : item.reverse()
642 # -----------------------------------------------------------------------------
643 def homo_support_nappe(l_f):
644 """Homogénéisation du support d'une liste de nappes.
645 Retourne une liste de nappes.
647 if type(l_f) not in (list, tuple):
651 assert isinstance(f, t_nappe), 'Erreur : homo_support_nappe est réservé aux nappes !'
654 __ff=__ff.homo_support(g)
658 # -----------------------------------------------------------------------------
659 def func_union(func,l_f) :
660 """Retourne la fonction x : func(y0=l_f[0](x), y1=l_f[1](x), ...)
661 sur la liste d'abscisses union de celles de self et de other.
663 para = copy.copy(l_f[0].para)
664 # Pour les prolongements et l'interpolation, c'est la première fonction qui prime
667 vale_x = vale_x + f.vale_x.tolist()
668 # on ote les abscisses doublons
669 vale_x = list(Set(vale_x))
671 vale_x = array(vale_x)
672 # interpolation des fonctions sur l'union des abscisses
673 vale_y = [map(f,vale_x) for f in l_f]
674 # applique la fonction
675 vale_y = map(func, *vale_y)
676 return t_fonction(vale_x, vale_y, para)
678 def enveloppe(l_f, crit):
679 """renvoie l'enveloppe supérieure ou inférieure de self et other.
681 if crit.upper() == 'SUP':
682 env = func_union(max, l_f)
683 elif crit.upper() == 'INF':
684 env = func_union(min, l_f)
686 raise FonctionError, 'enveloppe : le critère doit etre SUP ou INF !'
689 def fractile(l_f, fract):
690 """renvoie l'enveloppe supérieure ou inférieure de self et other.
692 para = copy.copy(l_f[0].para)
693 # Pour les prolongements et l'interpolation, c'est la première fonction qui prime
696 vale_x = vale_x + f.vale_x.tolist()
697 # on ote les abscisses doublons
698 vale_x = list(Set(vale_x))
700 vale_x = array(vale_x)
704 vale_y = map(f,vale_x)
705 l_vale_y.append(vale_y)
706 tab_val=transpose(array(l_vale_y))
707 tab_val=tab_val.tolist()
708 for l in tab_val : l.sort()
711 for l_val in tab_val :
712 vale_y.append(l_val[-1])
714 indice=int((len(tab_val[0])-1)*fract)
715 reste =(len(tab_val[0])-1)*fract-indice
716 for l_val in tab_val :
717 vale_y.append(l_val[indice]*(1-reste)+l_val[indice+1]*reste)
718 return t_fonction(vale_x, vale_y, para)