Sprout from BR_V1_11a1 2006-11-27 16:25:02 UTC Pascale Noyret <pascale.noyret@edf.fr> 'PN pour module'
Delete:
Aster/Cata/Utilitai/Graph.py
Aster/Cata/Utilitai/Table.py
Aster/Cata/Utilitai/UniteAster.py
Aster/Cata/Utilitai/Utmess.py
Aster/Cata/Utilitai/__init__.py
Aster/Cata/Utilitai/courbes.py
Aster/Cata/Utilitai/funct_root.py
Aster/Cata/Utilitai/partition.py
Aster/Cata/Utilitai/sup_gmsh.py
Aster/Cata/Utilitai/t_fonction.py
Aster/Cata/Utilitai/transpose.py
+++ /dev/null
-#@ MODIF Graph Utilitai DATE 24/05/2005 AUTEUR MCOURTOI M.COURTOIS
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-# RESPONSABLE MCOURTOI M.COURTOIS
-
-import sys
-import os
-import os.path
-import string
-import re
-import types
-import time
-import Numeric
-
-# try/except pour utiliser hors aster
-try:
- import aster
-except ImportError:
- class fake_aster:
- def repout(self): return '/opt/aster/outils'
- aster=fake_aster()
-
-try:
- from Utilitai.Utmess import UTMESS
-except ImportError:
- def UTMESS(code,sprg,texte):
- fmt='\n <%s> <%s> %s\n\n'
- print fmt % (code,sprg,texte)
-
-if not sys.modules.has_key('Table'):
- try:
- from Utilitai import Table
- except ImportError:
- import Table
-
-
-# ------------------------------------------------------------------------------
-class Graph:
- """Cette classe définit l'objet Graph pour Code_Aster.
-
- Important : Utiliser les méthodes dédiées à la manipulation des données
- (AjoutCourbe, ...) car elles tiennent à jour les attributs "privés"
- relatifs aux données : NbCourbe, les extrema...
-
- Attributs :
- - Données de chaque courbe :
- .Valeurs : liste des valeurs de chaque courbe, pour chaque courbe :
- (paramètres, parties réelles [, parties imaginaires])
- .Legendes : liste des noms de chaque courbe
- .Labels : liste des noms des colonnes de chaque courbe
- .Styles : liste des infices de styles de ligne
- .Couleurs : liste des indices de couleurs
- .Marqueurs : liste des indices de symboles/marqueurs
- .FreqMarq : liste des fréquences des marqueurs
- .Tri : liste du tri à effectuer sur les données ('N', 'X', 'Y',
- 'XY' ou 'YX')
- Pour Lignes, Couleurs, Marqueurs, FreqMarq, -1 signifie valeur par défaut
- du traceur.
-
- - Propriétés :
- .Titre : titre du graphique
- .SousTitre : sous-titre (appelé commentaire dans agraf)
- - Axes :
- .Min_X, .Max_X, .Min_Y, .Max_Y : bornes du tracé (une méthode permet de
- renseigner automatiquement ces valeurs avec les extréma globaux)
- .Legende_X, .Legende_Y : légende des axes
- .Echelle_X, .Echelle_Y : type d'échelle (LIN, LOG)
- .Grille_X, .Grille_Y : paramètre de la grille (pas ou fréquence au choix
- de l'utilisateur en fonction du traceur qu'il veut utiliser)
-
- Attributs privés (modifiés uniquement par les méthodes de la classe) :
- .NbCourbe : nombre de courbes
- .BBXmin, BBXmax, BBYmin, BBYmax : extrema globaux (bounding box)
- .LastTraceArgs, LastTraceFormat : données utilisées lors du dernier tracé
- """
-# ------------------------------------------------------------------------------
- def __init__(self):
- """Construction + valeurs par défaut des attributs
- """
- self.Valeurs = []
- self.Legendes = []
- self.Labels = []
- self.Styles = []
- self.Couleurs = []
- self.Marqueurs = []
- self.FreqMarq = []
- self.Tri = []
- self.Titre = ''
- self.SousTitre = ''
- self.Min_X = 1.e+99
- self.Max_X = -1.e+99
- self.Min_Y = 1.e+99
- self.MinP_X = 1.e+99 # minimum > 0 pour les échelles LOG
- self.MinP_Y = 1.e+99
- self.Max_Y = -1.e+99
- self.Legende_X = ''
- self.Legende_Y = ''
- self.Echelle_X = 'LIN'
- self.Echelle_Y = 'LIN'
- self.Grille_X = -1
- self.Grille_Y = -1
- # attributs que l'utilisateur ne doit pas modifier
- self.NbCourbe = len(self.Valeurs)
- self.BBXmin = self.Min_X
- self.BBXmax = self.Max_X
- self.BBYmin = self.Min_Y
- self.BBYmax = self.Max_Y
- # pour conserver les paramètres du dernier tracé
- self.LastTraceArgs = {}
- self.LastTraceFormat = ''
- return
-# ------------------------------------------------------------------------------
- def SetExtrema(self,marge=0., x0=None, x1=None, y0=None, y1=None):
- """Remplit les limites du tracé (Min/Max_X/Y) avec les valeurs de la
- bounding box +/- avec une 'marge'*(Max-Min)/2.
- x0,x1,y0,y1 permettent de modifier la bb.
- """
- if x0<>None: self.BBXmin=min([self.BBXmin, x0])
- if x1<>None: self.BBXmax=max([self.BBXmax, x1])
- if y0<>None: self.BBYmin=min([self.BBYmin, y0])
- if y1<>None: self.BBYmax=max([self.BBYmax, y1])
-
- dx=max(self.BBXmax-self.BBXmin,0.01*self.BBXmax)
- self.Min_X = self.BBXmin - marge*dx/2.
- self.Max_X = self.BBXmax + marge*dx/2.
- dy=max(self.BBYmax-self.BBYmin,0.01*self.BBYmax)
- self.Min_Y = self.BBYmin - marge*dy/2.
- self.Max_Y = self.BBYmax + marge*dy/2.
- return
-# ------------------------------------------------------------------------------
- def AutoBB(self,debut=-1):
- """Met à jour automatiquement la "bounding box"
- (extrema toutes courbes confondues)
- Appelé par les méthodes de manipulation des données
- """
- if debut == -1:
- debut=self.NbCourbe-1
- if debut == 0:
- X0 = 1.e+99
- X1 = -1.e+99
- Y0 = 1.e+99
- Y1 = -1.e+99
- else:
- X0 = self.BBXmin
- X1 = self.BBXmax
- Y0 = self.BBYmin
- Y1 = self.BBYmax
-
- for i in range(debut,self.NbCourbe):
- X0 = min([X0,]+list(self.Valeurs[i][0]))
- X1 = max([X1,]+list(self.Valeurs[i][0]))
- self.MinP_X = min([self.MinP_X,]+[x for x \
- in list(self.Valeurs[i][0]) if x>0])
- for ny in range(1,len(self.Valeurs[i])):
- Y0 = min([Y0,]+list(self.Valeurs[i][ny]))
- Y1 = max([Y1,]+list(self.Valeurs[i][ny]))
- self.MinP_Y = min([self.MinP_Y,]+[y for y \
- in list(self.Valeurs[i][0]) if y>0])
- self.BBXmin = X0
- self.BBXmax = X1
- self.BBYmin = Y0
- self.BBYmax = Y1
- return
-# ------------------------------------------------------------------------------
- def AjoutCourbe(self,Val,Lab,Leg='',Sty=-1,Coul=-1,Marq=-1,FreqM=-1,Tri='N'):
- """Ajoute une courbe dans les données
- Val : liste de 2 listes (ou 3 si complexe) : abs, ord[, imag]
- Leg : une chaine
- Lab : liste de 2 chaines (ou 3 si complexe)
- Sty : un entier
- Coul : un entier
- Marq : un entier
- FreqM : un entier
- Tri : chaine de caractères : N, X, Y, XY ou YX
- Met à jour les attributs : NbCourbe, BBXmin/Xmax/Ymin/Ymax
- """
- nbc = len(Val) # nombre de colonnes : 2 ou 3
-
- # verifications : "if not (conditions requises)"
- if not ( 2 <= nbc <= 3 and \
- type(Val[0]) in (types.ListType, types.TupleType) and \
- type(Val[1]) in (types.ListType, types.TupleType) and \
- (nbc==2 or type(Val[2]) in (types.ListType, types.TupleType)) and \
- len(Val[0]) == len(Val[1]) and (nbc==2 or len(Val[0]) == len(Val[2])) ):
- UTMESS('S','Graph','"Val" doit etre une liste de 2 ou 3 listes de rééls de meme longueur')
-
- if len(Lab) <> nbc:
- UTMESS('S','Graph','"Lab" doit etre une liste de 2 ou 3 chaines')
-
- # ajout dans les données
- self.Legendes.append(str(Leg))
- self.Labels.append([str(L) for L in Lab])
- self.Valeurs.append(Val)
- self.Styles.append(Sty)
- self.Couleurs.append(Coul)
- self.Marqueurs.append(Marq)
- self.FreqMarq.append(FreqM)
- self.Tri.append(Tri)
-
- self.NbCourbe = self.NbCourbe + 1
- self.AutoBB()
- return
-# ------------------------------------------------------------------------------
- def Courbe(self,n):
- """Permet de récupérer les données de la courbe d'indice n sous forme
- d'un dictionnaire.
- """
- dico={
- 'Leg' : self.Legendes[n], # légende de la courbe
- 'LabAbs' : self.Labels[n][0], # labels des abscisses
- 'LabOrd' : [self.Labels[n][1],], # labels des ordonnées
- 'NbCol' : len(self.Valeurs[n]), # nombre de colonnes
- 'NbPts' : len(self.Valeurs[n][0]), # nombre de points
- 'Abs' : self.Valeurs[n][0], # liste des abscisses
- 'Ord' : [self.Valeurs[n][1],], # liste des ordonnées
- 'Sty' : self.Styles[n], # style de la ligne
- 'Coul' : self.Couleurs[n], # couleur
- 'Marq' : self.Marqueurs[n], # marqueur
- 'FreqM' : self.FreqMarq[n], # fréquence du marqueur
- 'Tri' : self.Tri[n], # ordre de tri des données
- }
- if(dico['NbCol'] == 3):
- dico['LabOrd'].append(self.Labels[n][2]) # labels de la partie imaginaire
- dico['Ord'].append(self.Valeurs[n][2]) # liste des ordonnées partie imaginaire
- return dico
-# ------------------------------------------------------------------------------
- def Trace(self,FICHIER=None,FORMAT=None,dform=None,**opts):
- """Tracé du Graph selon le format spécifié.
- FICHIER : nom du(des) fichier(s). Si None, on dirige vers stdout
- dform : dictionnaire de formats d'impression (format des réels,
- commentaires, saut de ligne...)
- opts : voir TraceGraph.
- """
- para={
- 'TABLEAU' : { 'mode' : 'a', 'driver' : TraceTableau, },
- 'XMGRACE' : { 'mode' : 'a', 'driver' : TraceXmgrace, },
- 'AGRAF' : { 'mode' : 'a', 'driver' : TraceAgraf, },
- }
- kargs={}
- if self.LastTraceArgs=={}:
- kargs['FICHIER']=FICHIER
- kargs['dform']=dform
- kargs['opts']=opts
- else:
- kargs=self.LastTraceArgs.copy()
- if FORMAT==None:
- FORMAT=self.LastTraceFormat
- if FICHIER<>None:
- kargs['FICHIER']=FICHIER
- if dform<>None:
- kargs['dform']=dform
- if opts<>{}:
- kargs['opts']=opts
- if not FORMAT in para.keys():
- print ' <A> <Objet Graph> Format inconnu : %s' % FORMAT
- else:
- kargs['fmod']=para[FORMAT]['mode']
- self.LastTraceArgs = kargs.copy()
- self.LastTraceFormat = FORMAT
- # call the associated driver
- para[FORMAT]['driver'](self,**kargs)
-# ------------------------------------------------------------------------------
- def __repr__(self):
- """Affichage du contenu d'un Graph"""
- srep=''
- for attr in ['NbCourbe','Legendes','Labels','Valeurs','Min_X','Max_X','Min_Y','Max_Y','BBXmax','BBXmin','BBYmax','BBYmin','Legende_X','Legende_Y','Echelle_X','Echelle_Y','Grille_X','Grille_Y','Tri']:
- srep=srep + '%-10s : %s\n' % (attr,str(getattr(self,attr)))
- return srep
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class TraceGraph:
- """
- Cette classe définit le tracé d'un objet Graph dans un fichier.
-
- Attributs :
- .NomFich : liste de noms de fichier de sortie
-
- Attributs privés (modifiés uniquement par les méthodes de la classe) :
- .Fich : liste des objets 'fichier'
- .Graph : objet Graph que l'on veut tracer
- .DicForm : dictionnaire des formats de base (séparateur, format des réels...)
-
- Les méthodes Entete, DescrCourbe, Trace (définition de l'entete, partie descriptive
- d'une courbe, méthode de tracé/impression) sont définies dans une classe dérivée.
- """
-# ------------------------------------------------------------------------------
- def __init__(self,graph,FICHIER,fmod='w',dform=None,opts={}):
- """Construction, ouverture du fichier, surcharge éventuelle du formatage
- (dform), mode d'ouverture du fichier (fmod).
- opts : dictionnaire dont les valeurs seront affectées comme attributs
- de l'objet (A utiliser pour les propriétés spécifiques
- à un format, exemple 'PILOTE' pour Xmgrace).
- """
- # attributs optionnels (au début pour éviter un écrasement maladroit !)
- for k,v in opts.items():
- setattr(self,k,v)
-
- # Ouverture du(des) fichier(s)
- self.NomFich=[]
- if type(FICHIER) is types.StringType:
- self.NomFich.append(FICHIER)
- elif type(FICHIER) in (types.ListType, types.TupleType):
- self.NomFich=FICHIER[:]
- else:
- # dans ce cas, on écrira sur stdout (augmenter le 2 éventuellement)
- self.NomFich=[None]*2
- self.Fich=[]
- for ff in self.NomFich:
- if ff<>None:
- self.Fich.append(open(ff,fmod))
- else:
- self.Fich.append(sys.stdout)
-
- # objet Graph sous-jacent
- self.Graph=graph
- # si Min/Max incohérents
- if graph.Min_X > graph.Max_X or graph.Min_Y > graph.Max_Y:
- graph.SetExtrema(marge=0.05)
- if graph.Min_X < 0. and graph.Echelle_X=='LOG':
- graph.Min_X=graph.MinP_X
- if graph.Min_Y < 0. and graph.Echelle_Y=='LOG':
- graph.Min_Y=graph.MinP_Y
-
- # formats de base (identiques à ceux du module Table)
- self.DicForm={
- 'csep' : ' ', # séparateur
- 'ccom' : '#', # commentaire
- 'cdeb' : '', # début de ligne
- 'cfin' : '\n', # fin de ligne
- 'formK' : '%-12s', # chaines
- 'formR' : '%12.5E', # réels
- 'formI' : '%12d' # entiers
- }
- if dform<>None and type(dform)==types.DictType:
- self.DicForm.update(dform)
-
- # let's go
- self.Trace()
- return
-# ------------------------------------------------------------------------------
- def __del__(self):
- """Fermeture du(des) fichier(s) à la destruction"""
- if hasattr(self,'Fich'):
- self._FermFich()
-# ------------------------------------------------------------------------------
- def _FermFich(self):
- """Fermeture du(des) fichier(s)"""
- for fp in self.Fich:
- if fp<>sys.stdout:
- fp.close()
-# ------------------------------------------------------------------------------
- def _OuvrFich(self):
- """Les fichiers sont ouverts par le constructeur. S'ils ont été fermés,
- par un appel au Tracé, _OuvrFich ouvre de nouveau les fichiers dans le
- meme mode"""
- n=len(self.NomFich)
- for i in range(n):
- if self.Fich[i].closed:
- self.Fich[i]=open(self.NomFich[i],self.Fich[i].mode)
-
-# ------------------------------------------------------------------------------
- def Entete(self):
- """Retourne l'entete"""
- raise StandardError, "Cette méthode doit etre définie par la classe fille."
-# ------------------------------------------------------------------------------
- def DescrCourbe(self,**args):
- """Retourne la chaine de caractères décrivant les paramètres de la courbe.
- """
- raise StandardError, "Cette méthode doit etre définie par la classe fille."
-# ------------------------------------------------------------------------------
- def Trace(self):
- """Méthode pour 'tracer' l'objet Graph dans un fichier.
- Met en page l'entete, la description des courbes et les valeurs selon
- le format et ferme le fichier.
- """
- raise StandardError, "Cette méthode doit etre définie par la classe fille."
-
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class TraceTableau(TraceGraph):
- """
- Impression d'un objet Graph sous forme d'un tableau de colonnes,
- on suppose que les courbes partagent la meme liste d'abscisse à 'EPSILON'
- près, sinon on alarme.
- """
- EPSILON=1.e-4
-# ------------------------------------------------------------------------------
- def Trace(self):
- """Méthode pour 'tracer' l'objet Graph dans un fichier.
- Met en page l'entete, la description des courbes et les valeurs selon
- le format et ferme le fichier.
- L'ouverture et la fermeture du fichier sont gérées par l'objet Table.
- """
- g=self.Graph
- msg=[]
- if g.NbCourbe > 0:
- # validité des données (abscisses identiques)
- t0=Numeric.array(g.Courbe(0)['Abs'])
- max0=max(abs(t0))
- for i in range(1,g.NbCourbe):
- if g.Courbe(i)['NbPts']<>g.Courbe(0)['NbPts']:
- msg.append(" <A> <TraceTableau> La courbe %d n'a pas le meme " \
- "nombre de points que la 1ère." % i)
- else:
- ti=Numeric.array(g.Courbe(i)['Abs'])
- if max(abs((ti-t0).flat)) > self.EPSILON*max0:
- msg.append(" <A> <TraceTableau> Courbe %d : écart entre les "\
- "abscisses supérieur à %9.2E" % (i+1,self.EPSILON))
- msg.append(" Utilisez IMPR_FONCTION pour interpoler " \
- "les valeurs sur la première liste d'abscisses.")
- # objet Table
- Tab=Table.Table()
- # titre / sous-titre
- tit=[]
- tit.append(self.DicForm['ccom']+' '+g.Titre)
- tit.append(self.DicForm['ccom']+' '+g.SousTitre)
- # legendes
- for i in range(g.NbCourbe):
- tit.append(self.DicForm['ccom']+' Courbe '+str(i)+' '+g.Legendes[i])
- Tab.titr=self.DicForm['cfin'].join(tit)
- # noms des paramètres/colonnes
- Tab.para.append(g.Labels[0][0])
- for i in range(g.NbCourbe):
- for lab in g.Labels[i][1:]:
- Tab.para.append(lab)
- # types
- Tab.type=['R']*len(Tab.para)
- # lignes de la Table
- dC0=g.Courbe(0)
- for j in range(dC0['NbPts']):
- row={}
- row[dC0['LabAbs']]=dC0['Abs'][j]
- for i in range(g.NbCourbe):
- dCi=g.Courbe(i)
- for k in range(dCi['NbCol']-1):
- try:
- row[dCi['LabOrd'][k]]=dCi['Ord'][k][j]
- except IndexError:
- row[dCi['LabOrd'][k]]=None
- Tab.append(row)
- Tab.Impr(FICHIER=self.NomFich[0], FORMAT='TABLEAU')
- # erreurs ?
- if msg:
- print '\n'.join(msg)
- return
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class TraceXmgrace(TraceGraph):
- """
- Impression d'un objet Graph au format XMGRACE.
- Attribut supplémentaire : .PILOTE
- """
- PILOTE=''
-# ------------------------------------------------------------------------------
- def Entete(self):
- """Retourne l'entete du fichier .agr correspondant à la mise en forme
- """
- dic_ech={ 'LIN' : 'Normal', 'LOG' : 'Logarithmic' }
- g=self.Graph
- entete=[]
- entete.append("""
-# Grace project file
-#
-@version 50100
-@page size 842, 595
-@page scroll 5%
-@page inout 5%
-@link page off
-@map font 0 to "Times-Roman", "Times-Roman"
-@map font 1 to "Times-Italic", "Times-Italic"
-@map font 2 to "Times-Bold", "Times-Bold"
-@map font 3 to "Times-BoldItalic", "Times-BoldItalic"
-@map font 4 to "Helvetica", "Helvetica"
-@map font 5 to "Helvetica-Oblique", "Helvetica-Oblique"
-@map font 6 to "Helvetica-Bold", "Helvetica-Bold"
-@map font 7 to "Helvetica-BoldOblique", "Helvetica-BoldOblique"
-@map font 8 to "Courier", "Courier"
-@map font 9 to "Courier-Oblique", "Courier-Oblique"
-@map font 10 to "Courier-Bold", "Courier-Bold"
-@map font 11 to "Courier-BoldOblique", "Courier-BoldOblique"
-@map font 12 to "Symbol", "Symbol"
-@map font 13 to "ZapfDingbats", "ZapfDingbats"
-@map color 0 to (255, 255, 255), "white"
-@map color 1 to (0, 0, 0), "black"
-@map color 2 to (255, 0, 0), "red"
-@map color 3 to (0, 255, 0), "green"
-@map color 4 to (0, 0, 255), "blue"
-@map color 5 to (255, 255, 0), "yellow"
-@map color 6 to (188, 143, 143), "brown"
-@map color 7 to (220, 220, 220), "grey"
-@map color 8 to (148, 0, 211), "violet"
-@map color 9 to (0, 255, 255), "cyan"
-@map color 10 to (255, 0, 255), "magenta"
-@map color 11 to (255, 165, 0), "orange"
-@map color 12 to (114, 33, 188), "indigo"
-@map color 13 to (103, 7, 72), "maroon"
-@map color 14 to (64, 224, 208), "turquoise"
-@map color 15 to (0, 139, 0), "green4"
-@reference date 0
-@date wrap off
-@date wrap year 1950
-@timestamp off
-@default linewidth 1.0
-@default linestyle 1
-@default color 1
-@default pattern 1
-@default font 0
-@default char size 1.000000
-@default symbol size 1.000000
-@default sformat "%.8g"
-@background color 0
-@page background fill on
-@r0 off
-@link r0 to g0
-@r0 type above
-@r0 linestyle 1
-@r0 linewidth 1.0
-@r0 color 1
-@r0 line 0, 0, 0, 0
-@r1 off
-@link r1 to g0
-@r1 type above
-@r1 linestyle 1
-@r1 linewidth 1.0
-@r1 color 1
-@r1 line 0, 0, 0, 0
-@r2 off
-@link r2 to g0
-@r2 type above
-@r2 linestyle 1
-@r2 linewidth 1.0
-@r2 color 1
-@r2 line 0, 0, 0, 0
-@r3 off
-@link r3 to g0
-@r3 type above
-@r3 linestyle 1
-@r3 linewidth 1.0
-@r3 color 1
-@r3 line 0, 0, 0, 0
-@r4 off
-@link r4 to g0
-@r4 type above
-@r4 linestyle 1
-@r4 linewidth 1.0
-@r4 color 1
-@r4 line 0, 0, 0, 0
-@g0 on
-@g0 hidden false
-@g0 type XY
-@g0 stacked false
-@g0 bar hgap 0.000000
-@with g0
-@ stack world 0, 0, 0, 0
-@ znorm 1
-@ view xmin 0.150000
-@ view xmax 1.150000
-@ view ymin 0.150000
-@ view ymax 0.850000
-@ title font 0
-@ title size 1.500000
-@ title color 1
-@ subtitle font 0
-@ subtitle size 1.000000
-@ subtitle color 1
-@ xaxes invert off
-@ yaxes invert off
-@ xaxis on
-@ xaxis type zero false
-@ xaxis offset 0.000000 , 0.000000
-@ xaxis bar on
-@ xaxis bar color 1
-@ xaxis bar linestyle 1
-@ xaxis bar linewidth 1.0
-@ xaxis label layout para
-@ xaxis label place auto
-@ xaxis label char size 1.000000
-@ xaxis label font 0
-@ xaxis label color 1
-@ xaxis label place normal
-@ xaxis tick on
-@ xaxis tick minor ticks 1
-@ xaxis tick default 6
-@ xaxis tick place rounded true
-@ xaxis tick in
-@ xaxis tick major size 1.000000
-@ xaxis tick major color 1
-@ xaxis tick major linewidth 1.0
-@ xaxis tick major linestyle 2
-@ xaxis tick major grid on
-@ xaxis tick minor color 1
-@ xaxis tick minor linewidth 1.0
-@ xaxis tick minor linestyle 2
-@ xaxis tick minor grid off
-@ xaxis tick minor size 0.500000
-@ xaxis ticklabel on
-@ xaxis ticklabel format general
-@ xaxis ticklabel prec 5
-@ xaxis ticklabel angle 0
-@ xaxis ticklabel skip 0
-@ xaxis ticklabel stagger 0
-@ xaxis ticklabel place normal
-@ xaxis ticklabel offset auto
-@ xaxis ticklabel offset 0.000000 , 0.010000
-@ xaxis ticklabel start type auto
-@ xaxis ticklabel start 0.000000
-@ xaxis ticklabel stop type auto
-@ xaxis ticklabel stop 0.000000
-@ xaxis ticklabel char size 0.800000
-@ xaxis ticklabel font 0
-@ xaxis ticklabel color 1
-@ xaxis ticklabel formula ""
-@ xaxis ticklabel append ""
-@ xaxis ticklabel prepend ""
-@ xaxis tick place both
-@ xaxis tick spec type none
-@ yaxis on
-@ yaxis type zero false
-@ yaxis offset 0.000000 , 0.000000
-@ yaxis bar on
-@ yaxis bar color 1
-@ yaxis bar linestyle 1
-@ yaxis bar linewidth 1.0
-@ yaxis label layout para
-@ yaxis label place auto
-@ yaxis label char size 1.000000
-@ yaxis label font 0
-@ yaxis label color 1
-@ yaxis label place normal
-@ yaxis tick on
-@ yaxis tick minor ticks 1
-@ yaxis tick default 6
-@ yaxis tick place rounded true
-@ yaxis tick in
-@ yaxis tick major size 1.000000
-@ yaxis tick major color 1
-@ yaxis tick major linewidth 1.0
-@ yaxis tick major linestyle 2
-@ yaxis tick major grid on
-@ yaxis tick minor color 1
-@ yaxis tick minor linewidth 1.0
-@ yaxis tick minor linestyle 1
-@ yaxis tick minor grid off
-@ yaxis tick minor size 0.500000
-@ yaxis ticklabel on
-@ yaxis ticklabel format general
-@ yaxis ticklabel prec 5
-@ yaxis ticklabel angle 0
-@ yaxis ticklabel skip 0
-@ yaxis ticklabel stagger 0
-@ yaxis ticklabel place normal
-@ yaxis ticklabel offset auto
-@ yaxis ticklabel offset 0.000000 , 0.010000
-@ yaxis ticklabel start type auto
-@ yaxis ticklabel start 0.000000
-@ yaxis ticklabel stop type auto
-@ yaxis ticklabel stop 0.000000
-@ yaxis ticklabel char size 0.800000
-@ yaxis ticklabel font 0
-@ yaxis ticklabel color 1
-@ yaxis ticklabel formula ""
-@ yaxis ticklabel append ""
-@ yaxis ticklabel prepend ""
-@ yaxis tick place both
-@ yaxis tick spec type none
-@ altxaxis off
-@ altyaxis off
-@ legend on
-@ legend loctype view
-@ legend 0.85, 0.8
-@ legend box color 1
-@ legend box pattern 1
-@ legend box linewidth 1.0
-@ legend box linestyle 1
-@ legend box fill color 0
-@ legend box fill pattern 1
-@ legend font 0
-@ legend char size 0.750000
-@ legend color 1
-@ legend length 4
-@ legend vgap 1
-@ legend hgap 1
-@ legend invert false
-@ frame type 0
-@ frame linestyle 1
-@ frame linewidth 1.0
-@ frame color 1
-@ frame pattern 1
-@ frame background color 0
-@ frame background pattern 0
-""")
- entete.append('@ title "'+g.Titre+'"')
- entete.append('@ subtitle "'+g.SousTitre+'"')
- entete.append('@ xaxis label "'+g.Legende_X+'"')
- entete.append('@ yaxis label "'+g.Legende_Y+'"')
- entete.append('@ xaxes scale '+dic_ech[g.Echelle_X])
- entete.append('@ yaxes scale '+dic_ech[g.Echelle_Y])
- entete.append('@ xaxis tick major '+str(g.Grille_X))
- entete.append('@ yaxis tick major '+str(g.Grille_Y))
- entete.append('@ world xmin '+str(g.Min_X))
- entete.append('@ world xmax '+str(g.Max_X))
- entete.append('@ world ymin '+str(g.Min_Y))
- entete.append('@ world ymax '+str(g.Max_Y))
- return entete
-# ------------------------------------------------------------------------------
- def DescrCourbe(self,**args):
- """Retourne la chaine de caractères décrivant les paramètres de la courbe.
- """
- # valeurs par défaut
- sty = str(ValCycl(args['Sty'],0,8,1))
- color = str(ValCycl(args['Coul'],1,15,args['NumSet']+1))
- symbol= str(ValCycl(args['Marq'],0,10,args['NumSet']))
- freqm = str(ValCycl(args['FreqM'],0,-1,0))
-
- sn=str(args['NumSet'])
- descr=[]
- descr.append(string.replace("""
-@ s0 hidden false
-@ s0 type xy
-@ s0 symbol size 1.000000
-@ s0 symbol pattern 1
-@ s0 symbol linestyle 1
-@ s0 symbol fill pattern 0
-@ s0 symbol linewidth 1.0
-@ s0 symbol char 65
-@ s0 symbol char font 0
-@ s0 line type 1
-@ s0 line linewidth 1.0
-@ s0 line pattern 1
-@ s0 baseline type 0
-@ s0 baseline off
-@ s0 dropline off
-@ s0 fill type 0
-@ s0 fill rule 0
-@ s0 fill pattern 1
-@ s0 avalue off
-@ s0 avalue type 2
-@ s0 avalue char size 1.000000
-@ s0 avalue font 0
-@ s0 avalue rot 0
-@ s0 avalue format general
-@ s0 avalue prec 3
-@ s0 avalue prepend ""
-@ s0 avalue append ""
-@ s0 avalue offset 0.000000 , 0.000000
-@ s0 errorbar on
-@ s0 errorbar place both
-@ s0 errorbar pattern 1
-@ s0 errorbar size 1.000000
-@ s0 errorbar linewidth 1.0
-@ s0 errorbar linestyle 1
-@ s0 errorbar riser linewidth 1.0
-@ s0 errorbar riser linestyle 1
-@ s0 errorbar riser clip off
-@ s0 errorbar riser clip length 0.100000
-
-@ s0 comment ""
-""",' s0 ',' s'+sn+' '))
- descr.append('@ s'+sn+' symbol '+symbol)
- descr.append('@ s'+sn+' symbol color '+color)
- descr.append('@ s'+sn+' symbol skip '+freqm)
- descr.append('@ s'+sn+' symbol fill color '+color)
- descr.append('@ s'+sn+' line linestyle '+sty)
- descr.append('@ s'+sn+' line color '+color)
- descr.append('@ s'+sn+' fill color '+color)
- descr.append('@ s'+sn+' avalue color '+color)
- descr.append('@ s'+sn+' errorbar color '+color)
- descr.append('@ s'+sn+' legend "'+args['Leg']+'"')
- return descr
-# ------------------------------------------------------------------------------
- def Trace(self):
- """Méthode pour 'tracer' l'objet Graph dans un fichier.
- Met en page l'entete, la description des courbes et les valeurs selon
- le format et ferme le fichier.
- """
- g=self.Graph
- if self.PILOTE=='INTERACTIF':
- self.NomFich[0]='Trace_'+time.strftime('%y%m%d%H%M%S',time.localtime())+'.dat'
- self.Fich[0]=open(self.NomFich[0],'w')
- # initialise le graph
- self._FermFich()
- nbsets, x0, x1, y0, y1 = IniGrace(self.NomFich[0])
- NumSetIni = nbsets+1
- g.SetExtrema(0.05, x0, x1, y0, y1)
- # si Min/Max incohérents
- if g.Min_X < 0. and g.Echelle_X=='LOG':
- g.Min_X=g.MinP_X
- if g.Min_Y < 0. and g.Echelle_Y=='LOG':
- g.Min_Y=g.MinP_Y
-
- self._OuvrFich()
- fich=self.Fich[0]
- if g.NbCourbe < 1:
- self._FermFich()
- return
- # cohérence des valeurs par défaut
- if g.Grille_X<0 or g.Grille_Y<0:
- deltaX=g.Max_X-g.Min_X
- deltaY=g.Max_Y-g.Min_Y
- g.Grille_X=deltaX/5.
- g.Grille_Y=deltaY/5.
- if deltaX>4:
- g.Grille_X=int(round(g.Grille_X))
- if deltaY>4:
- g.Grille_Y=int(round(g.Grille_Y))
- # entete
- fich.write('\n'.join(self.Entete()))
- fich.write('\n')
- # valeurs
- it=-1
- for i in range(g.NbCourbe):
- dCi=g.Courbe(i)
- for k in range(dCi['NbCol']-1):
- it=it+1
- dCi['NumSet'] = NumSetIni + it
- fich.write('\n'.join(self.DescrCourbe(**dCi)))
- fich.write('\n')
- # partie données (.dat)
- lig=[]
- it=-1
- for i in range(g.NbCourbe):
- dCi=g.Courbe(i)
- for k in range(dCi['NbCol']-1):
- it=it+1
- lig.append('@target g0.s%d' % (NumSetIni + it))
- lig.append('@type xy')
- listX, listY = Tri(g.Tri, lx=dCi['Abs'], ly=dCi['Ord'][k])
- for j in range(dCi['NbPts']):
- svX=self.DicForm['formR'] % listX[j]
- svY=self.DicForm['formR'] % listY[j]
- lig.append(self.DicForm['formR'] % listX[j] + \
- ' ' + self.DicForm['formR'] % listY[j])
- lig.append('&')
- fich.write('\n'.join(lig))
- fich.write('\n')
- self._FermFich()
-
- # Production du fichier postscript, jpeg ou lancement interactif
- pilo=self.PILOTE
- if self.PILOTE<>'':
- xmgr=os.path.join(aster.repout(),'xmgrace')
- nfhard=self.NomFich[0]+'.hardcopy'
- # nom exact du pilote
- if pilo=='POSTSCRIPT':
- pilo='PostScript'
- elif pilo=='INTERACTIF':
- pilo='X11'
- # ligne de commande
- if pilo=='X11':
- lcmde=xmgr+' '+self.NomFich[0]
- if not os.environ.has_key('DISPLAY') or os.environ['DISPLAY']=='':
- os.environ['DISPLAY']=':0.0'
- UTMESS('A','TraceXmgrace','Variable DISPLAY non définie')
- UTMESS('I','TraceXmgrace','on fixe le DISPLAY à %s' % os.environ['DISPLAY'])
- else:
- if os.path.exists(os.path.join(aster.repout(),'gracebat')):
- xmgr=os.path.join(aster.repout(),'gracebat')
- lcmde=xmgr+' -hdevice '+pilo+' -hardcopy -printfile '+nfhard+' '+self.NomFich[0]
- # appel xmgrace
- UTMESS('I','TraceXmgrace','Lancement de : '+lcmde)
- if not os.path.exists(xmgr):
- UTMESS('S','TraceXmgrace','Fichier inexistant : '+xmgr)
- iret=os.system(lcmde)
- if iret==0 or os.path.exists(nfhard):
- if pilo not in ['','X11']:
- os.remove(self.NomFich[0]) # necessaire sous windows
- os.rename(nfhard,self.NomFich[0])
- else:
- UTMESS('A','TraceXmgrace',"Erreur lors de l'utilisation du filtre "+pilo+"\nLe fichier retourné est le fichier '.agr'")
- # menage
- if self.PILOTE=='INTERACTIF':
- os.remove(self.NomFich[0])
- return
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class TraceAgraf(TraceGraph):
- """
- Impression d'un objet Graph au format AGRAF.
- """
-# ------------------------------------------------------------------------------
- def Entete(self):
- """Retourne l'entete des directives Agraf"""
- dic_ech={ 'LIN' : '0', 'LOG' : '1' }
- g=self.Graph
- entete=[]
- entete.append("""
-ASPECT_GRAPHIQUE:
- En-tete :Departement Analyses Mecaniques et Acoustique
- Aspect :0
- Nombre de vues :1
- Cesure commentaire :40
- MinMax :0
- Fonte Titre :%helvetica-14
- Fonte Axes :%courier-12
- Fonte Autre :%times-12
-
- DEFAUT_COURBE:
- Couleur (rvb) : 0 0 0
-
- DEFAUT_COURBE:
- Couleur (rvb) : 65535 0 0
-
- DEFAUT_COURBE:
- Couleur (rvb) : 11822 35723 22359
-
- DEFAUT_COURBE:
- Couleur (rvb) : 0 0 65535
-
- DEFAUT_COURBE:
- Couleur (rvb) : 65535 0 65535
-
- DEFAUT_COURBE:
- Couleur (rvb) : 0 65535 65535
-
- DEFAUT_COURBE:
- Couleur (rvb) : 0 65535 0
-
- DEFAUT_COURBE:
- Couleur (rvb) : 41120 21074 11565
-
- DEFAUT_COURBE:
- Couleur (rvb) : 65535 42405 0
-
- DEFAUT_COURBE:
- Couleur (rvb) : 41120 8224 61680
-
- DEFAUT_COURBE:
- Couleur (rvb) : 65535 65535 0
-
- DEFAUT_COURBE:
- Couleur (rvb) : 53970 46260 35980
-
-GRAPHIQUE:
-""")
- if g.Titre=='':
- g.Titre='GRAPHIQUE CODE_ASTER'
- entete.append('Titre :'+g.Titre+'\n')
- if g.SousTitre<>'':
- entete.append('Commentaire :'+g.SousTitre+'\n')
- entete.append('Frequence Grille X :'+str(int(g.Grille_X))+'\n')
- entete.append('Frequence Grille Y :'+str(int(g.Grille_Y))+'\n')
- entete.append('Echelle X :'+dic_ech[g.Echelle_X]+'\n')
- entete.append('Echelle Y :'+dic_ech[g.Echelle_Y]+'\n')
- if g.Legende_X<>'':
- entete.append('Legende X :'+g.Legende_X+'\n')
- if g.Legende_Y<>'':
- entete.append('Legende Y :'+g.Legende_Y+'\n')
- entete.append('Min X : '+str(g.Min_X)+'\n')
- entete.append('Max X : '+str(g.Max_X)+'\n')
- entete.append('Min Y : '+str(g.Min_Y)+'\n')
- entete.append('Max Y : '+str(g.Max_Y)+'\n')
-
- return entete
-# ------------------------------------------------------------------------------
- def DescrCourbe(self,**args):
- """Retourne la chaine de caractères décrivant les paramètres de la courbe.
- """
- # valeurs par défaut
- sty = str(ValCycl(args['Sty'],0,2,0))
- color = str(ValCycl(args['Coul'],0,12,args['NumSet']))
- symbol= str(ValCycl(args['Marq'],0,12,args['NumSet']))
- freqm = str(ValCycl(args['FreqM'],0,-1,0))
-
- descr=[]
- descr.append(' COURBE:\n')
- descr.append(' Trait :'+sty+'\n')
- descr.append(' Couleur :'+color+'\n')
- descr.append(' Marqueur :'+symbol+'\n')
- descr.append(' Frequence Marqueur :'+freqm+'\n')
- if args['Leg']<>'':
- descr.append(' Legende :'+args['Leg']+'\n')
- descr.append(' Tri :'+args['Tri']+'\n')
- descr.append(' Abscisses : [ '+str(args['Bloc'])+', '+str(args['ColX'])+']\n')
- descr.append(' Ordonnees : [ '+str(args['Bloc'])+', '+str(args['ColY'])+']\n')
- return descr
-# ------------------------------------------------------------------------------
- def Trace(self):
- """Méthode pour 'tracer' l'objet Graph dans un fichier.
- Met en page l'entete, la description des courbes et les valeurs selon
- le format et ferme le fichier.
- """
- self._OuvrFich()
- fdogr=self.Fich[0]
- fdigr=self.Fich[1]
- g=self.Graph
- if g.NbCourbe > 0:
- # cohérence des valeurs par défaut
- if g.Grille_X<0 or g.Grille_Y<0:
- g.Grille_X=0
- g.Grille_Y=0
- # entete
- for lig in self.Entete():
- fdigr.write(lig)
- # valeurs
- for i in range(g.NbCourbe):
- dCi=g.Courbe(i)
- dCi['NumSet']=i
- # partie directives (.digr)
- for k in range(dCi['NbCol']-1):
- dCi['Bloc']=i+1
- dCi['ColX']=1
- dCi['ColY']=k+2
- for lig in self.DescrCourbe(**dCi):
- fdigr.write(lig)
- # partie données (.dogr)
- if dCi['Leg']<>'':
- leg=dCi['Leg']
- else:
- leg='COURBE_'+str(i)
- fdogr.write('#NOM DE LA FONCTION: '+leg+'\n')
- for j in range(dCi['NbPts']):
- for k in range(dCi['NbCol']):
- sv=self.DicForm['formR'] % g.Valeurs[i][k][j]
- fdogr.write(' '+sv)
- fdogr.write('\n')
- fdogr.write('\n')
- self._FermFich()
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-def ValCycl(val,vmin,vmax,vdef):
- """
- Retourne une valeur entre vmin et vmax (bornes incluses) :
- - si val<vmin, on utilise val=vdef,
- - si val>vmax, on cycle tel que val=vmax+1 retourne vmin, etc.
- - si vmax<vmin, il n'y a pas de max
- """
- v = val
- if v < vmin:
- v = vdef
- if vmax < vmin:
- return v
- else:
- return (((v-vmin) % (vmax+1-vmin))+vmin)
-
-# ------------------------------------------------------------------------------
-def Tri(tri, lx, ly):
- """Retourne les listes triées selon la valeur de tri ('X','Y','XY','YX').
- """
- dNumCol={ 'X' : 0, 'Y' : 1 }
- tab=Numeric.array((lx,ly),Numeric.Float64)
- tab=Numeric.transpose(tab)
- li=range(len(tri))
- li.reverse()
- for i in li:
- if tri[-i] in dNumCol.keys():
- icol=dNumCol[tri[-i]]
- tab = Numeric.take(tab, Numeric.argsort(tab[:,icol]))
- return [ tab[:,0].tolist(), tab[:,1].tolist() ]
-
-# ------------------------------------------------------------------------------
-def AjoutParaCourbe(dCourbe, args):
- """Ajoute les arguments fournis par l'utilisateur (args) dans le dictionnaire
- décrivant la courbe (dCourbe).
- """
- # correspondance : mot-clé Aster / clé du dico de l'objet Graph
- keys={
- 'LEGENDE' : 'Leg',
- 'STYLE' : 'Sty',
- 'COULEUR' : 'Coul',
- 'MARQUEUR' : 'Marq',
- 'FREQ_MARQUEUR' : 'FreqM',
- 'TRI' : 'Tri',
- }
- for mc, key in keys.items():
- if args.has_key(mc):
- dCourbe[key]=args[mc]
-
-# ------------------------------------------------------------------------------
-def IniGrace(fich):
- """Retourne le numéro de la dernière courbe d'un fichier xmgrace (sinon 0).
- """
- ns=0
- x0=None
- x1=None
- y0=None
- y1=None
- if os.path.exists(fich) and os.stat(fich).st_size<>0:
- os.rename(fich, fich+'.prev')
- fpre=open(fich+'.prev', 'r')
- fnew=open(fich, 'w')
- for line in fpre:
- ikeep=True
- mat=re.search('@target g[0-9]+\.s([0-9]+)', line)
- if mat<>None and int(mat.group(1))>ns:
- ns=int(mat.group(1))
- mat=re.search('@[ ]+world[ ]+xmin[ ]+([\-\+\.0-9eEdD]+)', line)
- if mat<>None:
- try:
- x0=float(mat.group(1))
- ikeep=False
- except ValueError:
- pass
- mat=re.search('@[ ]+world[ ]+xmax[ ]+([\-\+\.0-9eEdD]+)', line)
- if mat<>None:
- try:
- x1=float(mat.group(1))
- ikeep=False
- except ValueError:
- pass
- mat=re.search('@[ ]+world[ ]+ymin[ ]+([\-\+\.0-9eEdD]+)', line)
- if mat<>None:
- try:
- y0=float(mat.group(1))
- ikeep=False
- except ValueError:
- pass
- mat=re.search('@[ ]+world[ ]+ymax[ ]+([\-\+\.0-9eEdD]+)', line)
- if mat<>None:
- try:
- y1=float(mat.group(1))
- ikeep=False
- except ValueError:
- pass
- if ikeep:
- fnew.write(line)
- fpre.close()
- fnew.close()
- print """
- <I> Informations sur le fichier '%s' :
- Nombre de courbes : %3d
- Bornes des abscisses : [ %13.6G , %13.6G ]
- Bornes des ordonnées : [ %13.6G , %13.6G ]
-""" % (fich, ns, x0, x1, y0, y1)
- return ns, x0, x1, y0, y1
+++ /dev/null
-#@ MODIF Table Utilitai DATE 17/05/2005 AUTEUR DURAND C.DURAND
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-# RESPONSABLE MCOURTOI M.COURTOIS
-
-import sys
-import string
-import re
-
-from types import *
-EnumTypes=(ListType, TupleType)
-NumberTypes=(IntType, LongType, FloatType, ComplexType)
-
-# try/except pour utiliser hors aster
-try:
- from Utilitai.Utmess import UTMESS
-except ImportError:
- def UTMESS(code,sprg,texte):
- fmt='\n <%s> <%s> %s\n\n'
- print fmt % (code,sprg,texte)
-
-if not sys.modules.has_key('Graph'):
- try:
- from Utilitai import Graph
- except ImportError:
- import Graph
-
-# formats de base (identiques à ceux du module Graph)
-DicForm={
- 'csep' : ' ', # séparateur
- 'ccom' : '#', # commentaire
- 'cdeb' : '', # début de ligne
- 'cfin' : '\n', # fin de ligne
- 'formK' : '%-8s', # chaines
- 'formR' : '%12.5E', # réels
- 'formI' : '%8d' # entiers
-}
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class TableBase(object):
- """Classe pour partager les méthodes d'impression entre Table et Colonne
- (c'est surtout utile pour vérifier que l'extraction et les filtres sur les
- colonnes sont corrects).
- """
- def __repr__(self):
- return self.ReprTable()
- def Croise(self,**kargs):
- raise StandardError, 'Must be defined in a derived class'
-
-# ------------------------------------------------------------------------------
- def Impr(self,FICHIER=None,FORMAT='TABLEAU',dform=None,**opts):
- """Impresssion de la Table selon le format spécifié.
- FICHIER : nom du(des) fichier(s). Si None, on dirige vers stdout
- dform : dictionnaire de formats d'impression (format des réels,
- commentaires, saut de ligne...)
- opts : selon FORMAT.
- """
- para={
- 'TABLEAU' : { 'mode' : 'a', 'driver' : self.ImprTableau, },
- 'ASTER' : { 'mode' : 'a', 'driver' : self.ImprTableau, },
- 'XMGRACE' : { 'mode' : 'a', 'driver' : self.ImprGraph, },
- 'AGRAF' : { 'mode' : 'a', 'driver' : self.ImprTableau, },
- 'TABLEAU_CROISE' : { 'mode' : 'a', 'driver' : self.ImprTabCroise, },
- }
- kargs={
- 'FICHIER' : FICHIER,
- 'FORMAT' : FORMAT,
- 'dform' : DicForm.copy(),
- 'mode' : para[FORMAT]['mode'],
- }
- if dform<>None and type(dform)==DictType:
- kargs['dform'].update(dform)
- # ajout des options
- kargs.update(opts)
-
- if not kargs.get('PAGINATION'):
- # call the associated driver
- para[FORMAT]['driver'](**kargs)
-
- else:
- if not type(kargs['PAGINATION']) in EnumTypes:
- ppag=[kargs['PAGINATION'],]
- else:
- ppag=list(kargs['PAGINATION'])
- del kargs['PAGINATION']
- npag=len(ppag)
- # paramètres hors ceux de la pagination
- lkeep=[p for p in self.para if ppag.count(p)==0]
- # création des listes des valeurs distinctes
- lvd=[]
- for p in ppag:
- lvp=getattr(self,p).values()
- lvn=[]
- for it in lvp:
- if it<>None and lvn.count(it)==0:
- lvn.append(it)
- lvn.sort()
- lvd.append(lvn)
- # création des n-uplets
- s = '[['+','.join(['x'+str(i) for i in range(npag)])+'] '
- s+= ' '.join(['for x'+str(i)+' in lvd['+str(i)+']' for i in range(npag)])+']'
- try:
- lnup=eval(s)
- except SyntaxError, s:
- UTMESS('F','Table','Erreur lors de la construction des n-uplets')
- # pour chaque n-uplet, on imprime la sous-table
- for nup in lnup:
- tab=self
- for i in range(npag):
- tab = tab & (getattr(tab,ppag[i]) == nup[i])
- sl=''
- if tab.titr: sl='\n'
- tab.titr += sl+ppag[i]+': '+str(nup[i])
- tab[lkeep].Impr(**kargs)
-
-# ------------------------------------------------------------------------------
- def ImprTableau(self,**kargs):
- """Impression au format TABLEAU ou ASTER
- """
- # fichier ou stdout
- if kargs.get('FICHIER')<>None:
- f=open(kargs['FICHIER'],kargs['mode'])
- else:
- f=sys.stdout
- # ecriture
- f.write(self.ReprTable(**kargs) + '\n')
- # fermeture
- if kargs.get('FICHIER')<>None:
- f.close()
-
- def ReprTable(self,FORMAT='TABLEAU',dform=DicForm,**ignore):
- """Représentation d'une Table ou d'une Colonne sous forme d'un tableau.
- """
- rows=self.rows
- para=self.para
- typ =self.type
- if not type(para) in EnumTypes:
- para=[self.para,]
- typ =[self.type,]
- # est-ce que l'attribut .type est renseigné ?
- typdef=typ<>[None]*len(typ)
- txt=[]
- # ['']+ pour ajouter un séparateur en début de ligne
- lspa=['',]
- # lmax : largeur max des colonnes = max(form{K,R,I},len(parametre))
- lmax=[]
- for p in para:
- t=typ[para.index(p)]
- larg_max=max([len(str(p))] + \
- [len(FMT(dform,k,t) % 0) for k in ('formK','formR','formI')])
- lspa.append(FMT(dform,'formK',t,larg_max,str(p)) % p)
- lmax.append(larg_max)
- if typdef:
- stype=dform['csep'].join([''] + \
- [FMT(dform,'formK',typ[i],lmax[i]) % typ[i] for i in range(len(para))])
- txt.append('')
- txt.append('-'*80)
- txt.append('')
- ASTER=(FORMAT=='ASTER')
- if ASTER:
- txt.append('#DEBUT_TABLE')
- if self.titr:
- if ASTER:
- txt.extend(['#TITRE '+lig for lig in self.titr.split('\n')])
- else:
- txt.append(self.titr)
- txt.append(dform['csep'].join(lspa))
- if ASTER and typdef:
- txt.append(stype)
- for r in rows:
- lig=['']
- empty=True
- for v in para:
- i=para.index(v)
- t=typ[i]
- rep=r.get(v,None)
- if type(rep) is FloatType:
- lig.append(FMT(dform,'formR',t,lmax[i]) % rep)
- empty=False
- elif type(rep) is IntType:
- lig.append(FMT(dform,'formI',t,lmax[i]) % rep)
- empty=False
- else:
- if rep==None:
- rep='-'
- else:
- empty=False
- s=FMT(dform,'formK',t,lmax[i],rep) % str(rep)
- # format AGRAF = TABLEAU + '\' devant les chaines de caractères !
- if FORMAT=='AGRAF':
- s='\\'+s
- lig.append(s)
- if not empty:
- txt.append(dform['csep'].join(lig))
- if ASTER:
- txt.append('#FIN_TABLE')
- return dform['cfin'].join(txt)
-# ------------------------------------------------------------------------------
- def ImprTabCroise(self,**kargs):
- """Impression au format TABLEAU_CROISE d'une table ayant 3 paramètres.
- """
- # création du tableau croisé et impression au format TABLEAU
- tabc=self.Croise()
- kargs['FORMAT']='TABLEAU'
- tabc.Impr(**kargs)
-# ------------------------------------------------------------------------------
- def ImprGraph(self,**kargs):
- """Impression au format XMGRACE : via le module Graph
- """
- args=kargs.copy()
- if len(self.para)<>2:
- UTMESS('A','Table','La table doit avoir exactement deux paramètres.')
- return
- lx, ly = [[v for v in getattr(self,p).values() if v<>None] for p in self.para]
- # objet Graph
- graph=Graph.Graph()
- dicC={
- 'Val' : [lx, ly],
- 'Lab' : self.para,
- }
- if args['LEGENDE']==None: del args['LEGENDE']
- Graph.AjoutParaCourbe(dicC, args)
- graph.AjoutCourbe(**dicC)
-
- # Surcharge des propriétés du graphique et des axes
- # (bloc quasiment identique dans impr_fonction_ops)
- if args.get('TITRE'): graph.Titre=args['TITRE']
- if args.get('BORNE_X'):
- graph.Min_X=args['BORNE_X'][0]
- graph.Max_X=args['BORNE_X'][1]
- if args.get('BORNE_Y'):
- graph.Min_Y=args['BORNE_Y'][0]
- graph.Max_Y=args['BORNE_Y'][1]
- if args.get('LEGENDE_X'): graph.Legende_X=args['LEGENDE_X']
- if args.get('LEGENDE_Y'): graph.Legende_Y=args['LEGENDE_Y']
- if args.get('ECHELLE_X'): graph.Echelle_X=args['ECHELLE_X']
- if args.get('ECHELLE_Y'): graph.Echelle_Y=args['ECHELLE_Y']
- if args.get('GRILLE_X'): graph.Grille_X=args['GRILLE_X']
- if args.get('GRILLE_Y'): graph.Grille_Y=args['GRILLE_Y']
-
- try:
- graph.Trace(**args)
- except TypeError:
- UTMESS('A','Table','Les cellules ne doivent contenir que des nombres réels')
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class Table(TableBase):
- """Une table est construite comme une liste de lignes, chaque ligne est
- un dictionnaire.
- On crée puis on ajoute les lignes avec la méthode append :
- t=Table()
- t.append(dict(a=1,b=2))
- t.append(dict(a=3,b=4))
- La méthode __iter__ définit un itérateur sur les lignes de la table,
- __repr__ retourne une représentation de la table, utilisée par "print t".
- Grace à la classe Colonne et à sa méthode _extract, il est possible
- de construire une sous-table qui satisfait un critère donné.
- Le critère est donné par une fonction Python qui retourne vrai
- ou faux si la valeur d'une colonne respecte le critère ou non.
- Exemple:
- def critere(valeur):
- return valeur < 10
- soustable = t.a._extract(critere)
- t.a retourne un objet intermédiaire de la classe Colonne qui mémorise
- le nom de la colonne demandée (a, ici).
- """
- def __init__(self, rows=[], para=[], typ=[], titr=''):
- """Constructeur de la Table :
- rows : liste des lignes (dict)
- para : liste des paramètres
- type : liste des types des paramètres
- titr : titre de la table
- """
- self.rows=[r for r in rows if r.values()<>[None]*len(r.values())]
- self.para=list(para)
- if len(typ)==len(self.para):
- self.type=list(typ)
- else:
- self.type=[None]*len(self.para)
- self.titr=titr
-
- def append(self, obj):
- """Ajoute une ligne (type dict) à la Table"""
- self.rows.append(obj)
-
- def __iter__(self):
- """Itère sur les lignes de la Table"""
- return iter(self.rows)
-
- def __getattr__(self, column):
- """Construit un objet intermediaire (couple table, colonne)"""
- typ=None
- if not column in self.para:
- column=''
- else:
- typ=self.type[self.para.index(column)]
- return Colonne(self, column, typ)
-
- def sort(self, CLES=None, ORDRE='CROISSANT'):
- """Tri de la table.
- CLES : liste des clés de tri
- ORDRE : CROISSANT ou DECROISSANT (de longueur 1 ou len(keys))
- """
- # par défaut, on prend tous les paramètres
- if CLES==None:
- CLES=self.para[:]
- if not type(CLES) in EnumTypes:
- CLES=[CLES,]
- else:
- CLES=list(CLES)
- self.rows=sort_table(self.rows, self.para, CLES, (ORDRE=='DECROISSANT'))
-# if not type(order) in EnumTypes:
-# order=[order,]
-# print 'TRI clés=%s, order=%s' % (keys,order)
-# # on ne garde que le premier si les longueurs sont différentes
-# if len(order)<>len(keys):
-# order=[order[0],]
-# else:
-# # si toutes les valeurs sont identiques, on peut ne garder que la 1ère
-# d={}
-# for o in order: d[o]=None
-# if len(order)<>len(keys) or len(d.keys())==1:
-# order=[order[0],]
-# if len(order)==1:
-# self.rows=sort_table(self.rows, self.para, keys, (order[0]=='DECROISSANT'))
-# else:
-# # de la dernière clé à la première
-# for k,o in [(keys[i],order[i]) for i in range(len(keys)-1,-1,-1)]:
-# print 'TRI : clé=%s, order=%s' % (k,o)
-# self.rows=sort_table(self.rows, self.para, [k], (o=='DECROISSANT'))
-
- def __delitem__(self, args):
- """Supprime les colonnes correspondantes aux éléments de args """
- if not type(args) in EnumTypes:
- args=[args,]
- new_rows=self.rows
- new_para=self.para
- new_type=self.type
- for item in args:
- del new_type[new_para.index(item)]
- new_para.remove(item)
- for line in new_rows : del line[item]
- return Table(new_rows, new_para, new_type, self.titr)
-
- def __getitem__(self, args):
- """Extrait la sous table composée des colonnes dont les paramètres sont dans args """
- if not type(args) in EnumTypes:
- args=[args,]
- else:
- args=list(args)
- #print '<getitem> args=',args
- new_rows=[]
- new_para=args
- new_type=[]
- for item in new_para:
- if not item in self.para:
- return Table()
- new_type.append(self.type[self.para.index(item)])
- for line in self:
- new_line={}
- for item in new_para:
- new_line[item]=line.get(item)
- new_rows.append(new_line)
- return Table(new_rows, new_para, new_type, self.titr)
-
- def __and__(self, other):
- """Intersection de deux tables (opérateur &)"""
- if other.para<>self.para:
- UTMESS('A','Table','Les paramètres sont différents')
- return Table()
- else:
- tmp = [ r for r in self if r in other.rows ]
- return Table(tmp, self.para, self.type, self.titr)
-
- def __or__(self, other):
- """Union de deux tables (opérateur |)"""
- if other.para<>self.para:
- UTMESS('A','Table','Les paramètres sont différents')
- return Table()
- else:
- tmp = self.rows[:]
- tmp.extend([ r for r in other if r not in self ])
- return Table(tmp, self.para, self.type[:], self.titr)
-
- def values(self):
- """Renvoie la table sous la forme d'un dictionnaire de listes dont les
- clés sont les paramètres.
- """
- dico={}
- for column in self.para:
- dico[column]=Colonne(self, column).values()
- return dico
-
- def Array(self,Para,Champ):
- """Renvoie sous forme de NumArray le résultat d'une extraction dans une table
- méthode utile à macr_recal
- """
- import Numeric
- __Rep = self[Para,Champ].values()
- F=Numeric.zeros((len(__Rep[Para]),2),Numeric.Float)
- for i in range(len(__Rep[Para])):
- F[i][0] = __Rep[Para][i]
- F[i][1] = __Rep[Champ][i]
- del(__Rep)
- return F
-
- def Croise(self):
- """Retourne un tableau croisé P3(P1,P2) à partir d'une table ayant
- trois paramètres (P1, P2, P3).
- """
- if len(self.para)<>3:
- UTMESS('A','Table','La table doit avoir exactement trois paramètres.')
- return Table()
- py, px, pz = self.para
- ly, lx, lz = [getattr(self,p).values() for p in self.para]
- new_rows=[]
- #lpz='%s=f(%s,%s)' % (pz,px,py)
- lpz='%s/%s' % (px,py)
- new_para=[lpz,]
- # attention aux doublons dans lx et ly
- for it in ly:
- if it<>None and new_para.count(it)==0:
- new_para.append(it)
- newx=[]
- for it in lx:
- if it<>None and newx.count(it)==0:
- newx.append(it)
- for x in newx:
- if x<>None:
- d={ lpz : x, }
- taux = (getattr(self,px)==x)
- for dz in taux.rows:
- d[dz[py]]=dz[pz]
- new_rows.append(d)
- new_type=[self.type[0],] + [self.type[2]]*len(ly)
- new_titr=self.titr
- if new_titr<>'': new_titr+='\n'
- new_titr+=pz + ' FONCTION DE ' + px + ' ET ' + py
- return Table(new_rows, new_para, new_type, new_titr)
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-class Colonne(TableBase):
- """Classe intermédiaire pour mémoriser un couple (table, nom de colonne)
- et exprimer les critères d'extraction sous une forme naturelle en python
- en surchargeant les operateurs <, >, <> et =.
- Alors on peut écrire la requete simple :
- soustable=t.a<10
- Ainsi que des requetes plus complexes :
- soustable=t.a<10 & t.b <4
- ou
- soustable=t.a<10 | t.b <4
- Les "alias" EQ, NE, LE, LT, GE, GT permettent à la macro IMPR_TABLE
- d'utiliser directement le mot-clé utilisateur CRIT_COMP défini dans le
- catalogue : getattr(Table,CRIT_COMP).
- """
- def __init__(self, table, column, typ=None):
- """Constructeur (objet Table associé, paramètre de la colonne, type du
- paramètre).
- """
- self.Table=table
- self.rows=self.Table.rows
- self.para=column
- self.type=typ
- self.titr=''
-
- def _extract(self, fun):
- """Construit une table avec les lignes de self.Table
- dont l'élément de nom self.para satisfait le critère fun,
- fun est une fonction qui retourne vrai ou faux
- """
- return Table([row for row in self.Table if fun(row.get(self.para))], self.Table.para, self.Table.type, self.Table.titr)
-
- def __le__(self, VALE):
- return self._extract(lambda v: v<>None and v<=VALE)
-
- def __lt__(self, VALE):
- return self._extract(lambda v: v<>None and v<VALE)
-
- def __ge__(self, VALE):
- return self._extract(lambda v: v<>None and v>=VALE)
-
- def __gt__(self, VALE):
- return self._extract(lambda v: v<>None and v>VALE)
-
- def __eq__(self, VALE, CRITERE='RELATIF', PRECISION=0.):
- if type(VALE) in EnumTypes :
- return self._extract(lambda v: v in VALE)
- if PRECISION==0. or not type(VALE) in NumberTypes:
- if type(VALE) in StringTypes:
- return self._extract(lambda v: v<>None and str(v).strip()==VALE.strip())
- else:
- return self._extract(lambda v: v==VALE)
- else:
- if CRITERE=='ABSOLU':
- vmin=VALE-PRECISION
- vmax=VALE+PRECISION
- else:
- vmin=(1.-PRECISION)*VALE
- vmax=(1.+PRECISION)*VALE
- return self._extract(lambda v: v<>None and vmin<v<vmax)
-
- def __ne__(self, VALE, CRITERE='RELATIF', PRECISION=0.):
- if type(VALE) in EnumTypes :
- return self._extract(lambda v: v not in VALE)
- if PRECISION==0. or not type(VALE) in NumberTypes:
- if type(VALE) in StringTypes:
- return self._extract(lambda v: v<>None and str(v).strip()<>VALE.strip())
- else:
- return self._extract(lambda v: v<>VALE)
- else:
- if CRITERE=='ABSOLU':
- vmin=VALE-PRECISION
- vmax=VALE+PRECISION
- else:
- vmin=(1.-PRECISION)*VALE
- vmax=(1.+PRECISION)*VALE
- return self._extract(lambda v: v<>None and (v<vmin or vmax<v))
-
- def MAXI(self):
- # important pour les performances de récupérer le max une fois pour toutes
- maxi=max(self)
- return self._extract(lambda v: v==maxi)
-
- def MINI(self):
- # important pour les performances de récupérer le min une fois pour toutes
- mini=min(self)
- return self._extract(lambda v: v==mini)
-
- def ABS_MAXI(self):
- # important pour les performances de récupérer le max une fois pour toutes
- abs_maxi=max([abs(v) for v in self.values() if type(v) in NumberTypes])
- return self._extract(lambda v: v==abs_maxi or v==-abs_maxi)
-
- def ABS_MINI(self):
- # important pour les performances de récupérer le min une fois pour toutes
- abs_mini=min([abs(v) for v in self.values() if type(v) in NumberTypes])
- # tester le type de v est trop long donc pas de abs(v)
- return self._extract(lambda v: v==abs_mini or v==-abs_mini)
-
- def __iter__(self):
- """Itère sur les éléments de la colonne"""
- for row in self.Table:
- # si l'élément n'est pas présent on retourne None
- yield row.get(self.para)
- #yield row[self.para]
-
- def __getitem__(self, i):
- """Retourne la ième valeur d'une colonne"""
- return self.values()[i]
-
- def values(self):
- """Renvoie la liste des valeurs"""
- return [r[self.para] for r in self.Table]
-
- # équivalences avec les opérateurs dans Aster
- LE=__le__
- LT=__lt__
- GE=__ge__
- GT=__gt__
- EQ=__eq__
- NE=__ne__
- def VIDE(self) : return self.__eq__(None)
- def NON_VIDE(self): return self.__ne__(None)
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-def sort_table(rows,l_para,w_para,reverse=False):
- """Sort list of dict.
- rows : list of dict
- l_para : list of the keys of dict
- w_para : keys of the sort
- """
- c_para=[i for i in l_para if i not in w_para]
- new_rows=rows
- for i in w_para :
- new_key= '__'+str(w_para.index(i))+i
- for row in new_rows :
- row[new_key]=row[i]
- del row[i]
- for i in c_para :
- new_key= '___'+i
- for row in new_rows :
- row[new_key]=row[i]
- del row[i]
- new_rows.sort()
- if reverse:
- new_rows.reverse()
- for i in w_para :
- old_key= '__'+str(w_para.index(i))+i
- for row in new_rows :
- row[i]=row[old_key]
- del row[old_key]
- for i in c_para :
- old_key= '___'+i
- for row in new_rows :
- row[i]=row[old_key]
- del row[old_key]
- return new_rows
-
-# ------------------------------------------------------------------------------
-def FMT(dform, nform, typAster=None, larg=0, val=''):
- """Retourne un format d'impression Python à partir d'un type Aster ('R','I',
- 'K8', 'K16'...). Si typAster==None, retourne dform[nform].
- larg : largeur minimale du format (val permet de ne pas ajouter des blancs
- si la chaine à afficher est plus longue que le format, on prend le partie
- de ne pas tronquer les chaines)
- """
- if typAster==None:
- fmt=dform[nform]
- elif typAster in ('I', 'R'):
- if nform=='formK':
- # convertit %12.5E en %-12s
- fmt=re.sub('([0-9]+)[\.0-9]*[diueEfFgG]+','-\g<1>s',dform['form'+typAster])
- #print nform, typAster, fmt
- else:
- fmt=dform[nform]
- else:
- # typAster = Kn
- fmt='%-'+typAster[1:]+'s'
- # on ajoute éventuellement des blancs pour atteindre la largeur demandée
- if larg<>0:
- fmt=' '*max(min(larg-len(val),larg-len(fmt % 0)),0) + fmt
- return fmt
-
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-# ------------------------------------------------------------------------------
-if __name__ == "__main__":
- listdic = [
- {'NOEUD': 'N1' ,'NUME_ORDRE': 1 ,'INST': 0.5, 'DX': -0.00233, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N1' ,'NUME_ORDRE': 2 ,'INST': 1.0, 'DX': -0.00467, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N1' ,'NUME_ORDRE': 3 ,'INST': 1.5, 'DX': -0.00701, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N1' ,'NUME_ORDRE': 4 ,'INST': 2.0, 'DX': -0.00934, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N1' ,'NUME_ORDRE': 5 ,'INST': 2.5, 'DX': -0.01168, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N2' ,'NUME_ORDRE': 11,'INST': 5.5, 'DX': -0.00233, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N2' ,'NUME_ORDRE': 12,'INST': 6.0, 'DX': -0.00467, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N2' ,'NUME_ORDRE': 13,'INST': 6.5, 'DX': -0.00701, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N2' ,'NUME_ORDRE': 14,'INST': 7.0, 'DX': -0.00934, 'COOR_Y': 0.53033,},
- {'NOEUD': 'N2' ,'NUME_ORDRE': 15,'INST': 7.5, 'DX': -0.01168, 'COOR_Y': 0.53033,},
- ]
- import random
- random.shuffle(listdic)
- listpara=['NOEUD','NUME_ORDRE','INST','COOR_Y','DX']
- listtype=['K8','I','R','R','R']
- t=Table(listdic,listpara,listtype)
-
- tb=t[('NOEUD','DX')]
- print tb.para
- print tb.type
-
- print
- print "------Table initiale----"
- print t
- print
- print "--------- CRIT --------"
- print t.NUME_ORDRE <=5
- print
- print "------- CRIT & CRIT -----"
- print (t.NUME_ORDRE < 10) & (t.INST >=1.5)
- print
- print "----- EQ maxi / min(col), max(col) ------"
- print t.DX == max(t.DX)
- print min(t.DX)
- print max(t.DX)
- print "------ getitem sur 2 paramètres ------"
- print t.NUME_ORDRE
- print t.DX
- print t['DX','NUME_ORDRE']
- print "------ sort sur INST ------"
- t.sort('INST')
- print t
-
- print "------- TABLEAU_CROISE ------"
- tabc=t['NOEUD','INST','DX']
- tabc.Impr(FORMAT='TABLEAU_CROISE')
-
- N=5
- ldic=[]
- for i in range(N):
- ldic.append({'IND':float(i), 'VAL' : random.random()*i})
- para=['IND','VAL']
- t3=Table(ldic, para, titr='Table aléatoire')
- col=t3.VAL.ABS_MAXI()
- col=t3.VAL.MINI()
-
- t3.sort('VAL','IND')
-
- tg=tabc['INST','DX'].DX.NON_VIDE()
- #tg.Impr(FORMAT='XMGRACE')
-
- g=Graph.Graph()
- g.Titre="Tracé d'une fonction au format TABLEAU"
- g.AjoutCourbe(Val=[tg.INST.values(), tg.DX.values()], Lab=['INST','DX'])
- g.Trace(FORMAT='TABLEAU')
-
-# t.Impr(PAGINATION='NOEUD')
- t.Impr(PAGINATION=('NOEUD','INST'))
-
+++ /dev/null
-#@ MODIF UniteAster Utilitai DATE 11/05/2005 AUTEUR MCOURTOI M.COURTOIS
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2005 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-import types
-
-import aster
-from Cata.cata import _F, DEFI_FICHIER, INFO_EXEC_ASTER, DETRUIRE
-
-#-------------------------------------------------------------------------------
-#-------------------------------------------------------------------------------
-#-------------------------------------------------------------------------------
-class UniteAster:
- """Classe pour manipuler les fichiers en Python en accord avec les unités
- logiques utilisées en Fortran.
- De manière analogue au Fortran, les états possibles d'une unité sont :
- 'F' : fermé, 'O' : ouvert, 'R' : réservé.
-
- Méthodes :
- Nom : Retourne le nom du fichier associé à une unité,
- Etat : Retourne l'état d'une unité,
- Libre : Retourne un numéro d'unité libre,
- EtatInit : Remet une, plusieurs ou toutes les unités dans leur état initial.
-
- Méthode privée :
- _setinfo : pour remplir le dictionnaire des 'infos'
- Attribut privé :
- infos[numéro unité] = { 'nom' : x, 'etat' : x , 'etat_init' : x }
- """
-#-------------------------------------------------------------------------------
- def __init__(self):
- """Initialise le dictionnaire des unités.
- """
- self.infos = {}
-
-#-------------------------------------------------------------------------------
- def _setinfo(self, ul):
- """Remplit les infos de l'unité 'ul'.
- """
- # ul peut etre un entier Aster
- try:
- unit = ul.valeur
- except:
- unit = int(ul)
- # Si la clé n'existe pas
- ini = False
- if not self.infos.has_key(unit):
- self.infos[unit] = {}
- self.infos[unit]['nom'] = ''
- self.infos[unit]['etat'] = '?'
- self.infos[unit]['etat_init'] = '?'
- ini = True
-
- __tab=INFO_EXEC_ASTER(UNITE=unit, LISTE_INFO=('ETAT_UNITE'))
-
- # O:ouvert, F:fermé, R:réservé
- self.infos[unit]['etat'] = __tab['ETAT_UNITE',1].strip()[0]
- if ini:
- self.infos[unit]['etat_init'] = self.infos[unit]['etat']
-
- # nom du fichier
- if self.infos[unit]['etat'] in ['O', 'R']:
- nomfich=''.join([__tab['NOMFIC%d' % i,1] for i in range(1,5)]).strip()
- elif self.infos[unit]['etat'] == 'F':
- nomfich='fort.'+str(unit)
- else:
- message = "Etat de l'unité inconnu : %s" % self.infos[unit]['etat']
- print __tab.EXTR_TABLE()
- raise aster.FatalError,"<F> <UniteAster._setinfo> %s" % message
- self.infos[unit]['nom'] = nomfich
- DETRUIRE(CONCEPT=_F(NOM=__tab))
-
-#-------------------------------------------------------------------------------
- def Libre(self, nom=None):
- """Réserve et retourne une unité libre en y associant, s'il est fourni,
- le fichier 'nom'.
- """
- __tab=INFO_EXEC_ASTER(LISTE_INFO=('UNITE_LIBRE'))
- unit = __tab['UNITE_LIBRE',1]
- DETRUIRE(CONCEPT=_F(NOM=__tab))
- if nom==None:
- nom='fort.'+str(unit)
-
- # Si la clé existe, c'est que le fichier n'était pas libre
- if self.infos.has_key(unit):
- message = "Cette unité est déjà affectée au fichier %s" % \
- self.infos[unit]['nom']
- raise aster.FatalError,"<F> <UniteAster.Libre> %s" % message
-
- DEFI_FICHIER(ACTION='RESERVER', UNITE=unit , FICHIER=nom.strip())
- self.infos[unit] = {}
- self.infos[unit]['nom'] = nom.strip()
- self.infos[unit]['etat'] = 'R'
- self.infos[unit]['etat_init'] = 'F'
- return unit
-
-#-------------------------------------------------------------------------------
- def Nom(self, ul):
- """Retourne le nom du fichier associé à l'unité 'ul'.
- """
- # ul peut etre un entier Aster
- try:
- unit = ul.valeur
- except:
- unit = int(ul)
- # Si la clé n'existe pas
- if not self.infos.has_key(unit):
- self._setinfo(unit)
- return self.infos[unit]['nom']
-
-#-------------------------------------------------------------------------------
- def Etat(self, ul, **kargs):
- """Retourne l'état de l'unité si 'etat' n'est pas fourni
- et/ou change son état :
- kargs['etat'] : nouvel état,
- kargs['TYPE'] : type du fichier à ouvrir ASCII/BINARY/LIBRE,
- kargs['ACCES'] : type d'accès NEW/APPEND/OLD.
- """
- # ul peut etre un entier Aster
- try:
- unit = ul.valeur
- except:
- unit = int(ul)
- # Si la clé n'existe pas
- if not self.infos.has_key(unit):
- self._setinfo(unit)
- if not kargs.has_key('etat'):
- return self.infos[unit]['etat']
-
- # En fonction de la demande, on bascule son état ou pas
- new = kargs.get('etat')
- if not new in ['R', 'F', 'O']:
- message = "Nouvel état de l'unité incorrect : %s" % new
- raise aster.FatalError,"<F> <UniteAster.Etat> %s" % message
-
- if self.infos[unit]['etat'] == new:
- pass
- elif new == 'R':
- if self.infos[unit]['etat'] == 'O':
- DEFI_FICHIER(ACTION='LIBERER', UNITE=unit)
- DEFI_FICHIER(ACTION = 'RESERVER',
- UNITE = unit,
- FICHIER = self.infos[unit]['nom'])
- elif new == 'F':
- DEFI_FICHIER(ACTION='LIBERER', UNITE=unit)
- elif new == 'O':
- if self.infos[unit]['etat'] == 'R':
- DEFI_FICHIER(ACTION='LIBERER', UNITE=unit)
- DEFI_FICHIER(ACTION ='ASSOCIER',
- UNITE = unit,
- FICHIER = self.infos[unit]['nom'],
- TYPE = kargs.get('TYPE', 'ASCII'),
- ACCES = kargs.get('ACCES', 'APPEND'),)
- self.infos[unit]['etat'] = new
- return self.infos[unit]['etat']
-
-#-------------------------------------------------------------------------------
- def EtatInit(self, ul=None):
- """Remet l'unité 'ul' dans son état initial.
- Si 'ul' est omis, toutes les unités sont remises dans leur état initial.
- """
- if ul == None:
- for uli, vul in self.infos.items():
- self.Etat(uli, etat=vul['etat_init'])
- else:
- if not type(ul) in [types.ListType, types.TupleType]:
- ul=[ul,]
- for u in ul:
- # u peut etre un entier Aster
- try:
- unit = u.valeur
- except:
- unit = int(u)
- # Si la clé n'existe pas
- if not self.infos.has_key(unit):
- self._setinfo(unit)
- else:
- self.Etat(unit, etat=self.infos[unit]['etat_init'])
+++ /dev/null
-#@ MODIF Utmess Utilitai DATE 30/11/2004 AUTEUR MCOURTOI M.COURTOIS
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-import sys
-import aster
-
-def UTMESS(code, sprg, texte):
- """Utilitaire analogue à la routine fortran UTMESS.
- code : 'A', 'E', 'S', 'F'
- sprg : nom du module, classe ou fonction python où l'on se trouve
- texte : contenu du message
- """
- fmt='\n <%s> <%s> %s\n\n'
- UL={
- 'MESSAGE' : 6,
- 'RESULTAT' : 8,
- #'ERREUR' : 9,
- }
- # On importe la définition des commandes à utiliser dans la macro
-# if jdc:
-# DEFI_FICHIER = jdc.get_cmd('DEFI_FICHIER')
-# else:
-# # on se limite au print !
-# UL={ 'MESSAGE' : 6, }
- try:
- from Cata.cata import DEFI_FICHIER
- except ImportError:
- # on se limite au print !
- UL={ 'MESSAGE' : 6, }
-
- reason=fmt % (code, sprg, texte)
-
- for nom,ul in UL.items():
- if ul<>6:
- DEFI_FICHIER(ACTION='LIBERER', UNITE=ul, )
- f=open('fort.'+str(ul),'a')
- else:
- f=sys.stdout
- # écriture du message
- f.write(reason)
-
- if ul<>6:
- f.close()
- DEFI_FICHIER(ACTION='ASSOCIER', UNITE=ul, TYPE='ASCII', ACCES='APPEND')
-
- if code=='S':
- raise aster.error, reason
- elif code=='F':
- raise aster.FatalError, reason
+++ /dev/null
-#@ MODIF __init__ Utilitai DATE 20/09/2004 AUTEUR DURAND C.DURAND
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2003 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-
+++ /dev/null
-#@ MODIF courbes Utilitai DATE 14/09/2004 AUTEUR MCOURTOI M.COURTOIS
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-#==================================================
-# fonction "COURBES"
-# usage : permet de tracer des courbes en interactif
-# avec XMGRACE ou dans un fichier postscript
-#==================================================
-
-import Stanley
-from Stanley import xmgrace
-from Stanley import as_courbes
-
-def COURBES(listcourb,titre=' ',soustitre=' ',legx=' ',legy=' ',bornex=None,borney=None,fichier=None):
-
-# ARGUMENTS
-
-# listcourb : tuple de courbes, chaque courbe etant definie soit par
-# (TABLE1, NOM_PARA_X, TABLE2, NOM_PARA_Y, LEGENDE)
-# soit par :
-# (FONCTION,LEGENDE)
-# titre et sous_titre : facultatifs, titre et sous-tritre du graphique
-# legx, legy : facultatifs, legendes des axes
-# bornex, borney : facultatifs, bornes sur les axes
-# fichier : facultatif : sortie au format postscript si present
-#
-# exemples d'appel :
-#--------------------
-# courb1=(SYYPRT,'ABSC_CURV',SYYPRT,'SIYY','PRT')
-# courb2=(SYYMLC10,'ABSC_CURV',SYYMLC10,'SIYY','MLC10')
-# courb3=(SYYML100,'ABSC_CURV',SYYML100,'SIYY','MLC100')
-# listcourb=(courb1,courb2,courb3)
-# COURBES(listcourb,titre='Plaque trouee',legx='Abcisses curvilignes',legy='Contraintes (MPa)',bornex=(0,100),borney=(500,1000))
-# fonc1=(F1,'F_PRT')
-# fonc2=(F2,'F_MLC10')
-# fonc3=(F3,'F_MLC100')
-# listfonc=(fonc1,fonc2,fonc3)
-# COURBES(listfonc,titre='Fonctions')
-# postscript
-# COURBES(listfonc,titre='Plaque trouee',fichier='./fort.24')
-#--------------------------------------------------------------
-
-# initialisation du trace de courbes
-
- if (fichier!=None):
- graphe=xmgrace.Xmgr(10,' -hardcopy -nosafe')
- print "Nombre de courbes ",len(listcourb)," sur le fichier :",fichier
-
- else:
- graphe=xmgrace.Xmgr(10,' -noask')
- print "Nombre de courbes ",len(listcourb)
-
- graphe.Nouveau_graphe()
-
-# dimensionnement des axes
- if bornex != None :
- xmin=list(bornex)[0]
- xmax=list(bornex)[1]
- ctest1 = as_courbes.Courbe()
- ctest1.x=[xmin,xmax]
- ctest1.y=[0.0,0.0]
- graphe.Courbe(ctest1)
-
- if borney != None :
- ymin=list(borney)[0]
- ymax=list(borney)[1]
- ctest2 = as_courbes.Courbe()
- ctest2.x=[0.0,0.0]
- ctest2.y=[ymin,ymax]
- graphe.Courbe(ctest2)
-
- if titre != None :
- if soustitre != None :
- graphe.Titre(titre,soustitre)
- else :
- graphe.Titre(titre,' ')
-
- if legx != None :
- graphe.Axe_x(legx)
-
- if legy != None :
- graphe.Axe_y(legy)
-
- k = 0
-
- for courbi in listcourb:
- sigi = as_courbes.Courbe()
-
- try :
- # cas d une table
- sigi.Lire_x(courbi[0],courbi[1])
- sigi.Lire_y(courbi[2],courbi[3])
- legende=courbi[4]
- except :
- # cas d une fonction
- sigi.x,sigi.y=courbi[0].Valeurs()
- legende=courbi[1]
-
- graphe.Courbe(sigi,legende)
- graphe.Send('WITH G'+repr(graphe.gr_act))
- graphe.Send('S' + str(k) + ' SYMBOL ' + str(k+2))
- graphe.Send('S' + str(k) + ' SYMBOL SIZE 0.5')
- graphe.Send('S' + str(k) + ' SYMBOL COLOR '+str(k+2))
- graphe.Send('S' + str(k) + ' LINE COLOR '+str(k+2))
- k = k + 1
- graphe.Send('REDRAW')
-
- if (fichier!=None):
- graphe.Sortie_EPS(fichier)
- graphe.Fermer()
- else:
- graphe.Attendre()
-
- k=0
-
-#===========================================
-
-
+++ /dev/null
-#@ MODIF funct_root Utilitai DATE 14/09/2004 AUTEUR MCOURTOI M.COURTOIS
-# -*- coding: iso-8859-1 -*-
-################################################################################
-# Mathematical utility routines
-# Copyright (C) 1999, Wesley Phoa
-#
-# Reference: Numerical Recipes in C
-# [[[[extraits]]]]
-
-class BracketingException(Exception):
- pass
-
-class RootFindingException(Exception):
- pass
-
-class MinimizationException(Exception):
- pass
-
-GOLDEN = (1+5**.5)/2
-
-#
-# MISCELLANEOUS
-#
-
-def sgn(x):
- if x==0:
- return 0
- else:
- return x/abs(x)
-
-#
-# UNIVARIATE ROOT FINDING
-#
-
-def bracket_root(f, interval, max_iterations=50):
- """\
-Given a univariate function f and a tuple interval=(x1,x2),
-return a new tuple (bracket, fnvals) where bracket=(x1,x2)
-brackets a root of f and fnvals=(f(x1),f(x2)).
- """
- x1, x2 = interval
- if x1==x2:
- raise BracketingException("initial interval has zero width")
- elif x2<x1:
- x1, x2 = x2, x1
- f1, f2 = f(x1), f(x2)
- for j in range(max_iterations):
- while f1*f2 >= 0: # not currently bracketed
- if abs(f1)<abs(f2):
- x1 = x1 + GOLDEN*(x1-x2)
- else:
- x2 = x2 + GOLDEN*(x2-x1)
- f1, f2 = f(x1), f(x2)
- return (x1, x2), (f1, f2)
- raise BracketingException("too many iterations")
-
-def ridder_root(f, bracket, fnvals=None, accuracy=1e-6, max_iterations=50):
- """\
-Given a univariate function f and a tuple bracket=(x1,x2) bracketing a root,
-find a root x of f using Ridder s method. Parameter fnvals=(f(x1),f(x2)) is optional.
- """
- x1, x2 = bracket
- if fnvals==None:
- f1, f2 = f(x1), f(x2)
- else:
- f1, f2 = fnvals
- if f1==0:
- return x1
- elif f2==0:
- return x2
- elif f1*f2>=0:
- raise BracketingException("initial interval does not bracket a root")
- x4 = 123456789.
- for j in range(max_iterations):
- x3 = (x1+x2)/2
- f3 = f(x3)
- temp = f3*f3 - f1*f2
- x4, x4old = x3 + (x3-x1)*sgn(f1-f2)*f3/temp**.5, x4
- f4 = f(x4)
- if f1*f4<0: # x1 and x4 bracket root
- x2, f2 = x4, f4
- else: # x4 and x2 bracket root
- x1, f1 = x4, f4
- if min(abs(x1-x2),abs(x4-x4old))<accuracy or temp==0:
- return x4
- raise RootFindingException("too many iterations")
-
-def root(f, interval=(0.,1.), accuracy=1e-4, max_iterations=50):
- """\
-Given a univariate function f and an optional interval (x1,x2),
-find a root of f using bracket_root and ridder_root.
- """
- bracket, fnvals = bracket_root(f, interval, max_iterations)
- return ridder_root(f, bracket, fnvals, accuracy, max_iterations)
-
+++ /dev/null
-#@ MODIF partition Utilitai DATE 10/02/2005 AUTEUR ASSIRE A.ASSIRE
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-# RESPONSABLE ASSIRE A.ASSIRE
-
-import aster
-import string, os, time, sys, UserList, types
-import Numeric as N
-
-from Accas import _F
-from Noyau.N_utils import AsType
-
-# ------------------------------------------------------------------ #
-
-def enleve_doublons_liste(liste):
- """ A partir d une liste qui peut contenir des doublons, on renvoie une liste sans
- doublons (et qui reviendra triée)
- """
-
- """ Au cas ou ca ne serait pas deja un vecteur numpy """ # Exemple
- liste=N.sort(N.array(liste,copy=0)) # [1, 2, 2, 3, 3, 4, 5]
- liste_diff=liste[1:]-liste[:-1] # [1, 0, 1, 0, 1, 1]
- liste_temp=N.nonzero(liste_diff) # [0, 2, 4, 5]
- liste_indice=N.zeros(liste_temp.shape[0]+1)
- liste_indice[0]=0
- liste_indice[1:]=liste_temp+1 # [0, 1, 3, 5, 6]
- liste2=N.take(liste,liste_indice) # [1, 2, 3, 4, 5]
- return liste2
-
-
-# ============================================================================ #
-
-
-class CONNEC(UserList.UserList) :
- """
- Connectivite : sequence mutable de sequences de noeuds
- Pour l'instant, on s'appuie sur une liste de liste.
- """
-
- def __init__(self,nma) :
-
- UserList.UserList.__init__(self,[None]*nma)
-
-
- def Index(self) :
- """
- Retourne la connectivite sous forme de deux vecteurs numpy :
- ind -> tableau des index (y compris le n+1 eme)
- noe -> tableau des listes de noeuds
- """
-
- # Dimension des deux vecteurs
- nma = len(self)
- ltot = N.reduce(lambda x,y : x+len(y), self,0)
- ind = N.zeros(nma+1,Int)
- noe = N.zeros(ltot,Int)
-
- # Construction des vecteurs
- ind[0] = 0
- for i in range(nma) :
- m = self[i]
- ind[i+1] = ind[i] + len(m)
- noe[ind[i]:ind[i+1]] = N.array(m)
-
- return ind,noe
-
-
-# ============================================================================ #
-# ============================================================================ #
-
-class MAIL_PY :
-
- """
- SD PYTHON MAILLAGE
- La numeration est 0..N-1 pour les noeuds et 0..M-1 pour les mailles
- """
-
- def __init__(self,nno=0,nma=0) :
-
- self.cn = N.zeros((nno,3),N.Float)
- self.tm = N.zeros(nma,N.Int)
- self.co = CONNEC(nma)
- self.gma = {}
- self.gno = {}
- self.indice_noeuds = N.arange(nno)
-
- self.correspondance_noeuds = []
- self.correspondance_mailles = []
-
- try:
- import aster
- nom_mailles = (None,) + string.strip(aster.getvectjev('&CATA.TM.NOMTM'))
- except:
- nom_mailles = (None,
- 'POI1', 'SEG2', 'SEG22', 'SEG3', 'SEG33', 'SEG4', 'TRIA3',
- 'TRIA33', 'TRIA6', 'TRIA66', 'TRIA7', 'QUAD4', 'QUAD44', 'QUAD8',
- 'QUAD88', 'QUAD9', 'QUAD99', 'TETRA4', 'TETRA10','PENTA6', 'PENTA15',
- 'PYRAM5', 'PYRAM13','HEXA8', 'HEXA20', 'HEXA27', 'TR3QU4', 'QU4TR3',
- 'TR6TR3', 'TR3TR6', 'TR6QU4', 'QU4TR6', 'TR6QU8', 'QU8TR6', 'TR6QU9',
- 'QU9TR6', 'QU8TR3', 'TR3QU8', 'QU8QU4', 'QU4QU8', 'QU8QU9', 'QU9QU8',
- 'QU9QU4', 'QU4QU9', 'QU9TR3', 'TR3QU9', 'SEG32', 'SEG23' )
-
- dic_mailles = {}
- for i in range(len(nom_mailles)) :
- dic_mailles[nom_mailles[i]] = i
-
- self.nom = nom_mailles
- self.dic = dic_mailles
-
- try:
- psyco.bind(self.FromAster)
- except: pass
-
-
-# -------------------------------------------------------------
-
- def get_connexite(self, nom, nma):
- co=CONNEC(nma)
- dico_connexite = aster.getcolljev(nom)
- for element in dico_connexite.keys() :
- co[int(element)-1] = (N.array(dico_connexite[element])-1)
- return co
-
-
- def get_coordonnees_noeuds(self, nom, nombre_noeuds):
- coordonnees_noeuds = aster.getvectjev(nom)
- coordonnees_noeuds = N.array(coordonnees_noeuds)
- coordonnees_noeuds.shape = (nombre_noeuds,3)
- cn = coordonnees_noeuds
- return cn
-
-
-# -------------------------------------------------------------
-
- def FromAster(self,nom) :
-
- # On accepte le concept Aster ou bien la chaine texte de son nom
- if type(nom)!=types.StringType:
- nom_maillage = nom.nom
- else:
- nom_maillage = nom
-
- nom_maillage=string.ljust(nom_maillage,8)
-
- # recuperation de la taille
- self.dime_maillage = aster.getvectjev(nom_maillage+'.DIME')
- nombre_noeuds = self.dime_maillage[0]
- nombre_mailles = self.dime_maillage[2]
-
- # coordonnees des noeuds
- self.cn = self.get_coordonnees_noeuds(nom_maillage+'.COORDO .VALE', nombre_noeuds)
-
- # type de maille
- self.tm = N.array(aster.getvectjev(nom_maillage+'.TYPMAIL'))
-
- # connexite
- self.co = self.get_connexite(nom_maillage+'.CONNEX', nombre_mailles)
-
- self.indice_noeuds=N.arange(nombre_noeuds)
-
- # groupe de noeuds
- Lno_groupno_tot = aster.getcolljev(nom_maillage+'.GROUPENO')
-
- Lno_groupno={}
- for key in Lno_groupno_tot :
- # Tolerance car parfois defi_group crée des groupes nuls à clé entiere
- try:
- Lno_groupno[key.strip()]=N.array(Lno_groupno_tot[key])-1
- except: pass
- self.gno=Lno_groupno
-
- # groupe de mailles
- Lma_groupma_tot = aster.getcolljev(nom_maillage+'.GROUPEMA')
- Lma_groupma={}
- for key in Lma_groupma_tot :
- Lma_groupma[key.strip()]=N.array(Lma_groupma_tot[key])-1
- self.gma=Lma_groupma
-
- del(Lma_groupma_tot)
-
- # listes de correspondance entre Aster et Mail-Py
- self.correspondance_noeuds = aster.getvectjev(nom_maillage+'.NOMNOE')
- self.correspondance_mailles = aster.getvectjev(nom_maillage+'.NOMMAI')
-
-
-# -------------------------------------------------------------
-
- def ToAster(self) :
-
- try:
- INFO_EXEC_ASTER = self.jdc.get_cmd('INFO_EXEC_ASTER')
- DETRUIRE = self.jdc.get_cmd('DETRUIRE')
- LIRE_MAILLAGE = self.jdc.get_cmd('LIRE_MAILLAGE')
- except:
- try:
- from Cata.cata import INFO_EXEC_ASTER
- from Cata.cata import DETRUIRE
- from Cata.cata import LIRE_MAILLAGE
- except:
- print "\n\nERREUR : il faut lancer ce programme depuis Aster pour pouvoir \ngénérer un maillage Aster.\n\n"
- sys.exit(1)
-
- # Recuperation d'une unité logique libre
- _UL=INFO_EXEC_ASTER(LISTE_INFO='UNITE_LIBRE')
- ul1=_UL['UNITE_LIBRE',1]
- DETRUIRE(CONCEPT=(_F(NOM=_UL),), INFO=1)
-
- # Sauvegarde du maillage dans un fichier .mail
- fichier = 'fort.'+str(ul1)
- f = open(fichier, 'w')
- f.write( self.Voir_Mail() )
- f.close()
-
- # Récupération des concepts Aster présents dans la base
-# Note AA : ne marche pas encore bien :
-# Limitation : ne gere qu'un seul maillage genere par ToAster
-#
-# nom='_MSH'
-# try:
-# self.jdc_recupere=CONTEXT.get_current_step()
-# t_maillage=[]
-# for i in self.jdc_recupere.sds_dict.keys( ):
-# if self.jdc_recupere.sds_dict[i].__class__.__name__ == 'maillage_sdaster':
-# if (mail[0:4]==nom):
-# t_maillage.append( i )
-#
-# num=len(_lst)+1
-# except:
-# num=0
-
- # Lecture de ce fichier .mail pour l'injecter dans l'espace Aster
- _SMESH_ = LIRE_MAILLAGE(UNITE = ul1)
-
- return(_SMESH_)
-
-# -------------------------------------------------------------
-
- def __str__(self) :
- return self.Voir_Mail()
-
-# -------------------------------------------------------------
-
- def Voir_Mail(self) :
- """
- Impression au format ASTER
- """
-
- l = []
- l.append('TITRE')
- l.append('% CLASSE PY_MAIL -> MAIL')
- l.append('FINSF')
- l.append('%')
-
- (nno,ndim) = self.cn.shape
-
-
- # Coordonnees des noeuds
- l.append('COOR_3D')
-
- # Si le maillage initial ne provient pas d'Aster
- if len(self.correspondance_noeuds) == 0:
- for i in range(nno) :
- ch = 'N'+repr(i)+' '+`self.cn[i,0]` + ' ' + `self.cn[i,1]` + ' ' + `self.cn[i,2]`
-# ch = 'N'+repr(i+1)+' '+`self.cn[i,0]` + ' ' + `self.cn[i,1]` + ' ' + `self.cn[i,2]`
- l.append(ch)
-
- # Si le maillage initial provient d'Aster
- else:
- for i in range(nno) :
- ch = self.correspondance_noeuds[i]+' '+`self.cn[i,0]` + ' ' + `self.cn[i,1]` + ' ' + `self.cn[i,2]`
-# ch = 'N'+repr(i+1)+' '+`self.cn[i,0]` + ' ' + `self.cn[i,1]` + ' ' + `self.cn[i,2]`
- l.append(ch)
-
- # Connectivité des mailles
- ind = N.argsort(self.tm)
- ty = 0
-
- # Si le maillage initial ne provient pas d'Aster
- if len(self.correspondance_mailles) == 0:
- for m in ind :
- if self.tm[m] <> ty :
- l.append('FINSF') ; l.append('%')
- ty = self.tm[m]
- l.append(self.nom[ty])
- ch = 'M'+`m`+' '
-# ch = 'M'+`m+1`+' '
- for n in self.co[m] :
- ch = ch + 'N'+`n` + ' '
-# ch = ch + 'N'+`n+1` + ' '
- l.append(ch)
-
- # Si le maillage initial provient d'Aster
- else:
- for m in ind :
- if self.tm[m] <> ty :
- l.append('FINSF') ; l.append('%')
- ty = self.tm[m]
- l.append(self.nom[ty])
- ch = self.correspondance_mailles[m]+' '
-# ch = 'M'+`m+1`+' '
- for n in self.co[m] :
- ch = ch + self.correspondance_noeuds[n] + ' '
-# ch = ch + 'N'+`n+1` + ' '
- l.append(ch)
-
- l.append('FINSF') ; l.append('%')
-
-
- # Group_ma et Group_no
- entete = ['GROUP_MA','GROUP_NO']
- d_gp = [self.gma,self.gno]
- pref = ['M','N']
-
- # Si le maillage initial ne provient pas d'Aster
- if len(self.correspondance_mailles) == 0:
-
- for (d_gp,entete,prefixe) in [(self.gma,'GROUP_MA','M'),(self.gno,'GROUP_NO','N')] :
- for gp in d_gp :
- if len(d_gp[gp])>0: # On ne prend en compte que les group_* non vides
- l.append(entete)
- l.append(' ' + gp)
- ch = ' '
- i=0
- for o in d_gp[gp]:
- i+=1 # on ne met que 8 mailles sur une meme ligne
- if (len(ch) > 60 or i>7):
- l.append(ch)
- ch = ' '
- i=0
- ch = ch + prefixe + `o` + ' '
- l.append(ch)
- l.append('FINSF') ; l.append('%')
-
- # Si le maillage initial provient d'Aster
- else:
-
- for (d_gp,entete,prefixe) in [(self.gma,'GROUP_MA','M'),(self.gno,'GROUP_NO','N')] :
- for gp in d_gp :
- if len(d_gp[gp])>0: # On ne prend en compte que les group_* non vides
- l.append(entete)
- l.append(' ' + gp)
- ch = ' '
- i=0
- for o in d_gp[gp]:
- i+=1 # on ne met que 8 mailles sur une meme ligne
- if (len(ch) > 60 or i>7):
- l.append(ch)
- ch = ' '
- i=0
-# ch = ch + prefixe + `o` + ' '
- if prefixe=='M':
- ch = ch + self.correspondance_mailles[o] + ' '
- else:
- ch = ch + self.correspondance_noeuds[o] + ' '
- l.append(ch)
- l.append('FINSF') ; l.append('%')
-
- # Fin
- l.append('FIN\n')
- return string.join(l,'\n')
-
-
-
-
-# ============================================================================ #
-# ============================================================================ #
-
-class PARTITION:
-
- def __init__(self, jdc=None ,nb=0):
-
- self.jdc = jdc
-
- self.fichier_out = ''
- self.liste_mailles = N.array( [] )
- self.liste_sd = N.array( [] )
- self.liste_mailles_bord = []
- self.liste_sd_bord = []
-
- self.MAILLAGE_Python = None
-
- self.RELATIONS = { 'C_plus' : None,
- 'C_moins': None,
- 'nr': 0 }
-
- self.ASTER = { 'MAILLAGE' : None,
- 'MODELE' : None,
- 'GROUP_MA' : None,
- 'GROUP_MA_BORD' : None,
- 'DICO_SD_MAILLES' : None,
- }
-
- self.OPTIONS = { 'NB_PART' : '',
- 'ALGO' : '',
- 'INFO' : '',
- 'rep_metis' : aster.repout(),
- 'exe_metis' : aster.repout() + 'pmetis',
- 'fichier_in' : 'fort.66',
- 'fichier_out' : 'fort.68',
- 'elimine_bords': 'OUI',
- }
-
- self.Creation_Dico_Correspondance_Type_Maille()
-
-
-
-# ---------------------------------------------------------------------------- #
-
- def __str__(self) :
- """
- Impression du contenu de la partition
- """
- l = []
- l.append( 'Contenu de la partition :' )
- l.append( '-------------------------' )
- try: l.append( '- Maillage : ' + str(self.ASTER['MAILLAGE'].nom) )
- except: pass
- try: l.append( '- Modele : ' + str(self.ASTER['MODELE'].nom) )
- except: pass
- l.append( '- Nb part : ' + str(self.OPTIONS['NB_PART']) )
- l.append( '- Niveau INFO : ' + str(self.OPTIONS['INFO']) )
- l.append( '- Liste group_ma : ' + str(self.ASTER['GROUP_MA']) )
-
- return string.join(l,'\n')
-
-# ---------------------------------------------------------------------------- #
-
- def Partitionne_Aster(self, MAILLAGE, NB_PART, MODELE=None, METHODE=None, LOGICIEL=None, INFO=1):
-
- self.t00 = time.clock()
-
- self.OPTIONS['INFO'] = INFO
-
- if MODELE:
- # Recuperation de la liste des mailles à perndre en compte
- self.ASTER['MODELE'] = MODELE
- self.ASTER['MAILLAGE'] = MAILLAGE
- _LST_MA = self.Modele_to_Liste_Mailles(MODELE)
-
- elif MAILLAGE:
- self.ASTER['MAILLAGE'] = MAILLAGE
- _LST_MA = None
-
-
- # Creation du maillage Python correspondant au maillage Aster
- MAILLAGE_Python = MAIL_PY()
- MAILLAGE_Python.FromAster(MAILLAGE.nom)
-
- # Partitionne le maillage Python avec la liste de mailles _LST_MA
- self.Partitionne_Maillage(MAILLAGE_Python, NB_PART, MAILLE=_LST_MA, METHODE=METHODE, LOGICIEL=LOGICIEL, INFO=INFO)
-
- return
-
-
-# ---------------------------------------------------------------------------- #
-
- def Partitionne_Maillage(self, MAILLAGE_Python, NB_PART, MAILLE=None, METHODE=None, LOGICIEL=None, INFO=1):
-
- self.t00 = time.clock()
-
- if METHODE:
- self.OPTIONS['exe_metis'] = aster.repout() + string.lower(METHODE)
- elif LOGICIEL:
- self.OPTIONS['exe_metis'] = LOGICIEL
-
- self.OPTIONS['NB_PART'] = NB_PART
- self.OPTIONS['INFO'] = INFO
- self.MAILLAGE_Python = MAILLAGE_Python
-
- exe_metis = self.OPTIONS['exe_metis']
- f_metis = self.OPTIONS['fichier_in']
- fw_metis = self.OPTIONS['fichier_out']
-
- _LST_MA = MAILLE
-
- # On initialise la connectivité et la connectivité inverse des aretes
- self.MAILLAGE_Python.ca = {}
- self.MAILLAGE_Python.cia = {}
-
- _DIM = self.MAILLAGE_Python.dime_maillage[5]
- _LST_TMA = self.MAILLAGE_Python.tm
-
- if self.OPTIONS['INFO']>=5:
- print 'cn=', self.MAILLAGE_Python.cn
- print 'tm=', self.MAILLAGE_Python.tm
- print 'co=', self.MAILLAGE_Python.co
- print 'gma=', self.MAILLAGE_Python.gma
- print 'gno=', self.MAILLAGE_Python.gno
- print 'dim=', self.MAILLAGE_Python.dime_maillage
- if self.OPTIONS['INFO']>=5: print '_LST_MA=', _LST_MA
-
-
- # Elimination des mailles de bords
- if self.OPTIONS['elimine_bords']!='NON':
-
- # Liste des mailles à prendre en compte : dimension _DIM
- _D_DIM_MAILLES = self.Creation_Listes_Mailles_Par_Dim(self.MAILLAGE_Python.tm, _LST_MA=_LST_MA)
-
- # Connectivité et connectivité inverse sur les bords
- self.Connectivite_Aretes()
-
- self.liste_mailles = _D_DIM_MAILLES[ _DIM ]
-
- # Pour prendre en compte des mélanges d'elements de dimension differente
- _LST, _LST_BD = self.Elimination_Mailles_de_bords(MAILLAGE_Python, _D_DIM_MAILLES, _DIM)
- self.liste_mailles = N.concatenate( (self.liste_mailles,N.array(_LST)) )
-
- if self.OPTIONS['INFO']>=5:
- print '_LST_BD=',_LST_BD
- print '_LST=',_LST
-
- else:
- self.liste_mailles = _LST_MA
-
-
- # Restriction des connectivités aux mailles à prendre en compte
- self.Connectivite_Aretes(OPTION='all', _LST_OK=self.liste_mailles)
-
- # Creation de l'arbre de connectivité des bords
- self.Creation_Graphe()
-
- # Reduction de l'arbre de connectivité des bords
- _nb = self.Reduction_Graphe(_DIM)
-
- # Ecriture du fichier pour Metis/Chaco/Jostle
- _D_CORRES = self.Ecrire_Graphe(f_metis, _nb)
-
- # Lancement de metis sur le fichier fort.UL (production de fort.UL.part.N)
- txt = exe_metis + ' ' + f_metis + ' ' + str(NB_PART)
- print 'Commande : ',txt
- os.system( txt )
-
- # Lecture du fichier resultant de Metis
- self.fichier_out = f_metis + '.part.' + str(NB_PART)
- self.liste_sd = self.Lecture_fichier_sdd(self.fichier_out, self.liste_mailles)
-
- # Traitement des mailles de bords (on les reinjecte dans un SD)
- if self.OPTIONS['elimine_bords']!='NON':
- self.Affectation_Mailles_de_bords(_LST_BD, _DIM)
-
- t1 = time.clock()
- print "--- FIN PARTITIONNEMENT : ", t1 - self.t00
-
- return
-
-
-
-# ---------------------------------------------------------------------------- #
-
- def Creation_Dico_Correspondance_Type_Maille(self):
-
- # TYPE_ELEM : CF. &CATA.TM
- # 1 - >POI1 <>SEG2 <>SEG22 <>SEG3 <>SEG33 <>SEG4 <>TRIA3 <
- # 8 - >TRIA33 <>TRIA6 <>TRIA66 <>TRIA7 <>QUAD4 <>QUAD44 <>QUAD8 <
- # 15 - >QUAD88 <>QUAD9 <>QUAD99 <>TETRA4 <>TETRA10 <>PENTA6 <>PENTA15 <
- # 22 - >PYRAM5 <>PYRAM13 <>HEXA8 <>HEXA20 <>HEXA27 <>TR3QU4 <>QU4TR3 <
- # 29 - >TR6TR3 <>TR3TR6 <>TR6QU4 <>QU4TR6 <>TR6QU8 <>QU8TR6 <>TR6QU9 <
- # 36 - >QU9TR6 <>QU8TR3 <>TR3QU8 <>QU8QU4 <>QU4QU8 <>QU8QU9 <>QU9QU8 <
- # 43 - >QU9QU4 <>QU4QU9 <>QU9TR3 <>TR3QU9 <>SEG32 <>SEG23 <
-
- # Creation du dictionnaire des correspondance type_maille -> liste des aretes
- maille2aretes={}
- # POI
- maille2aretes[1] = [ ]
- # SEG
- maille2aretes[2] = [ (0,1) ]
- maille2aretes[3] = maille2aretes[4] = maille2aretes[5] = maille2aretes[6] = maille2aretes[2]
- # TRIA
- maille2aretes[7] = [ (0,1),(1,2),(0,2) ]
- maille2aretes[8] = maille2aretes[9] = maille2aretes[10] = maille2aretes[11] = maille2aretes[7]
- # QUAD
- maille2aretes[12] = [ (0,1),(1,2),(2,3),(0,3) ]
- maille2aretes[13] = maille2aretes[14] = maille2aretes[15] = maille2aretes[16] = maille2aretes[17] = maille2aretes[12]
- # TETRA
- maille2aretes[18] = [ (0,1,2),(0,1,3),(0,2,3),(1,3,2) ]
- maille2aretes[19] = maille2aretes[18]
- # PENTA
- maille2aretes[20] = [ (0,1,2),(3,4,5),(0,2,5,3),(0,1,4,3),(2,1,4,5) ]
- maille2aretes[21] = maille2aretes[20]
- # PYRAM
- maille2aretes[22] = [ (0,1,4),(1,2,4),(2,3,4),(3,0,4),(0,1,2,3) ]
- maille2aretes[23] = maille2aretes[22]
- # HEXA
- maille2aretes[24] = [ (0,1,2,3), (4,5,6,7), (1,2,6,5), (2,3,7,6), (7,4,0,3), (4,5,1,0) ]
- maille2aretes[25] = maille2aretes[26] = maille2aretes[24]
-
-
- # dictionnaire de correspondance entre type_maille -> nb noeud (maille linéaire)
- maille2nb={}
- # POI
- maille2nb[1] = 1
- # SEG
- maille2nb[2] = 2
- maille2nb[3] = maille2nb[4] = maille2nb[5] = maille2nb[6] = maille2nb[2]
- # TRIA
- maille2nb[7] = 3
- maille2nb[8] = maille2nb[9] = maille2nb[10] = maille2nb[11] = maille2nb[7]
- # QUAD
- maille2nb[12] = 4
- maille2nb[13] = maille2nb[14] = maille2nb[15] = maille2nb[16] = maille2nb[17] = maille2nb[12]
- # TETRA
- maille2nb[18] = 4
- maille2nb[19] = maille2nb[18]
- # PENTA
- maille2nb[20] = 5
- maille2nb[21] = maille2nb[20]
- # PYRAM
- maille2nb[22] = 5
- maille2nb[23] = maille2nb[22]
- # HEXA
- maille2nb[24] = 6
- maille2nb[25] = maille2nb[26] = maille2nb[24]
-
-
- # dictionnaire de correspondance entre type_maille -> dimension
- maille2dim = {}
- # POI
- maille2dim[1] = 0
- # SEG
- maille2dim[2] = 1
- maille2dim[3] = maille2dim[4] = maille2dim[5] = maille2dim[6] = maille2dim[2]
- # TRIA
- maille2dim[7] = 2
- maille2dim[8] = maille2dim[9] = maille2dim[10] = maille2dim[11] = maille2dim[7]
- # QUAD
- maille2dim[12] = 2
- maille2dim[13] = maille2dim[14] = maille2dim[15] = maille2dim[16] = maille2dim[17] = maille2dim[12]
- # TETRA
- maille2dim[18] = 3
- maille2dim[19] = maille2dim[18]
- # PENTA
- maille2dim[20] = 3
- maille2dim[21] = maille2dim[20]
- # PYRAM
- maille2dim[22] = 3
- maille2dim[23] = maille2dim[22]
- # HEXA
- maille2dim[24] = 3
- maille2dim[25] = maille2dim[26] = maille2dim[24]
-
- # On stocke les dictionnaires
- self.maille2aretes = maille2aretes
- self.maille2nb = maille2nb
- self.maille2dim = maille2dim
-
- return
-
-
-# ---------------------------------------------------------------------------- #
-
- def Modele_to_Liste_Mailles(self, MODELE):
-
- nommod = string.ljust(MODELE.nom,8)
- _DIC_MA = aster.getcolljev(nommod.ljust(8)+'.MODELE .LIEL')
-
- # Creation de la liste des mailles
- ll = []
- for type_maille in _DIC_MA.keys():
- ll.extend( _DIC_MA[type_maille][0:-1] )
- _LST_MA = N.array( ll ) - 1
-
- if self.OPTIONS['INFO']>=5:
- print '\n# ----- MODELE ----- #\n'
- print '_LST_MA=',len(_LST_MA),_LST_MA
- print '_DIC_MA=',len(_DIC_MA),_DIC_MA
-
- return _LST_MA
-
-
-# ---------------------------------------------------------------------------- #
-
- def Creation_Listes_Mailles_Par_Dim(self, _LST_TMA, _LST_MA=None):
-
- t0 = time.clock()
-
- # Si _LST_MA est renseigné on extrait la liste des TMA correspondante aux mailles de _LST_MA
- if _LST_MA != None:
- _LST_TMA = N.take(_LST_TMA,_LST_MA)
- else:
- _LST_MA = N.arange(len(_LST_TMA))
-
- _D_DIM_MAILLES = {}
-
- # Liste des mailles 3D (type maille de 18 à 26)
- _lst = N.where( _LST_TMA>=18, -3, _LST_TMA )
- _tmp = N.where( _lst==-3, -1, 0 )
-# _D_DIM_MAILLES[3] = N.nonzero( _tmp )
- _D_DIM_MAILLES[3] = N.take(_LST_MA, N.nonzero( _tmp ) )
-
- # Liste des mailles 2D (type maille de 7 à 17)
- _lst = N.where( _lst>=7, -2, _lst )
- _tmp = N.where( _lst==-2, -1, 0 )
- _D_DIM_MAILLES[2] = N.take(_LST_MA, N.nonzero( _tmp ) )
-
- # Liste des mailles 1D (type maille de 2 à 6)
- _lst = N.where( _lst>=2, -1, _lst )
- _tmp = N.where( _lst==-1, -1, 0 )
- _D_DIM_MAILLES[1] = N.take(_LST_MA, N.nonzero( _tmp ) )
-
- # Liste des mailles 0D (type maille 1)
- _lst = N.where( _lst>=1, -4, _lst )
- _tmp = N.where( _lst==-4, -1, 0 )
- _D_DIM_MAILLES[0] = N.take(_LST_MA, N.nonzero( _tmp ) )
-
-
- if self.OPTIONS['INFO']>=5:
- for i in _D_DIM_MAILLES.keys():
- print "-----------------"
- print 'Dim:',i, _D_DIM_MAILLES[i]
- print "-----------------"
-
- print "--- FIN Creation des listes de mailles par Dim : ", time.clock() - t0
-
- return _D_DIM_MAILLES
-
-
-# ---------------------------------------------------------------------------- #
-
- def Connectivite_Aretes(self, OPTION=None, _LST_OK=None):
-
- t0 = time.clock()
-
- # Si _LST_OK n'est pas renseigné on prend toutes les mailles
- if not _LST_OK: _LST_OK = N.arange(len(self.MAILLAGE_Python.tm))
-
- if self.OPTIONS['INFO']>=5: print '_LST_OK (ca)=',_LST_OK
-
- maille2aretes = self.maille2aretes
-
- # Creation de la :
- # - connectivite des aretes (self.MAILLAGE_Python.ca) : m1 -> [ (a1, a2), .. ]
- # - connectivite inverse des aretes (self.MAILLAGE_Python.cia) : (a1, a2) -> [ m1, m2, ... ]
-
- self.MAILLAGE_Python.ca = {}
- self.MAILLAGE_Python.cia = {}
-
- for n in _LST_OK:
-
- n1 = self.MAILLAGE_Python.tm[n]
-
- l_aretes = maille2aretes[n1] # liste des aretes de la maille n
- l_noeuds = self.MAILLAGE_Python.co[n] # liste des noeuds de la maille n
-
- for arete in l_aretes:
- ll = []
- for i in arete:
- ll.append( l_noeuds[i] )
- ll.sort()
- ll = tuple(ll)
-
- # Table de connectivité des aretes
- if OPTION:
- if not self.MAILLAGE_Python.ca.has_key(n): self.MAILLAGE_Python.ca[n]=[]
- self.MAILLAGE_Python.ca[n].append(ll)
-# try:
-# self.MAILLAGE_Python.ca[n].append(ll)
-# except KeyError:
-# self.MAILLAGE_Python.ca[n]=[ll]
-
- # Table de connectivité inverse des aretes
- if not self.MAILLAGE_Python.cia.has_key(ll): self.MAILLAGE_Python.cia[ll]=[]
- self.MAILLAGE_Python.cia[ll].append(n)
-# try:
-# self.MAILLAGE_Python.cia[ll].append(n)
-# except KeyError:
-# self.MAILLAGE_Python.cia[ll]=[n]
-
-
- if self.OPTIONS['INFO']>=5:
- for k in self.MAILLAGE_Python.cia.keys():
- print 'cia:',k, ' ', self.MAILLAGE_Python.cia[k]
- if OPTION:
- for k in self.MAILLAGE_Python.ca.keys():
- print 'ca: ',k, ' ', self.MAILLAGE_Python.ca[k]
-
-
- print "--- FIN Creation de la connectivite simple et inverse des aretes : ", time.clock() - t0
-
- return
-
-
-# ---------------------------------------------------------------------------- #
-
- def Elimination_Mailles_de_bords(self, MAILLAGE_Python, _D_DIM_MAILLES, _DIM):
- """
- Extraction des mailles de bords (mailles incluses dans un bord d une autre maille)
- """
-
- t0 = time.clock()
-
- _LST_TMA = self.MAILLAGE_Python.tm
-
- if self.OPTIONS['INFO']>=5:
- MAILLAGE = self.ASTER['MAILLAGE']
- nommail = string.ljust(MAILLAGE.nom,8)
- _LST_MAI = aster.getvectjev(nommail.ljust(8)+'.NOMMAI')
-
- # Le dico maille2nb donne le nombre de noeuds definissant un bord (lineaire)
- maille2nb = self.maille2nb
-
-
- # construction des listes des mailles de dim N-1 :
- # _LST_OK : Mailles de dim N-i qui ne sont pas un bord des mailles de dim N
- # _LST_BD : Mailles de dim N-i qui sont un bord
- #
- if self.OPTIONS['INFO']>=5: print '\n\nElimination des mailles de bord de DIM', _DIM - 1
-
- _LST4 = _D_DIM_MAILLES[ _DIM - 1 ]
- _LST_IND = N.arange( len(_LST4) ) + 1 # on ajoute 1 pour eviter le premier 0 dans les test nonzero plus bas
-
- if self.OPTIONS['INFO']>=5: print ' Mailles concernées=',_LST4
-
- i=0
- for m in _LST4:
- if self.OPTIONS['INFO']>=5: print '\n Maille de dim N-1:',m, ' Aster:',string.strip(_LST_MAI[m]), ' TMA:',self.MAILLAGE_Python.tm[m], ' CO:',self.MAILLAGE_Python.co[m], '(noeuds de cette maille)'
- nb = maille2nb[ self.MAILLAGE_Python.tm[m] ]
- ll = self.MAILLAGE_Python.co[m][0:nb]
- ll = N.sort(ll)
- ll = ll.tolist()
- ll = tuple(ll)
- if self.OPTIONS['INFO']>=5: print ' Bord (lineaire)', ll, nb
-
- try:
- if self.OPTIONS['INFO']>=5: print ' CIA=', self.MAILLAGE_Python.cia[ ll ], '(mailles de dim N qui ont cette maille pour bord)'
- _tmp=[]
- for maille in self.MAILLAGE_Python.cia[ ll ]:
- if self.OPTIONS['INFO']>=5: print ' Maille N:', maille, 'Aster:', string.strip(_LST_MAI[maille]), ' TMA:', self.MAILLAGE_Python.tm[maille]
-# self.liste_mailles_bord.append(m)
- except:
- if self.OPTIONS['INFO']>=5: print ' Maille non-bord'
- _LST_IND[i] = 0
-
- i+=1
-
- # Recuperation des mailles de bords et non-bords
- _LST_BD = N.nonzero(_LST_IND)
- _LST_BD = N.take(_LST4,_LST_BD)
-
- _LST_OK = N.where( _LST_IND==0, 1 , 0 )
- _LST_OK = N.nonzero(_LST_OK)
- _LST_OK = N.take(_LST4,_LST_OK)
-
- if self.OPTIONS['INFO']>=5: print '\nListe Maille de bords de DIM', _DIM - 1,' :',_LST_BD
- if self.OPTIONS['INFO']>=5: print 'Liste Maille de DIM', _DIM - 1,'qui ne sont pas des bords :',_LST_OK
-
- print "--- FIN Maille de bords de DIM",_DIM - 1, " : ", time.clock() - t0
- t0 = time.clock()
-
-
- # On cherche à marier les mailles de dimension N-2, N-3
- # Peut etre lent car on utilise la connectivité ! Mais pour le moment on a rien d'autre.
-
- _LST_BD0 = []
- _LST_OK0 = []
- _D_BD = {}
- for d in range(_DIM-1):
- _LST4 = _D_DIM_MAILLES[ d ]
- if self.OPTIONS['INFO']>=5: print '\n\nElimination des mailles de bord de DIM', d
- if self.OPTIONS['INFO']>=5: print ' Mailles concernées=',_LST4
- for mai in _LST4:
- if self.OPTIONS['INFO']>=5: print '\n Maille:', mai, ' Aster:',string.strip(_LST_MAI[mai]), ' TMA:',self.MAILLAGE_Python.tm[mai], ' CO:',self.MAILLAGE_Python.co[mai], '(noeuds de cette maille)'
-
- nb = maille2nb[ self.MAILLAGE_Python.tm[mai] ]
- ll = self.MAILLAGE_Python.co[mai][0:nb]
- ll = N.sort(ll)
- ll = ll.tolist()
- _tmp = tuple(ll)
-# _tmp = self.MAILLAGE_Python.co[mai]
-
- if self.OPTIONS['INFO']>=5: print ' Bord (lineaire):', _tmp, nb
-
- ok=0
- for arete in self.MAILLAGE_Python.cia:
- _nb=0
- for noe in _tmp:
- if noe in arete: _nb+=1
- if _nb == len(_tmp):
- if self.OPTIONS['INFO']>=5: print ' Maille N+i:', self.MAILLAGE_Python.cia[arete], '- Arete:', arete
- _LST_BD0.append( mai )
- ok=1
-# if not _D_BD.has_key( mai ): _D_BD[mai] = []
-# _D_BD[mai].append( self.MAILLAGE_Python.cia[arete] )
- break
- if ok == 0:
- _LST_OK0.append( mai )
-
-# print 'Mai:',mai, '_D_BD[mai]=',_D_BD[mai]
-
-
- if self.OPTIONS['INFO']>=5: print '\nListe Maille de bords de DIM', d,' :',_LST_BD0
- if self.OPTIONS['INFO']>=5: print 'Liste Maille de DIM', d,'qui ne sont pas des bords :',_LST_OK0
-
-
- print '--- FIN Maille de bords de DIM', d, ' :',time.clock() - t0
- t0 = time.clock()
-
-
- _LST_OK = N.concatenate( (_LST_OK, N.array(_LST_OK0)) )
- _LST_BD = N.concatenate( (_LST_BD, N.array(_LST_BD0)) )
-
- if self.OPTIONS['INFO']>=5: print '\nTotal:\nListe Maille de bords=',_LST_BD
- if self.OPTIONS['INFO']>=5: print 'Liste Maille non-bords=',_LST_OK,'\n'
-
-# print "--- FIN Maille de bords 3 : ", time.clock() - t0
-
- return _LST_OK, _LST_BD
-
-
-# ---------------------------------------------------------------------------- #
-
- def Affectation_Mailles_de_bords(self, _LST_BD, _DIM):
- """
- Affectation a un SD des mailles de bords (mailles incluses dans un bord d une autre maille)
- """
-
- if self.OPTIONS['INFO']>=5:
- print 'liste_mailles_bord=', self.liste_mailles_bord
- print 'liste_sd_bord', self.liste_sd_bord
- print '_LST_BD=',_LST_BD
-
-
- MAILLAGE = self.ASTER['MAILLAGE']
- _LST_TMA = self.MAILLAGE_Python.tm
-
- if self.OPTIONS['INFO']>=5:
- nommail = string.ljust(MAILLAGE.nom,8)
- _LST_MAI = aster.getvectjev(nommail.ljust(8)+'.NOMMAI')
-
- t0 = time.clock()
-
- # Affectation des mailles de bords à chacun des SD
-
- # Le dico maille2nb donne le nombre de noeuds definissant un bord (lineaire)
- maille2nb = self.maille2nb
-
- i = 0
- for m in _LST_BD:
- if self.OPTIONS['INFO']>=5: print '\n Maille de dim N-1:',m, ' Aster:',string.strip(_LST_MAI[m]), ' TMA:',self.MAILLAGE_Python.tm[m], ' CO:',self.MAILLAGE_Python.co[m], '(noeuds de cette maille)'
- nb = maille2nb[ self.MAILLAGE_Python.tm[m] ]
- ll = self.MAILLAGE_Python.co[m][0:nb]
- ll = N.sort(ll)
- ll = ll.tolist()
- ll = tuple(ll)
- if self.OPTIONS['INFO']>=5: print ' Bord (lineaire)', ll, nb
-
- # Cas particulier des POI1 en 2D et 3D (ils ne peuvent etre des bords d'elements 2D ou 3D)
- if ( (nb==1) and (_DIM>=2) ):
- _tmp=[]
- for arete in self.MAILLAGE_Python.cia.keys():
- if ll[0] in arete:
- for maille in self.MAILLAGE_Python.cia[ arete ]:
- if self.OPTIONS['INFO']>=5: print ' Maille N+i:', maille, ' Aster:',string.strip(_LST_MAI[maille]), ' Arete:', arete
- _tmp.append( self.liste_sd[maille] )
-
- # Cas particulier des SEG en 3D (ils ne peuvent etre des bords d'elements 3D)
- elif ( (nb==2) and (_DIM==3) ):
- _tmp=[]
- for arete in self.MAILLAGE_Python.cia.keys():
- _nb=0
- for noe in ll:
- if noe in arete: _nb+=1
- if _nb == len(ll):
- for maille in self.MAILLAGE_Python.cia[arete]:
- if self.OPTIONS['INFO']>=5: print ' Mailles N+i:', maille, ' Aster:',string.strip(_LST_MAI[maille]), ' Arete:', arete
- _tmp.append( self.liste_sd[maille] )
-
- # Autres mailles de bord
- else:
- if self.OPTIONS['INFO']>=5: print ' CIA=', self.MAILLAGE_Python.cia[ ll ], '(mailles de dim N qui ont cette maille pour bord)'
- _tmp=[]
- for maille in self.MAILLAGE_Python.cia[ ll ]:
- if self.OPTIONS['INFO']>=5: print ' Maille N+i:', maille, 'Aster:', string.strip(_LST_MAI[maille]), ' SD:', self.liste_sd[maille], ' TMA:', self.MAILLAGE_Python.tm[maille]
- _tmp.append( self.liste_sd[maille] )
-
- # integre la maille au SD le plus faible (pour que des groupes de bords se retrouvent dans le meme SD)
- _tmp.sort()
- self.liste_mailles_bord.append(m)
- self.liste_sd_bord.append( _tmp[0] )
- i += 1
- if self.OPTIONS['INFO']>=5: print ' ---> Maille:',m,'integree au SD:', _tmp[0]
-
- if self.OPTIONS['INFO']>=5:
- print '\n\nliste_mailles_bord=', self.liste_mailles_bord
- print 'liste_sd_bord=', self.liste_sd_bord
-
-
- print "--- FIN Affectation des mailles de bords : ", time.clock() - t0
-
- return
-
-
-# ---------------------------------------------------------------------------- #
-
- def Creation_Graphe(self):
-
- t0 = time.clock()
-
- # Creation du graphe complet
- self.GRAPH = {}
-
- for mai in self.liste_mailles:
- _ll=[]
- for are in self.MAILLAGE_Python.ca[mai]:
- _ll.extend( self.MAILLAGE_Python.cia[are] )
- _mm = enleve_doublons_liste(_ll) # coute cher!
- _tmp = _mm.tolist()
- _tmp.remove(mai)
- self.GRAPH[mai] = _tmp
-
- if self.OPTIONS['INFO']>=5: print 'self.GRAPH['+str(mai)+']=', self.GRAPH[mai]
-
- print "--- FIN Creation du graphe complet : ", time.clock() - t0
-
- return
-
-
-# ---------------------------------------------------------------------------- #
-
- def Reduction_Graphe(self, _DIM):
-
- t0 = time.clock()
-
- # Elimination des connectivités à interface nulle
- maille2dim = self.maille2dim
- _lst2 = []
- for mai in self.liste_mailles:
- if self.OPTIONS['INFO']>=5: print '\nmai:', mai, 'co:', self.MAILLAGE_Python.co[mai], 'tm:', self.MAILLAGE_Python.tm[mai]
- _DIM1 = maille2dim[ self.MAILLAGE_Python.tm[mai] ]
- _tmp2 =[]
- for mai2 in self.GRAPH[mai]:
- if self.OPTIONS['INFO']>=5: print 'mai2:', mai2, 'co:', self.MAILLAGE_Python.co[mai2], 'tm:', self.MAILLAGE_Python.tm[mai2]
- # calcule le nombre de noeuds communs aux deux mailles
- _nb = 0
- for noe in self.MAILLAGE_Python.co[mai2]:
- if noe in self.MAILLAGE_Python.co[mai]: _nb += 1
- _DIM2 = maille2dim[ self.MAILLAGE_Python.tm[mai2] ]
- if _nb >= min(_DIM1, _DIM2): # le min permet de faire du collage 3D-coque par exemple
- _tmp2.append( mai2 )
- _tmp = [mai, mai2]
- _tmp.sort()
- _lst2.append(_tmp)
- self.GRAPH[mai] = _tmp2
-
- print "--- FIN Elimination des connectivités avec une interface nulle : ", time.clock() - t0
- t0 = time.clock()
-
-
- # Calcul du nombre d'aretes
- # A voir : normalement il n'y a rien a faire car nb0 = 2*nb (a verifier...)
- _lst2.sort()
- _v = _lst2[0]
- _nb = 1
- for i in _lst2:
- if i != _v:
- _v = i
- _nb += 1
-
-
- if self.OPTIONS['INFO']>=5:
- print '----------------------------------------------'
- for mai in self.liste_mailles:
- print 'self.GRAPH['+str(mai)+']=', self.GRAPH[mai]
- print '----------------------------------------------'
-
- return _nb
-
-
-# ------------------------------------------------------------------ #
-
- def Ecrire_Graphe(self, f_metis, _nb):
-
- t0 = time.clock()
-
- # On doit renumeroter les mailles qui arrivent dans self.liste_mailles pour avoir 0... N-1
- _D_CORRES = {}
- for i in N.arange(len(self.liste_mailles)):
- _D_CORRES[ self.liste_mailles[i] ] = i
-
- # Ecriture du fichier fort.UL pour metis
- fw = open(f_metis,'w')
- fw.write( str(len(self.liste_mailles)) + ' ' + str(_nb) + '\n')
- for l in self.liste_mailles:
-# try:
- _tmp = []
- for t in self.GRAPH[l]:
- try:
- t = _D_CORRES[t]
- _tmp.append( str(t+1) ) # Necessaire car metis numerote de 1 à N
- except:
- print 'on oublie le bord:', t
- fw.write( string.join(_tmp, ' ') + '\n' )
-# except:
-# print 'Probleme ecriture graphe! On continue..'
- fw.close()
-
- print "--- FIN Ecriture du fichier du graphe pour metis : ", time.clock() - t0
-
- return _D_CORRES
-
-
-# ---------------------------------------------------------------------------- #
-
- def Lecture_fichier_sdd(self, fichier, _LST_OK):
-
- t0 = time.clock()
-
- # Lecture du fichier produit par metis (partie a optimiser)
- try:
- f = open( fichier, 'r' )
- except:
- print "\n\n ERREUR: le fichier est introuvable! Le partitionneur \n ne s'est probablement pas lancé.\n\n"
- sys.exit(1)
- else:
- _tmp = []
- for l in f.readlines():
- _tmp.append( int(string.strip(l)) )
- f.close()
- _l_domaines = N.array(_tmp,copy=0)
-
- # Pour garder le fichier metis
- os.system( 'mv ' + fichier + ' REPE_OUT/' )
-
- if self.OPTIONS['INFO']>=5: print '_l_domaines=',_l_domaines
-
- print "--- FIN Lecture du fichier produit par metis : ", time.clock() - t0
-
- return _l_domaines
-
-
-# ---------------------------------------------------------------------------- #
-
- def Creation_Group_ma_Python_par_SD(self, NOM='SD', NOM2='B'):
-
- t0 = time.clock()
-
- NB_PART = self.OPTIONS['NB_PART']
-
- # Creation du dictionnaire des listes des mailles par SD
- # d_gma : { num sd -> [ liste mailles ] }
- d_gma = {}
- for i in range(NB_PART):
- d_gma[i] = []
-
- i=0
- for sdd in self.liste_sd:
- d_gma[sdd].append( self.liste_mailles[i] )
- i+=1
-
-
- # Creation du dictionnaire des listes des mailles de bord par SD
- # d_gma_bord : { num sd -> [ liste mailles ] }
- d_gma_bord = {}
- for i in range(NB_PART):
- d_gma_bord[i] = []
-
- i=0
- for sdd in self.liste_sd_bord:
- d_gma_bord[sdd].append( self.liste_mailles_bord[i] )
- i+=1
-
-
- # Generation des listes de noms de groupes
- _l_sd = []
- _l_bord = []
- for i in range(NB_PART):
- if d_gma[i] != []:
- _l_sd.append( NOM + str(i) )
- if d_gma_bord[i] != []:
- _l_bord.append( NOM2 + str(i) )
-
- # Stockage
- self.ASTER['GROUP_MA'] = _l_sd
- self.ASTER['GROUP_MA_BORD'] = _l_bord
-
-
- # Creation des groupes de mailles dans le Maillage Python
- for i in range(NB_PART):
- self.MAILLAGE_Python.gma[NOM+str(i)] = d_gma[i]
- self.MAILLAGE_Python.gma[NOM2+str(i)] = d_gma_bord[i]
-
- print "--- FIN creation du dictionnaire des listes des mailles par SD ", time.clock() - t0
-
- return
-
-
-# ---------------------------------------------------------------------------- #
-
- def Creation_Group_ma_Aster_par_SD(self, NOM='SD', NOM2='B', INCLUSE='NON'):
-
- t0 = time.clock()
-
- MAILLAGE = self.ASTER['MAILLAGE']
- NB_PART = self.OPTIONS['NB_PART']
-
- nommail = string.ljust(MAILLAGE.nom,8)
- _LST_MAI = aster.getvectjev(nommail.ljust(8)+'.NOMMAI')
-
-
- # Creation du dictionnaire des listes des mailles par SD
- # d_gma : { num sd -> [ liste mailles ] }
- d_gma = {}
- for i in range(NB_PART):
- d_gma[i] = []
-
- m=0
- for sdd in self.liste_sd:
- d_gma[sdd].append( string.strip(_LST_MAI[ self.liste_mailles[m] ]) ) # voir si le strip coute cher !
- m += 1
-
-
- # Creation du dictionnaire des listes des mailles de bord par SD
- # d_gma_bord : { num sd -> [ liste mailles ] }
- d_gma_bord = {}
- for i in range(NB_PART):
- d_gma_bord[i] = []
-
- # On inclus directement les mailles de bords dans les SD
- if INCLUSE=='OUI':
- m=0
- for sdd in self.liste_sd_bord:
- d_gma[sdd].append( string.strip(_LST_MAI[ self.liste_mailles_bord[m] ]) ) # voir si le strip coute cher !
- m+=1
-
- else:
- m=0
- for sdd in self.liste_sd_bord:
- d_gma_bord[sdd].append( string.strip(_LST_MAI[ self.liste_mailles_bord[m] ]) ) # voir si le strip coute cher !
- m+=1
-
-
- print "--- FIN creation du dictionnaire des listes des mailles par SD ", time.clock() - t0
- t0 = time.clock()
-
-
- # Creation et lancement de la commande DEFI_GROUP associée
- try:
- DEFI_GROUP = self.jdc.get_cmd('DEFI_GROUP')
- except:
- try:
- from Cata.cata import DEFI_GROUP
- except:
- print "\n\nERREUR : il faut lancer ce programme depuis Aster pour pouvoir \ngénérer les groupes de mailles Aster.\n\n"
- return
-
- _tmp = []
- _l_sd = []
- _l_bord = []
- for i in range(NB_PART):
- if d_gma[i] != []:
- _tmp.append( {'MAILLE': d_gma[i],'NOM': NOM + str(i)} )
- _l_sd.append( NOM + str(i) )
- if d_gma_bord[i] != []:
- _tmp.append( {'MAILLE': d_gma_bord[i],'NOM': NOM2 + str(i)} )
- _l_bord.append( NOM2 + str(i) )
-
- motscle2= {'CREA_GROUP_MA': _tmp }
-
- DEFI_GROUP( reuse=MAILLAGE,
- MAILLAGE=MAILLAGE,
- INFO=1,
- **motscle2
- ) ;
-
- # Stockage
- self.ASTER['DICO_SD_MAILLES'] = d_gma
- self.ASTER['GROUP_MA'] = _l_sd
- self.ASTER['GROUP_MA_BORD'] = _l_bord
-
- print "--- FIN Creation et lancement de la commande DEFI_GROUP associée : ", time.clock() - t0
-
- return
-
-# ---------------------------------------------------------------------------- #
+++ /dev/null
-#@ MODIF sup_gmsh Utilitai DATE 10/05/2005 AUTEUR GJBHHEL E.LORENTZ
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-import os.path, string, os, copy
-from Numeric import *
-
-try :
- from Cata.cata import *
- from Accas import _F
-except :
- print 'Fonctionnalites Aster indisponibles'
-
-
-_CARAC = {
- 'prec' : 1.E-8
- }
-
-
-def Class_search(class_name, target_class) :
-
- """
- Check that class_name inherits from target_class
- (run recursively through the inheritance lists)
- """
-
- if class_name == target_class : return 1
-
- for cl in class_name.__bases__ :
- if Class_search(cl, target_class) : return 1
- return 0
-
-
-def Progress(L,**para) :
-
- """
- Compute the unknown parameters for a geometric progression :
- r = ratio
- N = number of elements
- h = initial size
-
- So that :
- L = Sum(i=0:N-1, h(i)) where h(i+1) = h(i)*r, h(0)=h
-
- Usage :
- Progress(L,r=...,N=...) -> h
- Progress(L,r=...,h=...) -> N
- Progress(L,h=...,N=...) -> r
-
- """
-
- prec = 1.E-4
-
- # Calcul de N
- if 'N' not in para.keys() :
- r = float(para['r'])
- h = float(para['h'])
- N = log(1+(r-1)*L/h)/log(r)
- N = int(N+0.5)
- return N
-
- # Calcul de h
- elif 'h' not in para.keys() :
- r = float(para['r'])
- N = int(para['N']+0.5)
- h = L*(r-1)/(r**N-1)
- return h
-
- # Calcul de r
- elif 'r' not in para.keys() :
- h = float(para['h'])
- N = int(para['N']+0.5)
- a = L/h
- if N > a :
- x = 0
- else :
- x = a**(1./(N-1))
-
- for i in xrange(100) :
- res = x**N - a*x + a-1
- if abs(res) < prec*(x-1)*a :
- return x
- dx = res/(a-N*x**(N-1))
- x = x+dx
-
- raise 'Solution failure'
-
- else :
- raise 'Unknown parameters'
-
-
-
-def Is_Geometric(object) :
-
- """
- return true if the object inherits of the Geometric class
- """
-
- return Class_search(object.__class__, Geometric)
-
-
-
-# --------------------------------------------------------------------------------------
-
-class Geometric :
-
- """
- GENERIC CLASS FOR GEOMETRICAL OBJECTS
-
- private attribute
- parameters : dictionnary of the attributes (except relation and parameters itself)
- see __getattr__ and __setattr__
-
-
- Attributes
- num : index among gmsh objects
- md : mesh descriptor
- mesh : related mesh object
- relation : model object in case of coincidence
- type : type of the object (None, Point, Line, Circle, Surface, Volume)
-
- Public methods
- Is_point : return true is the object inherits of the Point class
-
- Is_line : return true is the object inherits of the Line class
-
- Is_surface : return true is the object inherits of the Surface class
-
- Is_volume : return true is the object inherits of the Volume class
-
- Base_class : return the name of the base class : Point, Line, Surface
- or Volume
-
- Is_same_dimension : return true is both objects are of the same dimension
- (point, line, surface or volume)
- in -> object to compare to self
-
- Duplicate : duplicate an object and base its mesh_descriptor
- on the mesh_descriptor of the model
-
- Coincide : assert that an object is coincident with a model one
- All the attributes are then automatically read from
- the model object (see __setattr__ and __getattr__).
- in -> model object
-
- Private method
-
- Root :
- Provides the root object of an object, ie the object itself if there is no relation
- or the deepest model in case of relation.
-
- Geometric_coincide : check if a geometrical coincidence is possible
- return information about the coincidence, false else.
- in -> model object
-
- Deep_coincide : proceed recursively to ensure coincidence of the relevant sub-objects
- in -> model object
- in -> correspond (information returned by Geometric_coincide)
-
- __setattr__ : distinguish two sets of attributes
- relation (to express a relation with a model object in case of coincidence)
- all the other attributes which are stored in the dictionnary parameters
- instead of the usual __dict__ if there is no relation (see Coincide)
- and in the model object if there is a coincidence
-
- __getattr__ : if the object is related (relation <> None) the attribute is read
- in the model object. Else, it is read in the current object, actually
- in the dictionnary parameters (see __setattr__)
-
- Thanks to these two overloaded methods, the access to the attributes is usual if
- there is no relation whereas the attributes of the model object are accessed
- transparently if there is a relation.
-
- __cmp__ :
- The comparison of two objects involves possible coincidence. It is no more the object ids
- that are compared but the object roots (.relation if any).
-
- Gmsh : produce the source code for Gmsh
- in -> mesh
-
- Gmsh_send : send a line code to the gmsh interpreter
- in -> line_code (string)
-
- Intermediate_meshing : produce the source code for the intermediate objects
- in -> mesh
-
- Object meshing : produce the source code for the current object
- var -> object number (modified if several objects are created)
-
- """
-
- def __init__(self) :
- self.relation = None
- self.parameters = {}
- self.num = 0
- self.md = Mesh_Descriptor()
-
- types = {}
- types[Geometric] = 'Geometric'
- types[Point] = 'Point'
- types[Line] = 'Line'
- types[Circle] = 'Circle'
- types[Surface] = 'Surface'
- types[Volume] = 'Volume'
- types[LineLoop] = 'LineLoop'
- types[SurfaceLoop] = 'SurfaceLoop'
- try :
- self.type = types[self.__class__]
- except KeyError :
- raise 'Unknown object type'
-
-
- def Is_point(self) :
- return Class_search(self.__class__, Point)
-
- def Is_line(self) :
- return Class_search(self.__class__, Line)
-
- def Is_surface(self) :
- return Class_search(self.__class__, Surface)
-
- def Is_volume(self) :
- return Class_search(self.__class__, Volume)
-
- def Base_class(self) :
- if self.Is_volume() : return 'Volume'
- if self.Is_surface() : return 'Surface'
- if self.Is_line() : return 'Line'
- if self.Is_point() : return 'Point'
-
- def Is_same_dimension(self, obj) :
-
- return (
- (self.Is_point() and obj.Is_point() ) or
- (self.Is_line() and obj.Is_line() ) or
- (self.Is_surface() and obj.Is_surface() ) or
- (self.Is_volume() and obj.Is_volume() )
- )
-
-
- def __setattr__(self, attr, value) :
-
- if attr in ['relation','parameters'] :
- self.__dict__[attr] = value
- else :
- if self.relation :
- setattr(self.relation,attr,value)
- else :
- self.parameters[attr] = value
-
-
- def __getattr__(self,attr) :
-
- if self.relation :
- return (getattr(self.relation,attr))
- else :
- if attr in self.parameters.keys() :
- return self.parameters[attr]
- else :
- raise AttributeError,attr
-
-
- def Root(self) :
-
- o = self
- while o.relation : o = o.relation
- return o
-
-
- def __cmp__(self,obj) :
-
- if self.Root() is obj.Root() :
- return 0
- else :
- return -1
-
-
- def Geometric_coincide(self,obj) : return 0
-
- def Deep_coincide(self,obj,correspond) : pass
-
- def Coincide(self, obj) :
-
- if self == obj : return # in that way recursive loops cannot exist
-
- if self.relation : # the root is put in coincidence, not the object itself
- self.Root().Coincide(obj)
- return
-
- if not self.Is_same_dimension(obj) :
- raise 'Coincidence impossible : objects are not of the same dimension'
-
- correspond = self.Geometric_coincide(obj)
- if not correspond :
- raise 'The objects are not geometrically coincident'
-
- self.Deep_coincide(obj,correspond)
- self.relation = obj
-
-
- def Duplicate(self) :
-
- return copy.deepcopy(self) # special deepcopy for the Mesh_Descriptor
-
-
- def Gmsh(self,mesh) :
-
- if self.num : return # already meshed object
- self.mesh = mesh # Storing the mesh
- self.Intermediate_meshing(mesh) # creation of the intermediate objects
- num = mesh.num+1 # New object number
- self.Object_meshing(num) # object meshing (with current number num)
- mesh.num = num # Updating the current gmsh pointer
- self.num = num # Storing the current object number
-
-
- def Gmsh_send(self, line_code) :
-
- self.mesh.command.append(line_code)
-
-
- def Intermediate_meshing(self,mesh) :
- pass
-
-
- def Object_meshing(self,num) :
- raise "Creation of the Gmsh source not implemented"
-
-
-
-
-# -------------------- POINT OBJECTS ---------------------
-
-
-class Point(Geometric) :
-
- """
- POINT OBJECT
-
- Public methods
- __init__ :
- in -> coordinates (the 3rd is zero by default)
-
- Size : set the size of the neighbouring elements
- in -> size
-
- Attractor : define the point as an attractor
- in -> scale_x : size amplification factor in the x-direction
- in -> scale_y : size amplification factor in the y-direction
- in -> distance: influence distance for the perturbation
-
- Translate : translation of the point
- in -> x,y,z translation vector (default : z=0)
-
- Attributes
- coor : coordinates
- size : neighbouring element size
- attractor : parameters of the attractor
- """
-
- def __init__(self,x,y,z=0) :
-
- Geometric.__init__(self)
- self.coor = array([x,y,z], Float)
- self.attractor = None
-
-
- def Geometric_coincide(self,obj) :
-
- global _CARAC
- prec = _CARAC['prec']
-
- d = VectorNorm(self.coor - obj.coor)
- if d < prec*self.md.size :
- return 1
- else :
- return None
-
-
- def Size(self,h) :
-
- self.md.size = float(h)
-
-
- def Attractor(self, scale_x, scale_y, distance) :
-
- self.attractor = (float(scale_x), float(scale_y), float(distance))
-
-
- def Translate(self,x,y,z=0) :
-
- tran = array([x,y,z]).astype(Float)
- self.coor = self.coor + tran
-
-
- def Object_meshing(self,num) :
-
- ch = (
- 'Point(' + `num` + ') = {'
- + `self.coor[0]` + ', '
- + `self.coor[1]` + ', '
- + `self.coor[2]` + ', '
- + `self.md.size` + '};'
- )
- self.Gmsh_send(ch)
-
- if self.attractor :
- ch = (
- 'Attractor Point{' + `num` + '} = {'
- + `self.attractor[0]`+','
- + `self.attractor[1]`+','
- + `1./self.attractor[2]` + '};'
- )
- self.Gmsh_send(ch)
-
-
-
-# -------------------- LINE OBJECTS ----------------------
-
-
-class Line(Geometric) :
-
- """
- LINE OBJECT
-
-
- Public methods
-
- Attractor : define the point as an attractor
- in -> scale_x : size amplification factor in the x-direction
- in -> scale_y : size amplification factor in the y-direction
- in -> distance: influence distance for the perturbation
-
- """
-
-
- def __init__(self,*points) :
-
- Geometric.__init__(self)
-
- if len(points) <=1 :
- raise "There should be at least two points"
-
- for point in points :
- if not point.Is_point() :
- raise "Arguments should be points"
-
- self.points = list(points)
- self.attractor = None
-
-
- def Geometric_coincide(self,obj) :
-
- nb_points = len(self.points)
- if nb_points <> len(obj.points) :
- raise 'To coincide, lines should have the same number of points'
-
- # same order of points
- info = range(nb_points)
- for i in range(nb_points) :
- p1 = self.points[i]
- p2 = obj.points[info[i]]
- if not p1.Geometric_coincide(p2) :
- break
- else :
- return info
-
- # reverse order of points
- info.reverse()
- for i in range(nb_points) :
- p1 = self.points[i]
- p2 = obj.points[info[i]]
- if not p1.Geometric_coincide(p2) :
- break
- else :
- return info
-
- return None
-
-
- def Deep_coincide(self,obj,info) :
-
- for i in range(len(info)) :
- p1 = self.points[i]
- p2 = obj.points[info[i]]
- p1.Coincide(p2)
-
-
-
- def Translate(self,x,y,z=0) :
-
- for point in self.points :
- point.Translate(x,y,z)
-
-
- def Transfinite(self,number,progression = 1) :
-
- self.md.number = int(number)
- self.md.progression = float(progression)
-
-
- def Attractor(self,scale_x, scale_y, distance) :
-
- self.attractor = (float(scale_x), float(scale_y), float(distance))
-
-
- def __rmul__(self,base) :
-
- if len(self.points) > 2 :
- raise "Support (right argument) should be a straight line"
-
- if self.points[0] in base.points :
- supp_orig = 0
- supp_extr = 1
- elif self.points[1] in base.points :
- supp_orig = 1
- supp_extr = 0
- else :
- raise "No common point"
-
- if self.points[supp_orig] == base.points[0] :
- base_orig = 0
- base_extr = -1
- else :
- base_orig = -1
- base_extr = 0
-
- # Translation vector
- ce = self.points[supp_extr].coor
- co = self.points[supp_orig].coor
- tran = ce-co
-
- # Definition of the edge opposite to the base
- opp_base = base.Duplicate()
- opp_base.Translate(tran[0],tran[1],tran[2])
- opp_base.points[base_orig] = self.points[supp_extr]
-
- # Definition of the edge opposite to the support
- opp_supp = self.Duplicate()
- opp_supp.points[0] = base.points[base_extr]
- opp_supp.points[1] = opp_base.points[base_extr]
-
- surf = Surface(base,self,opp_base,opp_supp)
-
- if len(base.points) > 2 : surf.Ruled()
-
- return surf
-
-
- def Intermediate_meshing(self,mesh) :
-
- for point in self.points :
- point.Gmsh(mesh)
-
-
- def Object_meshing(self,num) :
-
- ch = self.type + '(' + `num` + ') = {'
- for point in self.points :
- ch = ch + `point.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
- if self.md.transfinite :
- ch = (
- 'Transfinite Line{' + `num` + '} = ' +
- `self.md.number+1` +
- ' Using Progression ' + `self.md.progression` + ';'
- )
- self.Gmsh_send(ch)
-
- if self.attractor :
- ch = (
- 'Attractor Line{' + `num` + '} = {'
- + `self.attractor[0]`+','
- + `self.attractor[1]`+','
- + `1./self.attractor[2]` + '};'
- )
- self.Gmsh_send(ch)
-
-
-
-class Circle(Line) : pass
-
-# The class inherits everything from Line but its name (to tell Gmsh
-# that it is a circle).
-
-
-
-def Curve(l_x,l_y,l_z=None) :
-
- if not l_z :
- l_z = [0.] * len(l_x)
-
- l_P = []
- for x,y,z in map(None,l_x,l_y,l_z) :
- l_P.append(Point(x,y,z))
-
- line = apply(Line,l_P)
- return line
-
-
-
-
-
-# -------------------- SURFACE OBJECTS ---------------------
-
-
-class Surface(Geometric) :
-
- """
- SURFACE OBJECT (inherit from the Geometric class)
-
- Public methods
- __init__ :
- in -> lines : external bounday of the surface (lines should be connected)
-
- Holes : set the internal holes (surfaces)
- in -> holes : list of holes
-
- Boundary : checks that the boundary is a closed loop and returns the orientation of the edges
-
- Summit : returns the summit list sorted according to the orientation (see Boundary method)
-
- Ruled : declare the surface is a ruled one
-
- Translate : translate the surface
- in -> x,y,z translation vector (default : z=0)
-
- Recombine : recombine the surface (try to mesh with quadrangles instead of triangles)
-
- Transfinite : Declare the mesh to be transfinite
-
-
- Attributes
- lines : list of external boundary lines
- holes : list of internal holes (surfaces)
- ruled : indicates (false or true) if the surface is a ruled surface
- loops : list of boundary (external and internal) loops (computed when meshing)
- """
-
-
- def __init__(self,*lines) :
-
- Geometric.__init__(self)
- self.lines = list(lines)
- self.holes = []
- self.ruled = 0
-
- # Check Assumptions
- for line in lines :
- if not line.Is_line() :
- raise "Arguments should be lines : " + repr(line)
- if lines == 0 : raise "There should be at least one line"
- self.Boundary()
-
-
- def Boundary(self) :
-
- # checking the boundary is a loop
- orie = []
- tmp = list(self.lines) + [self.lines[0]]
- for i in range(len(self.lines)) :
- lb = tmp[i]
- la = tmp[i+1]
- if lb.points[-1] in [la.points[0], la.points[-1]] :
- orie.append(1)
- elif lb.points[0] in [la.points[0], la.points[-1]] :
- orie.append(-1)
- else :
- raise "This is not a loop"
-
- # checking the boundary is closed
- if orie[0] == 1 : pi = self.lines[0].points[0]
- if orie[0] == -1 : pi = self.lines[0].points[-1]
- if orie[-1] == 1 : pf = self.lines[-1].points[-1]
- if orie[-1] == -1 : pf = self.lines[-1].points[0]
- if pi <> pf : raise "The loop is not closed"
-
- return orie
-
-
- def Summit(self) :
-
- summits = []
- for line, orie in map(None,self.lines,self.Boundary()) :
- if orie == 1 :
- summits.append(line.points[0])
- else :
- summits.append(line.points[-1])
- return summits
-
-
- def Holes(self,*holes) :
-
- for hole in holes :
- if not hole.Is_surface() :
- raise "Holes should be surfaces"
- self.holes = list(holes)
-
-
- def Geometric_coincide(self,obj) :
-
- """
- return (line_order, hole_order) :
- line_order : list of the coupled lines ith line of self with line_order[i]th line of obj
- hole_order : same as line_order but with the internal holes
- """
-
- if len(self.lines) <> len(obj.lines) :
- raise 'To coincide, surfaces should have the same number of border lines'
-
- if len(self.holes) <> len(obj.holes) :
- raise 'To coincide, surfaces should have the same number of internal holes'
-
- # Coincidence of the surface holes
- hole_order = []
- nb_holes = len(self.holes)
- for hole_1 in self.holes :
- for i in xrange(nb_holes) :
- if i in hole_order :
- continue
- hole_2 = obj.holes[i]
- if hole_1.Geometric_coincide(hole_2) :
- hole_order.append(i)
- break
- else :
- return None
-
- # Coincidence of the external boundary lines
- line_order = []
- nb_lines = len(self.lines)
- for line_1 in self.lines :
- for i in xrange(nb_lines) :
- if i in line_order :
- continue
- line_2 = obj.lines[i]
- if line_1.Geometric_coincide(line_2) :
- line_order.append(i)
- break
- else :
- return None
-
- return (line_order, hole_order)
-
-
- def Deep_coincide(self,obj,info) :
-
- line_order = info[0]
- hole_order = info[1]
-
- for i,j in map(None,xrange(len(line_order)),line_order) :
- l1 = self.lines[i]
- l2 = obj.lines[j]
- l1.Coincide(l2)
-
- for i,j in map(None,xrange(len(hole_order)),hole_order) :
- h1 = self.holes[i]
- h2 = obj.holes[j]
- h1.Coincide(h2)
-
-
-
- def Ruled(self) :
-
- self.ruled = 1
-
- if len(self.lines) not in [3,4] :
- raise "Ruled surfaces require 3 or 4 edges"
-
- if self.holes :
- raise "Holes are forbidden for ruled surfaces"
-
-
- def Translate(self,x,y,z=0) :
-
- l_points = []
- for surf in [self] + self.holes :
- for line in surf.lines :
- for point in line.points :
- if point not in l_points : l_points.append(point)
-
- for point in l_points :
- point.Translate(x,y,z)
-
-
- def Recombine(self,val=1) :
-
- self.md.recombine = val
-
-
- def Transfinite(self) :
-
- self.Ruled()
-
- if len(self.lines) == 4 :
- self.Recombine()
-
- self.md.transfinite = 1
-
- for line in self.lines :
- if not line.md.transfinite :
- raise "Transfinite surfaces require transfinite edges"
-
- if (
- self.lines[0].md.number <> self.lines[2].md.number or
- self.lines[1].md.number <> self.lines[3].md.number
- ) :
- raise "Coupled edges should have the same number of elements"
-
-
-
- def Intermediate_meshing(self,mesh) :
-
- self.loops = []
- for surf in [self]+self.holes :
- loop = LineLoop(surf)
- self.loops.append(loop)
- loop.Gmsh(mesh)
-
-
- def Object_meshing(self,num) :
-
- # Creation of the surface
- if self.ruled :
- ch = 'Ruled Surface(' + `num` + ') = {'
- else :
- ch = 'Plane Surface(' + `num` + ') = {'
- for loop in self.loops :
- ch = ch + `loop.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
- # Declaration of transfinite surface
- if self.md.transfinite :
- ch = 'Transfinite Surface {' + `num` + '} = {'
- for summit in self.Summit() :
- ch = ch + `summit.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
- # Recombine elements if requested
- if self.md.recombine :
- self.Gmsh_send('Recombine Surface {' + `num` + '} ;')
-
-
-
-class LineLoop(Geometric) : # Used only during the meshing phase
-
-
- def __init__(self,surface) :
-
- Geometric.__init__(self)
- self.surface = surface
-
-
- def Intermediate_meshing(self,mesh) :
-
- for line in self.surface.lines :
- line.Gmsh(mesh)
-
-
- def Object_meshing(self,num) :
-
- ch = 'Line Loop(' + `num` + ') = {'
- for line,orie in map(None,self.surface.lines,self.surface.Boundary()) :
- ch = ch + `orie*line.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
-
-
-
-class Volume(Geometric) :
-
-
- """
- VOLUME OBJECT (inherit from the Geometric class)
-
- Public methods
- __init__ :
- in -> surfaces : external bounday of the volume (surfaces should be connected)
-
- Edge : returns the list of edges
-
- Holes : set the internal holes (surfaces)
- in -> holes : list of holes
-
- Transfinite : Declare the mesh to be transfinite (force the surfaces to be transfinite too)
-
- Translate : translate the surface
- in -> x,y,z translation vector (default : z=0)
-
-
- Attributes
- surfaces : list of external boundary surfaces
- holes : list of internal holes (volumes)
- loops : list of boundary (external and internal) loops (computed when meshing)
-
-
- Private methods :
-
- Boundary : checks that the boundary is a closed loop and returns the orientation of the edges
-
- """
-
-
- def __init__(self,*surfaces) :
-
- Geometric.__init__(self)
- self.surfaces = list(surfaces)
- self.holes = []
- self.ruled = 0
-
- # Check Assumptions
- for surface in surfaces :
- if not surface.Is_surface() :
- raise "Arguments should be surfaces : " + repr(surface)
- if len(surfaces) < 2 : raise "There should be at least two surfaces"
- self.Boundary()
-
-
- def Boundary(self) :
-
- edges = []
- for surface in self.surfaces :
- edges = edges + surface.lines
-
- # each edge has to appear twice in the list of edges
- for edge in edges :
- if edges.count(edge) <> 2 :
- raise "The surface loop is not closed : each edge should appear twice"
-
-
- def Edge(self) :
-
- edges = []
- for surface in self.surfaces :
- for line in surface.lines :
- if line not in edges : edges.append(line)
- return edges
-
-
- def Holes(self,*holes) :
-
- for hole in holes :
- if not hole.Is_volume() :
- raise "Holes should be volumes"
- self.holes = list(holes)
-
-
- def Geometric_coincide(self,obj) :
-
- """
- return (surface_order, hole_order) :
- surface_order : list of the coupled surfaces ith surface of self with surface_order[i]th surface of obj
- hole_order : same as surface_order but with the internal holes
- """
-
- if len(self.surfaces) <> len(obj.surfaces) :
- raise 'To coincide, volumes should have the same number of border surfaces'
-
- if len(self.holes) <> len(obj.holes) :
- raise 'To coincide, volumes should have the same number of internal holes'
-
- # Coincidence of the surface holes
- hole_order = []
- nb_holes = len(self.holes)
- for hole_1 in self.holes :
- for i in xrange(nb_holes) :
- if i in hole_order :
- continue
- hole_2 = obj.holes[i]
- if hole_1.Geometric_coincide(hole_2) :
- hole_order.append(i)
- break
- else :
- return None
-
- # Coincidence of the external boundary lines
- surface_order = []
- nb_surfaces = len(self.surfaces)
- for surface_1 in self.surfaces :
- for i in xrange(nb_surfaces) :
- if i in surface_order :
- continue
- line_2 = obj.surfaces[i]
- if surface_1.Geometric_coincide(surface_2) :
- surface_order.append(i)
- break
- else :
- return None
-
- return (surface_order, hole_order)
-
-
- def Deep_coincide(self,obj,info) :
-
- surface_order = info[0]
- hole_order = info[1]
-
- for i,j in map(None,xrange(len(surface_order)),surface_order) :
- l1 = self.surfaces[i]
- l2 = obj.surfaces[j]
- l1.Coincide(l2)
-
- for i,j in map(None,xrange(len(hole_order)),hole_order) :
- h1 = self.holes[i]
- h2 = obj.holes[j]
- h1.Coincide(h2)
-
-
- def Transfinite(self) :
-
- if len(self.surfaces) == 5 :
- raise "Not implemented"
-
- if len(self.surfaces) not in [5,6] :
- raise "Transfinite volumes require 5 or 6 faces"
-
- if self.holes :
- raise "Holes are forbidden for transfinite volumes"
-
- self.md.transfinite = 1
-
- for surface in self.surfaces :
- if not surface.md.transfinite :
- surface.Transfinite() # attention : ce n'est pas vrai dans le cas des prismes
-
-# ATTENTION : ICI, IL FAUDRAIT VERIFIER QUE LES SURFACES PEUVENT ETRE MISES EN VIS A VIS
-
-
-
- def Translate(self,x,y,z=0) :
-
- l_points = []
- for volu in [self] + self.holes :
- for surf in volu.surfaces :
- for line in surf.lines :
- for point in line.points :
- if point not in l_points : l_points.append(point)
-
- for point in l_points :
- point.Translate(x,y,z)
-
-
-
- def Intermediate_meshing(self,mesh) :
-
- self.loops = []
- for volume in [self]+self.holes :
- loop = SurfaceLoop(volume)
- self.loops.append(loop)
- loop.Gmsh(mesh)
-
-
- def Object_meshing(self,num) :
-
- # Creation of the volume
- ch = 'Volume(' + `num` + ') = {'
- for loop in self.loops :
- ch = ch + `loop.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
- # Declaration of transfinite surface
- if self.md.transfinite :
-
- bottom_summits = self.surfaces[0].Summit()
- edges = self.Edge()
- top_summits = []
- for summit in bottom_summits :
- for edge in edges :
- if summit == edge.points[0] and edge.points[-1] not in bottom_summits :
- top_summits.append(edge.points[-1])
- break
- elif summit == edge.points[-1] and edge.points[0] not in bottom_summits :
- top_summits.append(edge.points[0])
- break
-
- ch = 'Transfinite Volume {' + `num` + '} = {'
- for summit in bottom_summits + top_summits :
- ch = ch + `summit.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
-
-
-class SurfaceLoop(Geometric) : # Used only during the meshing phase
-
-
- def __init__(self,volume) :
-
- Geometric.__init__(self)
- self.volume = volume
-
-
- def Intermediate_meshing(self,mesh) :
-
- for surface in self.volume.surfaces :
- surface.Gmsh(mesh)
-
-
- def Object_meshing(self,num) :
-
- ch = 'Surface Loop(' + `num` + ') = {'
- for surface in self.volume.surfaces :
- ch = ch + `surface.num` + ','
- ch = ch[:-1] + '};'
- self.Gmsh_send(ch)
-
-
-
-
-# ------------------- GEOMETRICAL TRANSFORMATION --------------
-
-def VectorProduct(u,v) :
-
- return array([u[1]*v[2]-u[2]*v[1],u[2]*v[0]-u[0]*v[2],u[0]*v[1]-u[1]*v[0]])
-
-
-def VectorNorm(u) :
-
- return sqrt(dot(u,u))
-
-
-class Rotation :
-
- def __init__(self,A,C,B) :
-
- self.C = c
- self.a = A-C
- n = VectorProduct(self.a,B-C)
- self.n = n / VectorNorm(n)
-
-
- def Proj(self,M) :
-
- lbd = dot(M-self.C,self.n)
- H = self.C + lbd*self.n
- return H
-
-
-def Scaling_P2(p,t) : return (1.-p)*t*t+p*t
-
-def Scaling_P3(p,t) :
- q = 1./p
- a = p+q-2
- b = 3-2*p-q
- return a*t**3 + b*t*t +p*t
-
-
-# -------------------- MESHING OPERATIONS ---------------------
-
-class Mesh_Descriptor :
-
- """
- Attributes
- relation Another mesh descriptor provides the mesh parameters
- parameters dictionnary of the mesh parameters
- size Point size
- transfinite Transfinite mesh (0 or 1)
- number Number of elements along a line (transfinite)
- progression Progression of element size (transfinite)
- recombine Recombine mesh or not
-
- Specific access :
- md.parameter_name = xxx -> the relation is destroyed (set to None)
- xxx = md.parameter_name -> if there is a relation, the effective
- parameter is looked for recursively
-
- Deep copying : a relation is set to the model instead of a true copy
- """
-
- List_Attr = ['size','transfinite','number','progression','recombine']
-
-
- def __init__(self) :
-
- self.relation = None
- self.parameters = {
- 'size' : 1. , # Point size
- 'transfinite': 0 , # Transfinite mesh (0 or 1)
- 'recombine' : 0 # Recombine mesh or not
- }
-
-
- def __setattr__(self, attr, value) :
-
- if attr in Mesh_Descriptor.List_Attr :
- self.relation = None
- self.parameters[attr] = value
-
- if attr == 'number' :
- self.transfinite = 1
-
- else :
- self.__dict__[attr] = value
-
-
- def __getattr__(self,attr) :
-
- if self.relation :
- return (getattr(self.relation,attr))
- else :
- if attr in self.parameters.keys() :
- return self.parameters[attr]
- else :
- raise AttributeError
-
-
- def __deepcopy__(self,visit) :
-
- md = copy.copy(self)
- md.parameters = copy.copy(self.parameters)
- md.relation = self
- return md
-
-
-
-class Mesh :
-
- """
-
- """
-
-
- def __init__(self, algo = 2, order = 1, gmsh='gmsh') :
-
- self.num_ph = 0
- self.num = 0
- self.order = order
- self.command = ['Mesh.Algorithm = ' + repr(algo) + ' ;']
- self.command += ['Mesh.ElementOrder = ' + repr(order) + ' ;']
- self.physicals = {}
- self.gmsh = gmsh
-
-
- def Physical(self, name, *l_lobj) :
-
- # Checking the name
- if type(name) <> type(' ') :
- raise 'First argument should be the name of the physical'
- if name in self.physicals.keys() :
- raise 'Physical '+name+' already exists'
-
- # treating the case of list of lists parameters
- l_obj = []
- for l in l_lobj :
- if type(l) == type([]) :
- l_obj = l_obj + l
- else :
- l_obj.append(l)
-
- # Checking all objects are geometric
- for obj in l_obj :
- if not Is_Geometric(obj) :
- raise "Non geometrical object : " + repr(obj) + " Physical = " + name
-
- cl = l_obj[0].Base_class()
- # Checking all objects are of the same dimension
- # ref_dim = l_obj[0]
- # for obj in l_obj[1:] :
- # if not ref_dim.Is_same_dimension(obj) :
- # raise "All objects are not of the same dimension : " + repr(obj)
-
- # Creation of the objects if necessary
- for obj in l_obj :
- obj.Gmsh(self)
-
- # Creation of the physical
- self.num_ph= self.num_ph + 1
- ch = name + '=' + `self.num_ph` + ';'
- self.command.append(ch)
- ch = 'Physical ' + cl + '(' + name + ') = {'
- for obj in l_obj :
- ch = ch + `obj.num` + ','
- ch = ch[:-1] + '};'
- self.command.append(ch)
-
- # Name of the physical
- name_gmsh = 'GM'+`self.num_ph`
- self.physicals[name] = name_gmsh
-
-
- def Save(self, file = 'fort.geo') :
-
- if os.path.isfile(file) :
- os.remove(file)
-
- f = open(file,'w')
- f.write(string.joinfields(self.command,'\n'))
- f.close()
-
-
- def View(self) :
-
- self.Save('fort.geo')
-# os.system('gmsh fort.geo')
- os.system(self.gmsh + ' fort.geo')
- os.remove('fort.geo')
-
-
- def Create(self, file = 'fort.19') :
-
- self.Save()
-# os.system('gmsh -3 fort.geo')
- os.system(self.gmsh + ' -3 fort.geo')
- os.rename('fort.msh',file)
-
-
- def Name(self, MA, CREA_GROUP_NO) :
-
- l_gma = []
- l_mcf = []
- for gma in self.physicals.keys() :
- l_gma.append(self.physicals[gma])
- l_mcf.append(_F(GROUP_MA = self.physicals[gma],NOM=gma))
-
- DEFI_GROUP(reuse = MA,
- MAILLAGE = MA,
- CREA_GROUP_MA = tuple(l_mcf),
- )
-
- SMESH_02 = CREA_MAILLAGE(
- MAILLAGE = MA,
- DETR_GROUP_MA = _F(GROUP_MA = tuple(l_gma)),
- )
-
- DETRUIRE(CONCEPT = _F(NOM = MA), INFO=1)
-
- if CREA_GROUP_NO == 'OUI' :
- DEFI_GROUP(reuse = SMESH_02,
- MAILLAGE = SMESH_02,
- CREA_GROUP_NO = _F(TOUT_GROUP_MA = 'OUI'),
- )
-
- else :
-# Traitement des GROUP_NO qui sont des points
- info_gno = SMESH_02.LIST_GROUP_NO()
- l_gno = []
- for gno in info_gno :
- if gno[1] == 1 : l_gno.append(gno[0])
-
- l_gma = []
- for gma in self.physicals.keys() :
- nom_gmsh = self.physicals[gma]
- if nom_gmsh in l_gno :
- l_gma.append(gma)
-
- if l_gma :
- DEFI_GROUP(reuse = SMESH_02,
- MAILLAGE = SMESH_02,
- CREA_GROUP_NO = _F(GROUP_MA = tuple(l_gma)),
- )
-
- return SMESH_02
-
-
-
- def LIRE_GMSH(self,
- UNITE_GMSH = 19,
- UNITE_MAILLAGE = 20,
- MODI_QUAD = 'NON',
- CREA_GROUP_NO = 'OUI'
- ) :
-
- """
- Lecture du maillage (format Aster) a partir de sa definition
- (format sup_gmsh)
- UNITE_GMSH = Numero d'unite logique pour le fichier msh
- UNITE_MAILLAGE = Numero d'unite logique pour le fichier mail
- MODI_QUAD = 'OUI' si line->quad, 'NON' sinon
- CREA_GROUP_NO = 'OUI' si on cree les group_no, 'NON' sinon
- """
-
- nom_gmsh = 'fort.' + repr(UNITE_GMSH)
- self.Create(nom_gmsh)
-
- PRE_GMSH(UNITE_GMSH=UNITE_GMSH, UNITE_MAILLAGE=UNITE_MAILLAGE)
-
- SMESH_00 = LIRE_MAILLAGE(UNITE = UNITE_MAILLAGE)
- DEFI_FICHIER(ACTION='LIBERER',UNITE = UNITE_GMSH)
- DEFI_FICHIER(ACTION='LIBERER',UNITE = UNITE_MAILLAGE)
-
- if MODI_QUAD == 'OUI' and self.order == 2 :
- raise 'The finite elements are already of second order'
-
- if MODI_QUAD == 'OUI' and self.order <> 2 :
- SMESH_01 = CREA_MAILLAGE(
- MAILLAGE = SMESH_00,
- LINE_QUAD = _F(TOUT = 'OUI')
- )
- DETRUIRE(CONCEPT=_F(NOM=SMESH_00), INFO=1)
- SMESH_00 = SMESH_01
-
- SMESH_00 = self.Name(SMESH_00,CREA_GROUP_NO)
-
- return SMESH_00
+++ /dev/null
-#@ MODIF t_fonction Utilitai DATE 31/05/2005 AUTEUR DURAND C.DURAND
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2005 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-from Numeric import *
-import copy
-import types
-
-def interp(typ_i,val,x1,x2,y1,y2) :
- if typ_i==['LIN','LIN']: return y1+(y2-y1)*(val-x1)/(x2-x1)
- if typ_i==['LIN','LOG']: return exp(log(y1)+(val-x1)*(log(y2)-log(y1))/(x2-x1))
- if typ_i==['LOG','LOG']: return exp(log(y1)+(log(val)-log(x1))*(log(y2)-log(y1))/(log(x2)-log(x1)))
- if typ_i==['LOG','LIN']: return y1+(log(val)-log(x1))*(y2-y1)/(log(x2)-log(x1))
- if typ_i[0]=='NON' :
- if val==x1 : return y1
- elif val==x2 : return y2
- else : raise StandardError, 'fonction : interpolation NON'
-def is_ordo(liste) :
- listb=dict([(i,0) for i in liste]).keys()
- listb.sort()
- return liste==listb
-
-class t_fonction :
- ### Classe pour fonctions réelles, équivalent au type aster = fonction_sdaster
- def __init__(self,vale_x,vale_y,para) :
- # création d'un objet fonction
- # vale_x et vale_y sont des listes de réels de meme longueur
- # para est un dictionnaire contenant les entrées PROL_DROITE, PROL_GAUCHE et INTERPOL (cf sd ASTER)
- pk=para.keys()
- pk.sort()
- if pk!=['INTERPOL','NOM_PARA','NOM_RESU','PROL_DROITE','PROL_GAUCHE'] :
- raise StandardError, 'fonction : parametres incorrects'
- if para['INTERPOL'] not in [['NON','NON'],['LIN','LIN'],['LIN','LOG'],['LOG','LOG'],['LOG','LIN'],] :
- raise StandardError, 'fonction : parametre INTERPOL incorrect'
- if para['PROL_DROITE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
- raise StandardError, 'fonction : parametre PROL_DROITE incorrect'
- if para['PROL_GAUCHE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
- raise StandardError, 'fonction : parametre PROL_GAUCHE incorrect'
- self.vale_x = array(vale_x)
- self.vale_y = array(vale_y)
- self.para = para
- if len(self.vale_x)!=len(self.vale_y) :
- raise StandardError, 'fonction : longueur abscisse <> longueur ordonnées'
- if not is_ordo(self.vale_x) :
- raise StandardError, 'fonction : abscisses non strictement croissantes'
-
- def __add__(self,other) :
- # addition avec une autre fonction ou un nombre, par surcharge de l'opérateur +
- if isinstance(other,t_fonction):
- para=copy.copy(self.para)
- vale_x,para['PROL_GAUCHE'],para['PROL_DROITE']=self.homo_support(other)
- fff=self.evalfonc(vale_x)
- ggg=other.evalfonc(vale_x)
- if isinstance(self,t_fonction_c): return t_fonction_c(vale_x,fff.vale_y+ggg.vale_y,para)
- else : return t_fonction(vale_x,fff.vale_y+ggg.vale_y,para)
- elif type(other) in [types.FloatType,types.IntType,types.ComplexType] :
- if isinstance(self,t_fonction_c): return t_fonction_c(self.vale_x,self.vale_y+other,self.para)
- else : return t_fonction(self.vale_x,self.vale_y+other,self.para)
- else: raise StandardError, 'fonctions : erreur de type dans __add__'
-
- def __mul__(self,other) :
- # multiplication avec une autre fonction ou un nombre, par surcharge de l'opérateur *
- if isinstance(other,t_fonction):
- para=copy.copy(self.para)
- vale_x,para['PROL_GAUCHE'],para['PROL_DROITE']=self.homo_support(other)
- fff=self.evalfonc(vale_x)
- ggg=other.evalfonc(vale_x)
- if isinstance(self,t_fonction_c): return t_fonction_c(vale_x,fff.vale_y*ggg.vale_y,para)
- else : return t_fonction(vale_x,fff.vale_y*ggg.vale_y,para)
- elif type(other) in [types.FloatType,types.IntType] :
- return t_fonction(self.vale_x,self.vale_y*other,self.para)
- elif type(other) ==types.ComplexType :
- return t_fonction_c(self.vale_x,self.vale_y*other,self.para)
- else: raise StandardError, 'fonctions : erreur de type dans __mul__'
-
- def __repr__(self) :
- # affichage de la fonction en double colonne
- texte=[]
- for i in range(len(self.vale_x)) :
- texte.append('%f %f' % (self.vale_x[i],self.vale_y[i]))
- return '\n'.join(texte)
-
- def __getitem__(self,other) :
- # composition de deux fonction F[G]=FoG=F(G(x))
- para=copy.copy(self.para)
- if other.para['NOM_RESU']!=self.para['NOM_PARA'] :
- raise StandardError,'''composition de fonctions : NOM_RESU1 et NOM_PARA2 incohérents '''
- para['NOM_PARA']==other.para['NOM_PARA']
- return t_fonction(other.vale_x,map(self,other.vale_y),para)
-
- def __call__(self,val,tol=1.e-6):
- # méthode pour évaluer f(x)
- # tolérance, par défaut 1.e-6 en relatif sur la longueur de l'intervalle
- # adjacent, pour capter les erreurs d'arrondi en cas de prolongement exclu
- i=searchsorted(self.vale_x,val)
- n=len(self.vale_x)
- if i==0 :
- if self.para['PROL_GAUCHE']=='EXCLU' :
- eps_g=(val-self.vale_x[0] )/(self.vale_x[1] -self.vale_x[0])
- if abs(eps_g)<=tol : return self.vale_y[0]
- else : raise StandardError, 'fonction évaluée hors du domaine de définition'
- else :
- if self.para['PROL_GAUCHE']=='CONSTANT' : return self.vale_y[0]
- if self.para['PROL_GAUCHE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val,self.vale_x[0],
- self.vale_x[1],
- self.vale_y[0],
- self.vale_y[1])
- elif i==n :
- if self.para['PROL_DROITE']=='EXCLU' :
- eps_d=(val-self.vale_x[-1])/(self.vale_x[-1]-self.vale_x[-2])
- if abs(eps_d)<=tol : return self.vale_y[-1]
- else : raise StandardError, 'fonction évaluée hors du domaine de définition'
- else :
- if self.para['PROL_DROITE']=='CONSTANT' : return self.vale_y[-1]
- if self.para['PROL_DROITE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val,self.vale_x[-1],
- self.vale_x[-2],
- self.vale_y[-1],
- self.vale_y[-2])
- else :
- return interp(self.para['INTERPOL'],val,self.vale_x[i-1],
- self.vale_x[i],
- self.vale_y[i-1],
- self.vale_y[i])
-
- def homo_support(self,other) :
- # renvoie le support d'abscisses homogénéisé entre self et other
- # i.e. si prolongement exclu, on retient plus grand min ou plus petit max, selon
- # si prolongement autorisé, on conserve les abscisses d'une fonction, extrapolantes
- # sur l'autre.
- # Pour les points intermédiaires : union et tri des valeurs des vale_x réunis.
- if other.vale_x[0]>self.vale_x[0]:
- if other.para['PROL_GAUCHE']!='EXCLU' : f_g=self
- else : f_g=other
- else :
- if self.para['PROL_GAUCHE'] !='EXCLU' : f_g=other
- else : f_g=self
- val_min =f_g.vale_x[0]
- prol_gauche=f_g.para['PROL_GAUCHE']
- if self.vale_x[0]>other.vale_x[0]:
- if other.para['PROL_DROITE']!='EXCLU' : f_d=self
- else : f_d=other
- else :
- if self.para['PROL_DROITE'] !='EXCLU' : f_d=other
- else : f_d=self
- val_max =f_d.vale_x[-1]
- prol_droite=f_d.para['PROL_DROITE']
- vale_x=concatenate((self.vale_x,other.vale_x))
- vale_x=clip(vale_x,val_min,val_max)
- vale_x=sort(dict([(i,0) for i in vale_x]).keys())
- return vale_x,prol_gauche,prol_droite
-
- def cut(self,rinf,rsup,prec,crit='RELATIF') :
- # renvoie la fonction self dont on a 'coupé' les extrémités en x=rinf et x=rsup
- # pour la recherche de rinf et rsup dans la liste d'abscisses :
- # prec=precision crit='absolu' ou 'relatif'
- para=copy.copy(self.para)
- para['PROL_GAUCHE']='EXCLU'
- para['PROL_DROITE']='EXCLU'
- if crit=='ABSOLU' : rinf_tab=greater(abs(self.vale_x-rinf),prec)
- elif crit=='RELATIF': rinf_tab=greater(abs(self.vale_x-rinf),prec*rinf)
- else : raise StandardError, 'fonction : cut : critère absolu ou relatif'
- if crit=='ABSOLU' : rsup_tab=greater(abs(self.vale_x-rsup),prec)
- elif crit=='RELATIF': rsup_tab=greater(abs(self.vale_x-rsup),prec*rsup)
- else : raise StandardError, 'fonction : cut : critère absolu ou relatif'
- if alltrue(rinf_tab) : i=searchsorted(self.vale_x,rinf)
- else : i=rinf_tab.tolist().index(0)+1
- if alltrue(rsup_tab) : j=searchsorted(self.vale_x,rsup)
- else : j=rsup_tab.tolist().index(0)
- vale_x=array([rinf,]+self.vale_x.tolist()[i:j]+[rsup,])
- vale_y=array([self(rinf),]+self.vale_y.tolist()[i:j]+[self(rsup),])
- return t_fonction(vale_x,vale_y,para)
-
- def cat(self,other,surcharge) :
- # renvoie une fonction concaténée avec une autre, avec règles de surcharge
- para=copy.copy(self.para)
- if self.para['INTERPOL']!=other.para['INTERPOL'] : raise StandardError, 'concaténation de fonctions à interpolations différentes'
- if min(self.vale_x)<min(other.vale_x) :
- f1=self
- f2=other
- else :
- f1=other
- f2=self
- para['PROL_GAUCHE']=f1.para['PROL_GAUCHE']
- if surcharge=='GAUCHE' :
- i=searchsorted(f2.vale_x,f1.vale_x[-1])
- if i!=len(f2.vale_x) : para['PROL_DROITE']=f2.para['PROL_DROITE']
- else : para['PROL_DROITE']=f1.para['PROL_DROITE']
- vale_x=array(f1.vale_x.tolist()+f2.vale_x[i:].tolist())
- vale_y=array(f1.vale_y.tolist()+f2.vale_y[i:].tolist())
- elif surcharge=='DROITE' :
- i=searchsorted(f1.vale_x,f2.vale_x[0])
- if i!=len(f1.vale_x) : para['PROL_DROITE']=f2.para['PROL_DROITE']
- else : para['PROL_DROITE']=f1.para['PROL_DROITE']
- vale_x=array(f1.vale_x[:i].tolist()+f2.vale_x.tolist())
- vale_y=array(f1.vale_y[:i].tolist()+f2.vale_y.tolist())
- return t_fonction(vale_x,vale_y,para)
-
- def tabul(self) :
- # mise en forme de la fonction selon un vecteur unique (x1,y1,x2,y2,...)
- __tab=array([self.vale_x,self.vale_y])
- return ravel(transpose(__tab)).tolist()
-
- def extreme(self) :
- # renvoie un dictionnaire des valeurs d'ordonnées min et max
- val_min=min(self.vale_y)
- val_max=max(self.vale_y)
- vm={}
- vm['min']=[[self.vale_y[i],self.vale_x[i]] for i in range(len(self.vale_x))\
- if self.vale_y[i]==val_min]
- vm['max']=[[self.vale_y[i],self.vale_x[i]] for i in range(len(self.vale_x))\
- if self.vale_y[i]==val_max]
- vm['min'].sort()
- vm['max'].sort()
- for item in vm['min'] : item.reverse()
- for item in vm['max'] : item.reverse()
- return vm
-
- def trapeze(self,coef) :
- # renvoie la primitive de la fonction, calculée avec la constante d'intégration 'coef'
- trapz = zeros(len(self.vale_y),Float)
- trapz[0] = coef
- trapz[1:] = (self.vale_y[1:]+self.vale_y[:-1])/2*(self.vale_x[1:]-self.vale_x[:-1])
- prim_y=cumsum(trapz)
- para=copy.copy(self.para)
- para['PROL_GAUCHE']='EXCLU'
- para['PROL_DROITE']='EXCLU'
- if para['NOM_RESU'][:4]=='VITE' : para['NOM_RESU']='DEPL'
- elif para['NOM_RESU'][:4]=='ACCE' : para['NOM_RESU']='VITE'
- return t_fonction(self.vale_x,prim_y,para)
-
- def simpson(self,coef) :
- # renvoie la primitive de la fonction, calculée avec la constante d'intégration 'coef'
- para=copy.copy(self.para)
- para['PROL_GAUCHE']='EXCLU'
- para['PROL_DROITE']='EXCLU'
- if para['NOM_RESU'][:4]=='VITE' : para['NOM_RESU']='DEPL'
- elif para['NOM_RESU'][:4]=='ACCE' : para['NOM_RESU']='VITE'
- fm = self.vale_y[0]
- fb = self.vale_y[1]
- h2 = self.vale_x[1] - self.vale_x[0]
- tabl = [coef,coef +(fb+fm)*h2/2.]
- prim_y = copy.copy(tabl)
- iperm = 0
- ip = (1,0)
- eps = 1.E-4
- for i in range(2,len(self.vale_x)) :
- h1 = h2
- h2 = self.vale_x[i] - self.vale_x[i-1]
- bma = h1 + h2
- fa = fm
- fm = fb
- fb = self.vale_y[i]
- deltah = h2 - h1
- if h1==0. or h2==0. or abs( deltah / h1 ) <= eps :
- ct = (1.,4.,1.)
- else :
- epsi = deltah / (h1*h2)
- ct = (1.-epsi*h2,2.+epsi*deltah,1.+epsi*h1)
- tabl[iperm] = tabl[iperm] + (bma/6.)*(ct[0]*fa+ct[1]*fm+ct[2]*fb)
- prim_y.append(tabl[iperm])
- iperm = ip[iperm]
- return t_fonction(self.vale_x,prim_y,para)
-
- def derive(self) :
- # renvoie la dérivée de la fonction
- pas=self.vale_x[1:]-self.vale_x[:-1]
- pentes=(self.vale_y[1:]-self.vale_y[:-1])/(self.vale_x[1:]-self.vale_x[:-1])
- derive=(pentes[1:]*pas[1:]+pentes[:-1]*pas[:-1])/(pas[1:]+pas[:-1])
- derv_y=[pentes[0]]+derive.tolist()+[pentes[-1]]
- para=copy.copy(self.para)
- para['PROL_GAUCHE']='EXCLU'
- para['PROL_DROITE']='EXCLU'
- if para['NOM_RESU'][:4]=='DEPL' : para['NOM_RESU']='VITE'
- elif para['NOM_RESU'][:4]=='VITE' : para['NOM_RESU']='ACCE'
- return t_fonction(self.vale_x,derv_y,para)
-
- def inverse(self) :
- # renvoie l'inverse de la fonction
- # on intervertit vale_x et vale_y, on swape interpolation
- para=copy.copy(self.para)
- para['NOM_RESU']='TOUTRESU'
- para['NOM_PARA']=self.para['NOM_PARA']
- para['INTERPOL'].reverse()
- if para['PROL_GAUCHE']=='CONSTANT' : para['PROL_GAUCHE']='EXCLU'
- if para['PROL_DROITE']=='CONSTANT' : para['PROL_DROITE']='EXCLU'
- vale_x=self.vale_y
- vale_y=self.vale_x
- if not is_ordo(vale_x) :
- vale_x=vale_x[::-1]
- vale_y=vale_y[::-1]
- return t_fonction(vale_x,vale_y,para)
-
- def abs(self) :
- # renvoie la mm fonction avec valeur absolue des ordonnées
- para=copy.copy(self.para)
- if para['PROL_GAUCHE']=='LINEAIRE' : para['PROL_GAUCHE']='EXCLU'
- if para['PROL_DROITE']=='LINEAIRE' : para['PROL_DROITE']='EXCLU'
- return t_fonction(self.vale_x,absolute(self.vale_y),para)
-
- def evalfonc(self,liste_val) :
- # renvoie la mm fonction interpolée aux points définis par la liste 'liste_val'
- return t_fonction(liste_val,map(self,liste_val),self.para)
-
- def sup(self,other) :
- # renvoie l'enveloppe supérieure de self et other
- para=copy.copy(self.para)
-# commentaire : pour les prolongements et l'interpolation, c'est self
-# qui prime sur other
- vale_x=self.vale_x.tolist()+other.vale_x.tolist()
-# on ote les abscisses doublons
- vale_x=dict([(i,0) for i in vale_x]).keys()
- vale_x.sort()
- vale_x=array(vale_x)
- vale_y1=map(self, vale_x)
- vale_y2=map(other,vale_x)
- vale_y=map(max,vale_y1,vale_y2)
- return t_fonction(vale_x,vale_y,para)
-
- def inf(self,other) :
- # renvoie l'enveloppe inférieure de self et other
- para=copy.copy(self.para)
-# commentaire : pour les prolongements et l'interpolation, c'est self
-# qui prime sur other
- vale_x=self.vale_x.tolist()+other.vale_x.tolist()
-# on ote les abscisses doublons
- vale_x=dict([(i,0) for i in vale_x]).keys()
- vale_x.sort()
- vale_x=array(vale_x)
- vale_y1=map(self, vale_x)
- vale_y2=map(other,vale_x)
- vale_y=map(min,vale_y1,vale_y2)
- return t_fonction(vale_x,vale_y,para)
-
- def suppr_tend(self) :
- # pour les corrections d'accélérogrammes
- # suppression de la tendance moyenne d'une fonction
- para=copy.copy(self.para)
- xy=sum(self.vale_x*self.vale_y)
- x0=sum(self.vale_x)
- y0=sum(self.vale_y)
- xx=sum(self.vale_x*self.vale_x)
- n=len(self.vale_x)
- a1 = ( n*xy - x0*y0) / (n*xx - x0*x0)
- a0 = (xx*x0 - x0*xy) / (n*xx - x0*x0)
- return t_fonction(self.vale_x,self.vale_y-a1*self.vale_x-a0,self.para)
-
- def normel2(self) :
- # norme de la fonction
- __ex=self*self
- __ex=__ex.trapeze(0.)
- return sqrt(__ex.vale_y[-1])
-
- def fft(self,methode) :
- # renvoie la transformée de Fourier rapide FFT
- import FFT
- para=copy.copy(self.para)
- para['NOM_PARA']='FREQ'
- if self.para['NOM_PARA']!='INST' :
- raise StandardError, 'fonction réelle : FFT : NOM_PARA=INST pour une transformée directe'
- pas = self.vale_x[1]-self.vale_x[0]
- for i in range(1,len(self.vale_x)) :
- ecart = abs(((self.vale_x[i]-self.vale_x[i-1])-pas)/pas)
- if ecart>1.e-2 :
- raise StandardError, 'fonction réelle : FFT : la fonction doit etre à pas constant'
- n=int(log(len(self.vale_x))/log(2))
- if methode=='TRONCATURE' :
- vale_y=self.vale_y[:2**n]
- elif methode=='PROL_ZERO' :
- vale_y=self.vale_y.tolist()
- vale_y=vale_y+[0.]*(2**(n+1)-len(self.vale_x))
- vale_y=array(vale_y)
- vect=FFT.fft(vale_y)
- pasfreq =1./(pas*(len(vect)-1))
- vale_x =[pasfreq*i for i in range(len(vect))]
- vale_y =vect*pas
- return t_fonction_c(vale_x,vale_y,para)
-
-class t_fonction_c(t_fonction) :
- ### Classe pour fonctions complexes, équivalent au type aster = fonction_c
- def tabul(self) :
- # mise en forme de la fonction selon un vecteur unique (x1,yr1,yi1,x2,yr2,yr2,...)
- __tab=array([self.vale_x,self.vale_y.real,self.vale_y.imag])
- return ravel(transpose(__tab)).tolist()
-
- def __repr__(self) :
- # affichage de la fonction en double colonne
- texte=[]
- for i in range(len(self.vale_x)) :
- texte.append('%f %f + %f .j' % (self.vale_x[i],self.vale_y[i].real,self.vale_y[i].imag))
- return '\n'.join(texte)
-
- def fft(self,methode,syme) :
- # renvoie la transformée de Fourier rapide FFT (sens inverse)
- import FFT
- para=copy.copy(self.para)
- para['NOM_PARA']='INST'
- if self.para['NOM_PARA']!='FREQ' :
- raise StandardError, 'fonction complexe : FFT : NOM_PARA=FREQ pour une transformée directe'
- pas = self.vale_x[1]-self.vale_x[0]
- for i in range(1,len(self.vale_x)) :
- ecart = abs(((self.vale_x[i]-self.vale_x[i-1])-pas)/pas)
- if ecart>1.e-3 :
- raise StandardError, 'fonction complexe : FFT : la fonction doit etre à pas constant'
- n=int(log(len(self.vale_x))/log(2))
- if syme=='OUI' and len(self.vale_x)==2**n :
- vale_fonc=self.vale_y
- elif syme=='NON' and len(self.vale_x)==2**n :
- vale_fonc=self.vale_y.tolist()
- vale_fon1=vale_fonc[:]
- vale_fon1.reverse()
- vale_fonc=vale_fonc+vale_fon1
- vale_fonc=array(vale_fonc)
- elif syme=='NON' and len(self.vale_x)!=2**n and methode=='PROL_ZERO' :
- vale_fonc=self.vale_y.tolist()+[complex(0.)]*(2**(n+1)-len(self.vale_x))
- vale_fon1=vale_fonc[:]
- vale_fon1.reverse()
- vale_fonc=vale_fonc+vale_fon1
- vale_fonc=array(vale_fonc)
- elif syme=='NON' and len(self.vale_x)!=2**n and methode=='TRONCATURE' :
- vale_fonc=self.vale_y[:2**n]
- vale_fonc=vale_fonc.tolist()
- vale_fon1=vale_fonc[:]
- vale_fon1.reverse()
- vale_fonc=vale_fonc+vale_fon1
- vale_fonc=array(vale_fonc)
- if syme=='OUI' and len(self.vale_x)!=2**n :
- raise StandardError, 'fonction complexe : FFT : syme=OUI et nombre de points<>2**n'
- part1=vale_fonc[:len(vale_fonc)/2+1]
- part2=vale_fonc[1:len(vale_fonc)/2]
- part2=conjugate(part2)
- part2=part2.tolist()
- part2.reverse()
- vale_fonc=array(part1.tolist()+part2)
- vect=FFT.inverse_fft(vale_fonc)
- vect=vect.real
- pasfreq =1./(pas*(len(vect)-1))
- vale_x =[pasfreq*i for i in range(len(vect))]
- pas2 =(1./self.vale_x[-1])*((len(self.vale_x))/float(len(vect)))
- vale_y =vect/pas2
- return t_fonction(vale_x,vale_y,para)
-
-
-class t_nappe :
- ### Classe pour nappes, équivalent au type aster = nappe_sdaster
- def __init__(self,vale_para,l_fonc,para) :
- # création d'un objet nappe
- # vale_para est la liste de valeur des parametres (mot clé PARA dans DEFI_NAPPE)
- # para est un dictionnaire contenant les entrées PROL_DROITE, PROL_GAUCHE et INTERPOL (cf sd ASTER)
- # l_fonc est la liste des fonctions, de cardinal égal à celui de vale_para
- pk=para.keys()
- pk.sort()
- if pk!=['INTERPOL','NOM_PARA','NOM_PARA_FONC','NOM_RESU','PROL_DROITE','PROL_GAUCHE'] :
- raise StandardError, 'nappe : parametres incorrects'
- if para['INTERPOL'] not in [['NON','NON'],['LIN','LIN'],
- ['LIN','LOG'],['LOG','LOG'],['LOG','LIN'],] :
- raise StandardError, 'nappe : parametre INTERPOL incorrect'
- if para['PROL_DROITE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
- raise StandardError, 'nappe : parametre PROL_DROITE incorrect'
- if para['PROL_GAUCHE'] not in ['EXCLU','CONSTANT','LINEAIRE'] :
- raise StandardError, 'nappe : parametre PROL_GAUCHE incorrect'
- self.vale_para = array(vale_para)
- if type(l_fonc) not in (types.ListType,types.TupleType) :
- raise StandardError, 'nappe : la liste de fonctions fournie n est pas une liste'
- if len(l_fonc)!=len(vale_para) :
- raise StandardError, 'nappe : nombre de fonctions différent du nombre de valeurs du paramètre'
- for f in l_fonc :
- if not isinstance(f,t_fonction) and not isinstance(f,t_fonction_c) :
- raise StandardError, 'nappe : les fonctions fournies ne sont pas du bon type'
- self.l_fonc = l_fonc
- self.para = para
-
- def __call__(self,val1,val2):
- # méthode pour évaluer nappe(val1,val2)
- i=searchsorted(self.vale_para,val1)
- n=len(self.vale_para)
- if i==0 :
- if val1==self.vale_para[0] : return self.l_fonc[0](val2)
- if val1 <self.vale_para[0] :
- if self.para['PROL_GAUCHE']=='EXCLU' : raise StandardError, 'nappe évaluée hors du domaine de définition'
- if self.para['PROL_GAUCHE']=='CONSTANT' : return self.l_fonc[0](val2)
- if self.para['PROL_GAUCHE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val1,
- self.vale_para[0],
- self.vale_para[1],
- self.l_fonc[0](val2),
- self.l_fonc[1](val2))
- elif i==n :
- if val1==self.vale_para[-1] : return self.l_fonc[-1](val2)
- if val1 >self.vale_para[-1] :
- if self.para['PROL_DROITE']=='EXCLU' : raise StandardError, 'nappe évaluée hors du domaine de définition'
- if self.para['PROL_DROITE']=='CONSTANT' : return self.l_fonc[-1](val2)
- if self.para['PROL_DROITE']=='LINEAIRE' : return interp(self.para['INTERPOL'],val1,
- self.vale_para[-1],
- self.vale_para[-2],
- self.l_fonc[-1](val2),
- self.l_fonc[-2](val2))
- else :
- return interp(self.para['INTERPOL'],val1,self.vale_para[i-1],
- self.vale_para[i],
- self.l_fonc[i-1](val2),
- self.l_fonc[i](val2))
-
- def __add__(self,other) :
- # addition avec une autre nappe ou un nombre, par surcharge de l'opérateur +
- l_fonc=[]
- if isinstance(other,t_nappe):
- if self.vale_para!=other.vale_para : raise StandardError, 'nappes à valeurs de paramètres différentes'
- for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]+other.l_fonc[i])
- elif type(other) in [types.FloatType,types.IntType] :
- for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]+other)
- else: raise StandardError, 't_nappe : erreur de type dans __add__'
- return t_nappe(self.vale_para,l_fonc,self.para)
-
- def __mul__(self,other) :
- # multiplication avec une autre fonction ou un nombre, par surcharge de l'opérateur *
- l_fonc=[]
- if isinstance(other,t_nappe):
- if self.vale_para!=other.vale_para : raise StandardError, 'nappes à valeurs de paramètres différentes'
- for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]*other.l_fonc[i])
- elif type(other) in [types.FloatType,types.IntType] :
- for i in range(len(self.l_fonc)) : l_fonc.append(self.l_fonc[i]*other)
- else: raise StandardError, 't_nappe : erreur de type dans __mul__'
- return t_nappe(self.vale_para,l_fonc,self.para)
-
- def __repr__(self) :
- # affichage de la nappe en double colonne
- texte=[]
- for i in range(len(self.vale_para)) :
- texte.append('paramètre : %f' % self.vale_para[i])
- texte.append(repr(self.l_fonc[i]))
- return '\n'.join(texte)
-
- def homo_support(self,other) :
- # renvoie la nappe self avec un support union de celui de self et de other
- # le support est la discrétisation vale_para et les discrétisations des fonctions
- if self==other : return self
- if self.para!=other.para : raise StandardError, 'combinaison de nappes à caractéristiques interpolation et prolongement différentes'
- vale_para=self.vale_para.tolist()+other.vale_para.tolist()
- vale_para=dict([(i,0) for i in vale_para]).keys()
- vale_para.sort()
- vale_para=array(vale_para)
- l_fonc=[]
- for val in vale_para :
- if val in self.vale_para : l_fonc.append(self.l_fonc[searchsorted(self.vale_para,val)])
- elif val in other.vale_para :
- other_fonc=other.l_fonc[searchsorted(other.vale_para,val)]
- new_vale_x=other_fonc.vale_x
- new_para =other_fonc.para
- new_vale_y=[self(x) for x in new_vale_x]
- if isinstance(other_fonc,t_fonction) :
- l_fonc.append(t_fonction(new_vale_x,new_vale_y,new_para))
- if isinstance(other_fonc,t_fonction_c) :
- l_fonc.append(t_fonction_c(new_vale_x,new_vale_y,new_para))
- else : raise StandardError, 'combinaison de nappes : incohérence'
- return t_nappe(vale_para,l_fonc,self.para)
-
- def extreme(self) :
- # renvoie un dictionnaire des valeurs d'ordonnées min et max
- val_min=min([min(fonc.vale_y) for fonc in self.l_fonc])
- val_max=max([max(fonc.vale_y) for fonc in self.l_fonc])
- vm={'min':[],'max':[]}
- for i in range(len(self.vale_para)) :
- for j in range(len(self.l_fonc[i].vale_y)) :
- y = self.l_fonc[i].vale_y[j]
- if y==val_min : vm['min'].append([y,self.l_fonc[i].vale_x[j],self.vale_para[i]])
- if y==val_max : vm['max'].append([y,self.l_fonc[i].vale_x[j],self.vale_para[i]])
- vm['min'].sort()
- vm['max'].sort()
- for item in vm['min'] : item.reverse()
- for item in vm['max'] : item.reverse()
- return vm
+++ /dev/null
-#@ MODIF transpose Utilitai DATE 14/09/2004 AUTEUR MCOURTOI M.COURTOIS
-# -*- coding: iso-8859-1 -*-
-# CONFIGURATION MANAGEMENT OF EDF VERSION
-# ======================================================================
-# COPYRIGHT (C) 1991 - 2004 EDF R&D WWW.CODE-ASTER.ORG
-# THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
-# IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
-# THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
-# (AT YOUR OPTION) ANY LATER VERSION.
-#
-# THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
-# WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
-# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
-# GENERAL PUBLIC LICENSE FOR MORE DETAILS.
-#
-# YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
-# ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
-# 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
-# ======================================================================
-
-
-
-######################################################################
-#### méthode de transposition de double liste
-#### à résorber quand on aura numarray et possibilité
-#### d opérations sur des arrays de strings
-######################################################################
-def transpose(liste):
- n=range(len(liste[0]))
- m=range(len(liste))
- liste_t=[[] for i in n]
- for i in n :
- for j in m :
- liste_t[i].append(liste[j][i])
- return liste_t
