Salome HOME
Modif V6_4_°
[tools/eficas.git] / Aster / Cata / cataSTA8 / Macro / post_k1_k2_k3_ops.py
1 #@ MODIF post_k1_k2_k3_ops Macro  DATE 06/11/2006   AUTEUR GALENNE E.GALENNE 
2 # -*- coding: iso-8859-1 -*-
3 #            CONFIGURATION MANAGEMENT OF EDF VERSION
4 # ======================================================================
5 # COPYRIGHT (C) 1991 - 2006  EDF R&D                  WWW.CODE-ASTER.ORG
6 # THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY  
7 # IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY  
8 # THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR     
9 # (AT YOUR OPTION) ANY LATER VERSION.                                                  
10 #                                                                       
11 # THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT   
12 # WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF            
13 # MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU      
14 # GENERAL PUBLIC LICENSE FOR MORE DETAILS.                              
15 #                                                                       
16 # YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE     
17 # ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,         
18 #    1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.        
19 # ======================================================================
20
21 def veri_tab(tab,nom,ndim) :
22    from Utilitai.Utmess     import UTMESS
23    macro = 'POST_K1_K2_K3'
24    for label in ('DX','DY','COOR_X','COOR_Y','ABSC_CURV') :
25        if label not in tab.para :
26           message='le label '+label+' doit etre present dans la table : '+nom
27           UTMESS('F', macro, message)
28    if ndim==3 :
29       if 'DZ'     not in tab.para :
30           message='le label DZ doit etre present dans la table : '+nom
31           UTMESS('F', macro, message)
32       if 'COOR_Z' not in tab.para :
33           message='le label COOR_Z doit etre present dans la table : '+nom
34           UTMESS('F', macro, message)
35
36 def cross_product(a,b):
37     cross = [0]*3
38     cross[0] = a[1]*b[2]-a[2]*b[1]
39     cross[1] = a[2]*b[0]-a[0]*b[2]
40     cross[2] = a[0]*b[1]-a[1]*b[0]
41     return cross
42  
43 def moy(t):
44     m = 0
45     for value in t :
46       m += value
47     return (m/len(t))
48      
49 def post_k1_k2_k3_ops(self,MODELISATION,FOND_FISS,MATER,RESULTAT,
50                    TABL_DEPL_SUP,TABL_DEPL_INF,ABSC_CURV_MAXI,PREC_VIS_A_VIS,
51                    TOUT_ORDRE,NUME_ORDRE,LIST_ORDRE,INST,LIST_INST,SYME_CHAR,
52                    INFO,VECT_K1,TITRE,**args):
53    """
54    Macro POST_K1_K2_K3
55    Calcul des facteurs d'intensité de contraintes en 2D et en 3D
56    par extrapolation des sauts de déplacements sur les lèvres de
57    la fissure. Produit une table.
58    """
59    import aster
60    import string
61    import copy
62    import math
63    import Numeric
64    from Numeric import array,asarray,Float,sqrt,matrixmultiply,transpose,sign,resize,dot,multiply
65    from math import pi
66    from types import ListType, TupleType
67    from Accas import _F
68    from Utilitai.Table      import Table, merge
69    EnumTypes = (ListType, TupleType)
70
71    macro = 'POST_K1_K2_K3'
72    from Accas               import _F
73    from Utilitai.Utmess     import UTMESS
74
75    ier = 0
76    # La macro compte pour 1 dans la numerotation des commandes
77    self.set_icmd(1)
78
79    # Le concept sortant (de type table_sdaster ou dérivé) est tab
80    self.DeclareOut('tabout', self.sd)
81    
82    # On importe les definitions des commandes a utiliser dans la macro
83    # Le nom de la variable doit etre obligatoirement le nom de la commande
84    CREA_TABLE    = self.get_cmd('CREA_TABLE')
85    CALC_TABLE    = self.get_cmd('CALC_TABLE')
86    POST_RELEVE_T    = self.get_cmd('POST_RELEVE_T')
87    DETRUIRE      = self.get_cmd('DETRUIRE')
88    DEFI_GROUP      = self.get_cmd('DEFI_GROUP')
89    MACR_LIGN_COUPE      = self.get_cmd('MACR_LIGN_COUPE')
90
91 #   ------------------------------------------------------------------
92 #                         CARACTERISTIQUES MATERIAUX
93 #   ------------------------------------------------------------------
94    matph = aster.getvectjev( string.ljust(MATER.nom,8)+'.MATERIAU.NOMRC         ')
95    phenom=None
96    for cmpt in matph :
97        if cmpt[:4]=='ELAS' :
98           phenom=cmpt
99           break
100    if phenom==None : UTMESS('F', macro, 'IL FAUT DEFINIR ELAS DANS DEFI_MATERIAU')
101 #   --- RECHERCHE SI LE MATERIAU DEPEND DE LA TEMPERATURE:
102    valk = aster.getvectjev( string.ljust(MATER.nom,8)+'.'+phenom[:10]+'.VALK')
103    valk = [x.strip() for x in valk]
104    valr = aster.getvectjev( string.ljust(MATER.nom,8)+'.'+phenom[:10]+'.VALR')
105    dicmat=dict(zip(valk,valr))
106 #   --- PROPRIETES MATERIAUX DEPENDANTES DE LA TEMPERATURE
107    Tempe3D = False
108    if dicmat.has_key('TEMP_DEF') and FOND_FISS and RESULTAT : 
109 # on recupere juste le nom du resultat thermique
110       ndim   = 3
111       Lchar = aster.getvectjev(string.ljust(RESULTAT.nom,8)+'.0000.EXCIT.LCHA        ')
112       for i in range(len(Lchar)):
113          resuth = aster.getvectjev(Lchar[i][0:8]+'.CHME.TEMPE.TEMP        ')
114          if resuth !=None :
115             Tempe3D = True
116             break
117    elif dicmat.has_key('TEMP_DEF') and not Tempe3D :
118       message = 'LES PROPRIETES MATERIAUX, NECESSAIRES AUX CALCULS \n'
119       message = message +'DES COEFFICIENTS D INTENSITE DES CONTRAINTES, ONT ETE OBTENUES A LA\n'
120       message = message +'TEMPERATURE DE REFERENCE DU MATERIAU ET NON A LA TEMPERATURE CALCULEE.'
121       UTMESS('A', macro, message)
122       nompar = ('TEMP',)
123       valpar = (dicmat['TEMP_DEF'],)
124       nomres=['E','NU']
125       valres,codret = MATER.RCVALE('ELAS',nompar,valpar,nomres,'F')
126       e = valres[0]
127       nu = valres[1]
128       
129
130 #   --- PROPRIETES MATERIAUX INDEPENDANTES DE LA TEMPERATURE
131    else :
132       e  = dicmat['E']
133       nu = dicmat['NU']  
134    
135    if not Tempe3D :
136       coefd3 = 0.
137       coefd  = e * sqrt(2.*pi)
138       unmnu2 = 1. - nu**2
139       unpnu  = 1. + nu
140       if MODELISATION=='3D' :
141          UTMESS('I', macro, 'l operateur CALC_G -option CALC_K_G- calcule plus precisement les K1 K2 K3')
142          ndim   = 3
143          coefd  = coefd      / ( 8.0 * unmnu2 )
144          coefd3 = e*sqrt(2*pi) / ( 8.0 * unpnu )
145          coefg  = unmnu2 / e
146          coefg3 = unpnu  / e
147       elif MODELISATION=='AXIS' :
148          ndim   = 2
149          coefd  = coefd  / ( 8. * unmnu2 )
150          coefg  = unmnu2 / e
151          coefg3 = unpnu  / e
152       elif MODELISATION=='D_PLAN' :
153          UTMESS('I', macro, 'l operateur CALC_G -option CALC_K_G- calcule plus precisement les K1 K2')
154          ndim   = 2
155          coefd  = coefd / ( 8. * unmnu2 )
156          coefg  = unmnu2 / e
157          coefg3 = unpnu  / e
158       elif MODELISATION=='C_PLAN' :
159          UTMESS('I', macro, 'l operateur CALC_G -option CALC_K_G- calcule plus precisement les K1 K2')
160          ndim   = 2
161          coefd  = coefd / 8.
162          coefg  = 1. / e
163          coefg3 = unpnu / e
164       else :
165          UTMESS('F', macro, 'modélisation non implantée')
166
167
168 #   ------------------------------------------------------------------
169 #                        CAS FOND_FISS
170 #   ------------------------------------------------------------------
171    if FOND_FISS : 
172       MAILLAGE = args['MAILLAGE']
173       NOEUD = args['NOEUD']
174       SANS_NOEUD = args['SANS_NOEUD']
175       GROUP_NO = args['GROUP_NO']
176       SANS_GROUP_NO = args['SANS_GROUP_NO']
177       TOUT = args['TOUT']
178       TYPE_MAILLAGE = args['TYPE_MAILLAGE']
179       NB_NOEUD_COUPE = args['NB_NOEUD_COUPE']
180       LNOFO = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.FOND      .NOEU        ')
181       RECOL = False
182 # Cas double fond de fissure : par convention les noeuds sont ceux de fond_inf
183       if LNOFO==None :
184          RECOL = True
185          LNOFO = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.FOND_INF  .NOEU        ')
186          if LNOFO==None : UTMESS('F', macro, 'PROBLEME A LA RECUPERATION DES NOEUDS DU FOND DE FISSURE \n')
187       LNOFO = map(string.rstrip,LNOFO)
188       Nbfond = len(LNOFO)
189
190       if MODELISATION=='3D' :
191 #   ----------Mots cles TOUT, NOEUD, SANS_NOEUD -------------
192         Typ = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.FOND      .TYPE        ')
193         if (Typ[0]=='SEG2    ') or (Typ[0]=='SEG3    ' and TOUT=='OUI') :
194            pas = 1
195         elif (Typ[0]=='SEG3    ') : 
196            pas = 2
197         else :
198            UTMESS('F', macro, 'TYPE DE MAILLES DU FOND DE FISSURE NON DEFINI')
199 ####
200         NO_SANS = []
201         NO_AVEC = []
202         if GROUP_NO!=None :
203           collgrno=aster.getcolljev(string.ljust(MAILLAGE.nom,8)+'.GROUPENO')
204           cnom = aster.getvectjev(string.ljust(MAILLAGE.nom,8)+'.NOMNOE')
205           if type(GROUP_NO) not in EnumTypes : GROUP_NO = (GROUP_NO,)
206           for m in range(len(GROUP_NO)) :
207             ngrno=GROUP_NO[m].ljust(8).upper()
208             if ngrno not in collgrno.keys() :
209               UTMESS('F', macro, "LE GROUP_NO "+ngrno+" N EST PAS DANS LE MAILLAGE")
210             for i in range(len(collgrno[ngrno])) : NO_AVEC.append(cnom[collgrno[ngrno][i]-1])
211           NO_AVEC= map(string.rstrip,NO_AVEC)
212         if NOEUD!=None : 
213           if type(NOEUD) not in EnumTypes : NO_AVEC = (NOEUD,)
214           else : NO_AVEC = NOEUD
215         if SANS_GROUP_NO!=None :
216           collgrno=aster.getcolljev(string.ljust(MAILLAGE.nom,8)+'.GROUPENO')
217           cnom = aster.getvectjev(string.ljust(MAILLAGE.nom,8)+'.NOMNOE')
218           if type(SANS_GROUP_NO) not in EnumTypes : SANS_GROUP_NO = (SANS_GROUP_NO,)
219           for m in range(len(SANS_GROUP_NO)) :
220             ngrno=SANS_GROUP_NO[m].ljust(8).upper()
221             if ngrno not in collgrno.keys() :
222               UTMESS('F', macro, "LE GROUP_NO "+ngrno+" N EST PAS DANS LE MAILLAGE")
223             for i in range(len(collgrno[ngrno])) : NO_SANS.append(cnom[collgrno[ngrno][i]-1])
224           NO_SANS= map(string.rstrip,NO_SANS)
225         if SANS_NOEUD!=None : 
226           if type(SANS_NOEUD) not in EnumTypes : NO_SANS = (SANS_NOEUD,)
227           else : NO_SANS = SANS_NOEUD
228 # Creation de la liste des noeuds du fond a traiter : Lnf1
229         Lnf1 = []
230         Nbf1 = 0
231         if len(NO_AVEC)!=0 :
232           for i in range(len(NO_AVEC)) :
233             if NO_AVEC[i] in LNOFO : 
234               Lnf1.append(NO_AVEC[i])
235               Nbf1 = Nbf1 +1
236             else : 
237               UTMESS('F', macro, 'LE NOEUD %s N APPARTIENT PAS AU FOND DE FISSURE'%NO_AVEC[i])
238         else :
239            for i in range(0,Nbfond,pas) :
240               if not (LNOFO[i] in NO_SANS) :
241                  Lnf1.append(LNOFO[i])
242                  Nbf1 = Nbf1 +1
243       else :
244         Lnf1 = LNOFO
245         Nbf1 = 1
246         
247 ##### Cas maillage libre###########
248 # creation des directions normales et macr_lign_coup
249       if TYPE_MAILLAGE =='LIBRE':
250         if not RESULTAT : UTMESS('F', macro, 'MOT CLE RESULTAT OBLIGATOIRE POUR TYPE_MAILLAGE = LIBRE')
251         Lnofon = Lnf1
252         Nbnofo = Nbf1
253         ListmaS = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.LEVRESUP  .MAIL        ')
254         if SYME_CHAR=='SANS':
255           ListmaI = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.LEVREINF  .MAIL        ')
256         __NCOFON=POST_RELEVE_T(ACTION=_F(INTITULE='Tab pour coordonnees noeuds du fond',
257                                             NOEUD=LNOFO,
258                                             RESULTAT=RESULTAT,
259                                             NOM_CHAM='DEPL',NUME_ORDRE=1,NOM_CMP=('DX',),
260                                             OPERATION='EXTRACTION',),);
261         tcoorf=__NCOFON.EXTR_TABLE()
262         DETRUIRE(CONCEPT=_F(NOM=__NCOFON),INFO=1) 
263         nbt = len(tcoorf['NOEUD'].values()['NOEUD'])
264         xs=array(tcoorf['COOR_X'].values()['COOR_X'][:nbt],Float)
265         ys=array(tcoorf['COOR_Y'].values()['COOR_Y'][:nbt],Float)
266         if ndim==2 : zs=Numeric.zeros(nbt,Float)
267         elif ndim==3 : zs=array(tcoorf['COOR_Z'].values()['COOR_Z'][:nbt],Float)
268         ns = tcoorf['NOEUD'].values()['NOEUD'][:nbt]
269         ns = map(string.rstrip,ns)
270         l_coorf =  [[ns[i],xs[i],ys[i],zs[i]] for i in range(0,nbt)]
271         l_coorf = [(i[0],i[1:]) for i in l_coorf]
272         d_coorf = dict(l_coorf) 
273 # Coordonnee d un pt quelconque des levres pr determination sens de propagation
274         cmail=aster.getvectjev(string.ljust(MAILLAGE.nom,8)+'.NOMMAI')
275         for i in range(len(cmail)) :
276             if cmail[i] == ListmaS[0] : break
277         colcnx=aster.getcolljev(string.ljust(MAILLAGE.nom,8)+'.CONNEX')
278         cnom = aster.getvectjev(string.ljust(MAILLAGE.nom,8)+'.NOMNOE')
279         NO_TMP = []
280         for k in range(len(colcnx[i+1])) : NO_TMP.append(cnom[colcnx[i+1][k]-1])
281         __NCOLEV=POST_RELEVE_T(ACTION=_F(INTITULE='Tab pour coordonnees pt levre',
282                                           NOEUD = NO_TMP,
283                                            RESULTAT=RESULTAT,
284                                            NOM_CHAM='DEPL',NUME_ORDRE=1,NOM_CMP=('DX',),
285                                            OPERATION='EXTRACTION',),);
286         tcoorl=__NCOLEV.EXTR_TABLE()
287         DETRUIRE(CONCEPT=_F(NOM=__NCOLEV),INFO=1) 
288         nbt = len(tcoorl['NOEUD'].values()['NOEUD'])
289         xl=moy(tcoorl['COOR_X'].values()['COOR_X'][:nbt])
290         yl=moy(tcoorl['COOR_Y'].values()['COOR_Y'][:nbt])
291         zl=moy(tcoorl['COOR_Z'].values()['COOR_Z'][:nbt])
292         Plev = array([xl, yl, zl])
293 # Calcul des normales a chaque noeud du fond
294         v1 =  array(VECT_K1)
295         VN = [None]*Nbfond
296         absfon = [0,]
297         if MODELISATION=='3D' :
298           DTANOR = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.DTAN_ORIGINE')
299           Pfon2 = array([d_coorf[LNOFO[0]][0],d_coorf[LNOFO[0]][1],d_coorf[LNOFO[0]][2]])
300           VLori = Pfon2 - Plev
301           if DTANOR != None :
302             VN[0] = array(DTANOR)
303           else :
304             Pfon3 = array([d_coorf[LNOFO[1]][0],d_coorf[LNOFO[1]][1],d_coorf[LNOFO[1]][2]])
305             VT = (Pfon3 - Pfon2)/sqrt(dot(transpose(Pfon3-Pfon2),Pfon3-Pfon2))
306             VN[0] = array(cross_product(VT,v1))
307           for i in range(1,Nbfond-1):
308             Pfon1 = array([d_coorf[LNOFO[i-1]][0],d_coorf[LNOFO[i-1]][1],d_coorf[LNOFO[i-1]][2]])
309             Pfon2 = array([d_coorf[LNOFO[i]][0],d_coorf[LNOFO[i]][1],d_coorf[LNOFO[i]][2]])
310             Pfon3 = array([d_coorf[LNOFO[i+1]][0],d_coorf[LNOFO[i+1]][1],d_coorf[LNOFO[i+1]][2]])
311             absf = sqrt(dot(transpose(Pfon1-Pfon2),Pfon1-Pfon2)) + absfon[i-1]
312             absfon.append(absf)
313             VT = (Pfon3 - Pfon2)/sqrt(dot(transpose(Pfon3-Pfon2),Pfon3-Pfon2))
314             VT = VT+(Pfon2 - Pfon1)/sqrt(dot(transpose(Pfon2-Pfon1),Pfon2-Pfon1))
315             VN[i] = array(cross_product(VT,v1)) 
316             VN[i] = VN[i]/sqrt(dot(transpose(VN[i]),VN[i]))
317           i = Nbfond-1
318           Pfon1 = array([d_coorf[LNOFO[i-1]][0],d_coorf[LNOFO[i-1]][1],d_coorf[LNOFO[i-1]][2]])
319           Pfon2 = array([d_coorf[LNOFO[i]][0],d_coorf[LNOFO[i]][1],d_coorf[LNOFO[i]][2]])
320           VLextr = Pfon2 - Plev
321           absf = sqrt(dot(transpose(Pfon1-Pfon2),Pfon1-Pfon2)) + absfon[i-1]
322           absfon.append(absf)
323           DTANEX = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.DTAN_EXTREMITE')
324           if DTANEX != None :
325             VN[i] = array(DTANEX)
326           else :
327             VT = (Pfon2 - Pfon1)/sqrt(dot(transpose(Pfon2-Pfon1),Pfon2-Pfon1))
328             VN[i] = array(cross_product(VT,v1))
329           dicoF = dict([(LNOFO[i],absfon[i]) for i in range(Nbfond)])  
330           dicVN = dict([(LNOFO[i],VN[i]) for i in range(Nbfond)])
331 #Sens de la tangente       
332           v = cross_product(VLori,VLextr)
333           sens = sign(dot(transpose(v),v1))
334 #Cas 2D              
335         if MODELISATION!='3D' :
336           DTANOR = False
337           DTANEX = False
338           VT = array([0.,0.,1.])
339           VN = array(cross_product(v1,VT))
340           dicVN = dict([(LNOFO[0],VN)])
341           Pfon = array([d_coorf[LNOFO[0]][0],d_coorf[LNOFO[0]][1],d_coorf[LNOFO[0]][2]])
342           VLori = Pfon - Plev
343           sens = sign(dot(transpose(VN),VLori))
344 #Extraction dep sup/inf sur les normales          
345         TlibS = [None]*Nbf1
346         TlibI = [None]*Nbf1
347         if NB_NOEUD_COUPE < 3 : 
348           message = 'LE NOMBRE DE NOEUDS NB_NOEUD_COUPE DOIT ETRE SUPERIEUR A 3 : ON PREND LA VALEUR PAR DEFAUT'
349           UTMESS('A', macro, message)
350           NB_NOEUD_COUPE = 5
351         MOD = aster.getvectjev(string.ljust(RESULTAT.nom,19)+'.MODL        ')
352         if MOD==None : UTMESS('F', macro, 'PROBLEME A LA RECUPERATION DU MODELE DANS LA SD RESULTAT FOURNIE')
353         MOD = map(string.rstrip,MOD)
354         MODEL = self.jdc.sds_dict[MOD[0]]
355         for i in range(Nbf1):
356           Porig = array(d_coorf[Lnf1[i]] )
357           if Lnf1[i]==LNOFO[0] and DTANOR : Pextr = Porig - ABSC_CURV_MAXI*dicVN[Lnf1[i]]
358           elif Lnf1[i]==LNOFO[Nbfond-1] and DTANEX : Pextr = Porig - ABSC_CURV_MAXI*dicVN[Lnf1[i]]
359           else : Pextr = Porig - ABSC_CURV_MAXI*dicVN[Lnf1[i]]*sens
360           TlibS[i] = MACR_LIGN_COUPE(RESULTAT=RESULTAT,
361                 NOM_CHAM='DEPL',MODELE=MODEL, MAILLE = ListmaS,
362                 LIGN_COUPE=_F(NB_POINTS=NB_NOEUD_COUPE,COOR_ORIG=(Porig[0],Porig[1],Porig[2],),
363                                           COOR_EXTR=(Pextr[0],Pextr[1],Pextr[2]),),);
364           if SYME_CHAR=='SANS':
365             TlibI[i] = MACR_LIGN_COUPE(RESULTAT=RESULTAT,
366                   NOM_CHAM='DEPL',MODELE=MODEL, MAILLE = ListmaI,
367                 LIGN_COUPE=_F(NB_POINTS=NB_NOEUD_COUPE,COOR_ORIG=(Porig[0],Porig[1],Porig[2],),
368                                           COOR_EXTR=(Pextr[0],Pextr[1],Pextr[2]),),);
369
370
371 ##### Cas maillage regle###########
372       else:
373 #   ---------- Dictionnaires des levres  -------------  
374         NnormS = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.SUPNORM   .NOEU        ')
375         if NnormS==None : 
376           message= 'PROBLEME A LA RECUPERATION DES NOEUDS DE LA LEVRE SUP : VERIFIER '
377           message=message+'QUE LE MOT CLE LEVRE_SUP EST BIEN RENSEIGNE DANS DEFI_FOND_FISS\n'
378           UTMESS('F', macro, message)
379         NnormS = map(string.rstrip,NnormS)
380         if LNOFO[0]==LNOFO[-1] and MODELISATION=='3D' : Nbfond=Nbfond-1  # Cas fond de fissure ferme
381         NnormS = [[LNOFO[i],NnormS[i*20:(i+1)*20]] for i in range(0,Nbfond)]
382         NnormS = [(i[0],i[1][0:]) for i in NnormS]
383         dicoS = dict(NnormS)
384         if SYME_CHAR=='SANS':
385            NnormI = aster.getvectjev(string.ljust(FOND_FISS.nom,8)+'.INFNORM   .NOEU        ')
386            if NnormI==None : 
387              message= 'PROBLEME A LA RECUPERATION DES NOEUDS DE LA LEVRE INF : VERIFIER '
388              message=message+'QUE LE MOT CLE LEVRE_INF EST BIEN RENSEIGNE DANS DEFI_FOND_FISS\n'
389              UTMESS('F', macro, message)
390            NnormI = map(string.rstrip,NnormI)
391            NnormI = [[LNOFO[i],NnormI[i*20:(i+1)*20]] for i in range(0,Nbfond)]
392            NnormI = [(i[0],i[1][0:]) for i in NnormI]
393            dicoI = dict(NnormI)
394
395 #   ---------- Dictionnaire des coordonnees  -------------  
396         if RESULTAT :
397           Ltot = LNOFO
398           for i in range(Nbf1) :
399             for k in range(0,20) :
400               if dicoS[Lnf1[i]][k] !='': Ltot.append(dicoS[Lnf1[i]][k])
401           if SYME_CHAR=='SANS':
402             for i in range(Nbf1) :
403               for k in range(0,20) :
404                 if dicoI[Lnf1[i]][k] !='': Ltot.append(dicoI[Lnf1[i]][k])
405           Ltot=dict([(i,0) for i in Ltot]).keys()
406           __NCOOR=POST_RELEVE_T(ACTION=_F(INTITULE='Tab pour coordonnees noeuds des levres',
407                                             NOEUD=Ltot,
408                                             RESULTAT=RESULTAT,
409                                             NOM_CHAM='DEPL',NUME_ORDRE=1,NOM_CMP=('DX',),
410                                             OPERATION='EXTRACTION',),);
411           tcoor=__NCOOR.EXTR_TABLE()
412           DETRUIRE(CONCEPT=_F(NOM=__NCOOR),INFO=1)  
413         else :  
414           if SYME_CHAR=='SANS':
415             __NCOOR=CALC_TABLE(TABLE=TABL_DEPL_SUP,
416                         ACTION=_F(OPERATION = 'COMB',NOM_PARA='NOEUD',TABLE=TABL_DEPL_INF,))
417             tcoor=__NCOOR.EXTR_TABLE()
418             DETRUIRE(CONCEPT=_F(NOM=__NCOOR),INFO=1)  
419           else :
420             tcoor=TABL_DEPL_SUP.EXTR_TABLE()
421         nbt = len(tcoor['NOEUD'].values()['NOEUD'])
422         xs=array(tcoor['COOR_X'].values()['COOR_X'][:nbt],Float)
423         ys=array(tcoor['COOR_Y'].values()['COOR_Y'][:nbt],Float)
424         if ndim==2 : zs=Numeric.zeros(nbt,Float)
425         elif ndim==3 : zs=array(tcoor['COOR_Z'].values()['COOR_Z'][:nbt],Float)
426         ns = tcoor['NOEUD'].values()['NOEUD'][:nbt]
427         ns = map(string.rstrip,ns)
428         l_coor =  [[ns[i],xs[i],ys[i],zs[i]] for i in range(0,nbt)]
429         l_coor = [(i[0],i[1:]) for i in l_coor]
430         d_coor = dict(l_coor)
431
432 #   ---------- Abscisse curviligne du fond  -------------  
433         absfon = [0,]
434         for i in range(Nbfond-1) :
435           Pfon1 = array([d_coor[LNOFO[i]][0],d_coor[LNOFO[i]][1],d_coor[LNOFO[i]][2]])
436           Pfon2 = array([d_coor[LNOFO[i+1]][0],d_coor[LNOFO[i+1]][1],d_coor[LNOFO[i+1]][2]])
437           absf = sqrt(dot(transpose(Pfon1-Pfon2),Pfon1-Pfon2)) + absfon[i]
438           absfon.append(absf)
439         dicoF = dict([(LNOFO[i],absfon[i]) for i in range(Nbfond)])
440
441      
442 # ---Noeuds LEVRE_SUP et LEVRE_INF: ABSC_CURV_MAXI et PREC_VIS_A_VIS-----
443    
444         NBTRLS = 0
445         NBTRLI = 0
446         Lnosup = [None]*Nbf1
447         Lnoinf = [None]*Nbf1
448         Nbnofo = 0
449         Lnofon = []
450         precv = PREC_VIS_A_VIS
451         if ABSC_CURV_MAXI!=None : rmax = ABSC_CURV_MAXI
452         else                   : rmax = 100
453         precn = precv * rmax
454         rmprec= rmax*(1.+precv/10.)
455         for i in range(0,Nbf1) :
456            Pfon = array([d_coor[Lnf1[i]][0],d_coor[Lnf1[i]][1],d_coor[Lnf1[i]][2]])
457            Tmpsup = []
458            Tmpinf = []
459            itots = 0
460            itoti = 0
461            NBTRLS = 0
462            NBTRLI = 0
463            for k in range(0,20) :
464               if dicoS[Lnf1[i]][k] !='':
465                  itots = itots +1
466                  Nsup =  dicoS[Lnf1[i]][k]
467                  Psup = array([d_coor[Nsup][0],d_coor[Nsup][1],d_coor[Nsup][2]])
468                  abss = sqrt(dot(transpose(Pfon-Psup),Pfon-Psup))
469                  if abss<rmprec :
470                     NBTRLS = NBTRLS +1
471                     Tmpsup.append(dicoS[Lnf1[i]][k])
472               if SYME_CHAR=='SANS':
473                  if dicoI[Lnf1[i]][k] !='':
474                     itoti = itoti +1
475                     Ninf =  dicoI[Lnf1[i]][k]
476                     Pinf = array([d_coor[Ninf][0],d_coor[Ninf][1],d_coor[Ninf][2]])
477                     absi = sqrt(dot(transpose(Pfon-Pinf),Pfon-Pinf))
478 # On verifie que les noeuds sont en vis a vis
479                     if abss<rmprec :
480                       dist = sqrt(dot(transpose(Psup-Pinf),Psup-Pinf))
481                       if dist>precn : 
482                         message= 'LES NOEUDS NE SONT PAS EN VIS-A-VIS \n'
483                         message=message+'DANS LE PLAN PERPENDICULAIRE AU NOEUD %s \n'%Lnf1[i]
484                         UTMESS('A', macro, message)
485                       else :
486                         NBTRLI = NBTRLI +1
487                         Tmpinf.append(dicoI[Lnf1[i]][k])
488 # On verifie qu il y a assez de noeuds
489            if NBTRLS < 3 : 
490               message= 'IL MANQUE DES POINTS DANS LE PLAN DEFINI PAR LA LEVRE \n'
491               message=message+'SUPERIEURE ET PERPENDICULAIRE AU FOND %s :\n'%Lnf1[i]
492               if itots<3 : message=message+' Augmenter PREC_NORM dans DEFI_FOND_FISS \n'
493               else : message=message+' Augmenter ABSC_CURV_MAXI'
494               if Lnf1[i]==LNOFO[0] or Lnf1[i]==LNOFO[-1]: message=message+' OU VERIFIER LES TANGENTES EXTREMITES'
495               UTMESS('A',macro, message)
496            elif (SYME_CHAR=='SANS') and (NBTRLI < 3) :
497               message= 'IL MANQUE DES POINTS DANS LE PLAN DEFINI PAR LA LEVRE \n'
498               message=message+'INFERIEURE ET PERPENDICULAIRE AU FOND %s :\n'%Lnf1[i]
499               if itoti<3 : message=message+' Augmenter PREC_NORM dans DEFI_FOND_FISS \n'
500               else : message=message+' Augmenter ABSC_CURV_MAXI'
501               if Lnf1[i]==LNOFO[0] or Lnf1[i]==LNOFO[-1]: message=message+' OU VERIFIER LES TANGENTES EXTREMITES'
502               UTMESS('A',macro, message)
503            else :
504               Lnosup[Nbnofo] = Tmpsup
505               if SYME_CHAR=='SANS' : Lnoinf[Nbnofo] = Tmpinf
506               Lnofon.append(Lnf1[i])
507               Nbnofo = Nbnofo+1
508         if Nbnofo == 0 :
509           message= 'CALCUL POSSIBLE POUR AUCUN NOEUD DU FOND :'
510           message=message+' VERIFIER LES DONNEES'
511           UTMESS('F',macro, message)
512                  
513    else :
514      Nbnofo = 1
515  
516 #   ----------Recuperation de la temperature au fond -------------  
517    if Tempe3D :
518       resuth = map(string.rstrip,resuth)
519       Rth = self.jdc.sds_dict[resuth[0]]
520       __TEMP=POST_RELEVE_T(ACTION=_F(INTITULE='Temperature fond de fissure',
521                                        NOEUD=Lnofon,TOUT_CMP='OUI',
522                                        RESULTAT=Rth,NOM_CHAM='TEMP',TOUT_ORDRE='OUI',
523                                        OPERATION='EXTRACTION',),);
524       tabtemp=__TEMP.EXTR_TABLE()
525       DETRUIRE(CONCEPT=_F(NOM=__TEMP),INFO=1) 
526    
527
528 #   ------------------------------------------------------------------
529 #                         BOUCLE SUR NOEUDS DU FOND
530 #   ------------------------------------------------------------------
531    for ino in range(0,Nbnofo) :
532       if FOND_FISS and INFO==2 :
533             texte="\n\n--> TRAITEMENT DU NOEUD DU FOND DE FISSURE: %s"%Lnofon[ino]
534             aster.affiche('MESSAGE',texte)
535 #   ------------------------------------------------------------------
536 #                         TABLE 'DEPSUP'
537 #   ------------------------------------------------------------------
538       if FOND_FISS : 
539          if TYPE_MAILLAGE =='LIBRE':
540             tabsup=TlibS[ino].EXTR_TABLE()
541             DETRUIRE(CONCEPT=_F(NOM=TlibS[ino]),INFO=1)
542          elif RESULTAT :
543             __TSUP=POST_RELEVE_T(ACTION=_F(INTITULE='Deplacement SUP',
544                                           NOEUD=Lnosup[ino],
545                                           RESULTAT=RESULTAT,
546                                           NOM_CHAM='DEPL',
547                                           TOUT_ORDRE='OUI',
548                                           NOM_CMP=('DX','DY','DZ',),
549                                           OPERATION='EXTRACTION',),);
550             tabsup=__TSUP.EXTR_TABLE()
551             DETRUIRE(CONCEPT=_F(NOM=__TSUP),INFO=1)      
552          else :
553             tabsup=TABL_DEPL_SUP.EXTR_TABLE()
554             veri_tab(tabsup,TABL_DEPL_SUP.nom,ndim)
555             Ls = [string.ljust(Lnosup[ino][i],8) for i in range(len(Lnosup[ino]))]
556             tabsup=tabsup.NOEUD==Ls
557       else :
558          tabsup=TABL_DEPL_SUP.EXTR_TABLE()
559          veri_tab(tabsup,TABL_DEPL_SUP.nom,ndim)
560
561 #   ------------------------------------------------------------------
562 #                          TABLE 'DEPINF'
563 #   ------------------------------------------------------------------
564       if SYME_CHAR=='SANS': 
565          if FOND_FISS : 
566             if TYPE_MAILLAGE =='LIBRE':
567                tabinf=TlibI[ino].EXTR_TABLE()
568                DETRUIRE(CONCEPT=_F(NOM=TlibI[ino]),INFO=1)
569             elif RESULTAT :
570                __TINF=POST_RELEVE_T(ACTION=_F(INTITULE='Deplacement INF',
571                                           NOEUD=Lnoinf[ino],
572                                           RESULTAT=RESULTAT,
573                                           NOM_CHAM='DEPL',
574                                           TOUT_ORDRE='OUI',
575                                           NOM_CMP=('DX','DY','DZ',),
576                                           OPERATION='EXTRACTION',),);
577                tabinf=__TINF.EXTR_TABLE()   
578                DETRUIRE(CONCEPT=_F(NOM=__TINF),INFO=1)                 
579             else :
580                tabinf=TABL_DEPL_INF.EXTR_TABLE()
581                if TABL_DEPL_INF==None : UTMESS('F', macro, 'TABL_DEPL_SUP et TABL_DEPL_INF sont obligatoires si SYME_CHAR=SANS')
582                veri_tab(tabinf,TABL_DEPL_INF.nom,ndim)
583                Li = [string.ljust(Lnoinf[ino][i],8) for i in range(len(Lnoinf[ino]))]
584                tabinf=tabinf.NOEUD==Li
585          else :
586             if TABL_DEPL_INF==None : UTMESS('F', macro, 'TABL_DEPL_SUP et TABL_DEPL_INF sont obligatoires si SYME_CHAR=SANS')
587             tabinf=TABL_DEPL_INF.EXTR_TABLE()
588             veri_tab(tabinf,TABL_DEPL_INF.nom,ndim)
589
590
591 #   ------------------------------------------------------------------
592 #               LES INSTANTS DE POST-TRAITEMENT
593 #   ------------------------------------------------------------------
594       if 'INST' in tabsup.para : 
595          l_inst=None
596          l_inst_tab=tabsup['INST'].values()['INST']
597          l_inst_tab=dict([(i,0) for i in l_inst_tab]).keys() #elimine les doublons
598          l_inst_tab.sort()
599          if LIST_ORDRE !=None or NUME_ORDRE !=None :
600            l_ord_tab = tabsup['NUME_ORDRE'].values()['NUME_ORDRE']
601            l_ord_tab.sort()
602            l_ord_tab=dict([(i,0) for i in l_ord_tab]).keys() 
603            d_ord_tab= [[l_ord_tab[i],l_inst_tab[i]] for i in range(0,len(l_ord_tab))]
604            d_ord_tab= [(i[0],i[1]) for i in d_ord_tab]
605            d_ord_tab = dict(d_ord_tab)
606            if NUME_ORDRE !=None : 
607              if type(NUME_ORDRE) not in EnumTypes : NUME_ORDRE=(NUME_ORDRE,)
608              l_ord=list(NUME_ORDRE)
609            elif LIST_ORDRE !=None : 
610               l_ord= aster.getvectjev(string.ljust(LIST_ORDRE.nom,19)+'.VALE') 
611            l_inst = []
612            for ord in l_ord :
613              if ord in l_ord_tab : l_inst.append(d_ord_tab[ord])
614              else :  
615                message ='LE NUMERO D ORDRE %i N A PAS ETE ETE TROUVE DANS LA TABLE\n'%ord 
616                UTMESS('F', macro, message)
617          if INST !=None or LIST_INST !=None :
618             CRITERE = args['CRITERE']
619             PRECISION = args['PRECISION']
620          else :
621             l_inst=l_inst_tab
622             PRECISION = 1.E-6
623             CRITERE='ABSOLU'
624          if        INST !=None : 
625             if type(INST) not in EnumTypes : INST=(INST,)
626             l_inst=list(INST)
627          elif LIST_INST !=None : l_inst=LIST_INST.Valeurs()
628          if      l_inst !=None :
629            for inst in l_inst  :
630                if CRITERE=='RELATIF' and inst!=0.: match=[x for x in l_inst_tab if abs((inst-x)/inst)<PRECISION]
631                else                              : match=[x for x in l_inst_tab if abs(inst-x)<PRECISION]
632                if len(match)==0 : 
633                  message = 'PAS D INSTANT TROUVE DANS LA TABLE POUR L INSTANT %f\n'%inst
634                  UTMESS('F', macro, message)
635                if len(match)>=2 :
636                  message = 'PLUSIEURS INSTANTS TROUVES DANS LA TABLE POUR L INSTANT %f\n'%inst 
637                  UTMESS('F', macro, message)
638          else                  :
639            l_inst  = l_inst_tab
640       else :
641          l_inst  = [None,]
642    
643 #   ------------------------------------------------------------------
644 #                         BOUCLE SUR LES INSTANTS
645 #   ------------------------------------------------------------------
646       for iord in range(len(l_inst)) :
647         inst=l_inst[iord]
648         if INFO==2 and inst!=None:
649             texte="==> INSTANT: %f"%inst
650             aster.affiche('MESSAGE',texte)
651         if inst!=None:
652            if PRECISION == None : PRECISION = 1.E-6
653            if CRITERE == None : CRITERE='ABSOLU'
654            if inst==0. :
655              tabsupi=tabsup.INST.__eq__(VALE=inst,CRITERE='ABSOLU',PRECISION=PRECISION)
656              if SYME_CHAR=='SANS': tabinfi=tabinf.INST.__eq__(VALE=inst,CRITERE='ABSOLU',PRECISION=PRECISION)
657            else :
658              tabsupi=tabsup.INST.__eq__(VALE=inst,CRITERE=CRITERE,PRECISION=PRECISION)
659              if SYME_CHAR=='SANS': tabinfi=tabinf.INST.__eq__(VALE=inst,CRITERE=CRITERE,PRECISION=PRECISION)
660         else :
661            tabsupi=tabsup
662            if SYME_CHAR=='SANS': tabinfi=tabinf
663
664 #     --- LEVRE SUP :  "ABSC_CURV" CROISSANTES, < RMAX ET DEP ---
665         abscs = getattr(tabsupi,'ABSC_CURV').values()
666         if not FOND_FISS :
667           refs=copy.copy(abscs)
668           refs.sort()
669           if refs!=abscs : UTMESS('F', macro, 'ABSC_CURV NON CROISSANTS POUR TABL_DEPL_INF')
670           if ABSC_CURV_MAXI!=None : rmax = ABSC_CURV_MAXI
671           else                    : rmax = abscs[-1]
672           precv = PREC_VIS_A_VIS
673           rmprec= rmax*(1.+precv/10.)
674           refsc=[x for x in refs if x<rmprec]
675           nbval = len(refsc)
676         else :
677           nbval=len(abscs)
678         abscs=array(abscs[:nbval])
679         coxs=array(tabsupi['COOR_X'].values()['COOR_X'][:nbval],Float)
680         coys=array(tabsupi['COOR_Y'].values()['COOR_Y'][:nbval],Float)
681         if ndim==2 :  cozs=Numeric.zeros(nbval,Float)
682         elif ndim==3 :  cozs=array(tabsupi['COOR_Z'].values()['COOR_Z'][:nbval],Float)
683
684         if FOND_FISS and not RESULTAT : #tri des noeuds avec abscisse
685           Pfon = array([d_coor[Lnofon[ino]][0],d_coor[Lnofon[ino]][1],d_coor[Lnofon[ino]][2]])
686           abscs = sqrt((coxs-Pfon[0])**2+(coys-Pfon[1])**2+(cozs-Pfon[2])**2)
687           tabsupi['Abs_fo'] = abscs
688           tabsupi.sort('Abs_fo')
689           abscs = getattr(tabsupi,'Abs_fo').values()
690           abscs=array(abscs[:nbval])
691           coxs=array(tabsupi['COOR_X'].values()['COOR_X'][:nbval],Float)
692           coys=array(tabsupi['COOR_Y'].values()['COOR_Y'][:nbval],Float)
693           if ndim==2 :  cozs=Numeric.zeros(nbval,Float)
694           elif ndim==3 :  cozs=array(tabsupi['COOR_Z'].values()['COOR_Z'][:nbval],Float)
695           
696         if FOND_FISS and INFO==2 and iord==0 and not TYPE_MAILLAGE =='LIBRE':
697           for ks in range(0,nbval) :
698             texte="NOEUD RETENU POUR LA LEVRE SUP: %s  %f"%(Lnosup[ino][ks],abscs[ks])
699             aster.affiche('MESSAGE',texte)
700         dxs=array(tabsupi['DX'].values()['DX'][:nbval],Float)
701         dys=array(tabsupi['DY'].values()['DY'][:nbval],Float)
702         if ndim==2 : dzs=Numeric.zeros(nbval,Float)
703         elif ndim==3 : dzs=array(tabsupi['DZ'].values()['DZ'][:nbval],Float)
704         
705 #     --- LEVRE INF :  "ABSC_CURV" CROISSANTES et < RMAX ---
706         if SYME_CHAR=='SANS': 
707           absci = getattr(tabinfi,'ABSC_CURV').values()
708           if not FOND_FISS :
709             refi=copy.copy(absci)
710             refi.sort()
711             if refi!=absci : UTMESS('F', macro, 'ABSC_CURV NON CROISSANTS POUR TABL_DEPL_SUP')
712             refic=[x for x in refi if x<rmprec]
713             nbvali=len(refic)
714           else :
715             nbvali=len(absci)
716           if nbvali!=nbval :
717              message= 'DIFFERENCE DE POINTS ENTRE LA LEVRE SUPERIEURE ET LA LEVRE INFERIEURE'
718              if FOND_FISS : message=message+' POUR TRAITER LE NOEUD %.s \n'%Lnofon[i]
719              message=message+' Nombre de points - levre superieure : %i\n'%len(refsc)
720              message=message+' Nombre de points - levre inferieure : %i\n'%len(refic)
721              UTMESS('A',macro, message)
722           nbval=min(nbval,nbvali)
723           absci=array(absci[:nbval])
724           coxi=array(tabinfi['COOR_X'].values()['COOR_X'][:nbval],Float)
725           coyi=array(tabinfi['COOR_Y'].values()['COOR_Y'][:nbval],Float)
726           if ndim==2 : cozi=Numeric.zeros(nbval,Float)
727           elif ndim==3 : cozi=array(tabinfi['COOR_Z'].values()['COOR_Z'][:nbval],Float)
728 #     --- ON VERIFIE QUE LES NOEUDS SONT EN VIS_A_VIS  (SYME=SANS)   ---
729           if not FOND_FISS :
730             precn = precv * rmax
731             dist=(coxs-coxi)**2+(coys-coyi)**2+(cozs-cozi)**2
732             dist=sqrt(dist)
733             for d in dist :
734                if d>precn : UTMESS('F', macro, 'LES NOEUDS NE SONT PAS EN VIS_A_VIS')
735           
736           if FOND_FISS and not RESULTAT :#tri des noeuds avec abscisse
737             Pfon = array([d_coor[Lnofon[ino]][0],d_coor[Lnofon[ino]][1],d_coor[Lnofon[ino]][2]])
738             absci = sqrt((coxi-Pfon[0])**2+(coyi-Pfon[1])**2+(cozi-Pfon[2])**2)
739             tabinfi['Abs_fo'] = absci
740             tabinfi.sort('Abs_fo')
741             absci = getattr(tabinfi,'Abs_fo').values()
742             absci=array(abscs[:nbval])
743             coxi=array(tabinfi['COOR_X'].values()['COOR_X'][:nbval],Float)
744             coyi=array(tabinfi['COOR_Y'].values()['COOR_Y'][:nbval],Float)
745             if ndim==2 :  cozi=Numeric.zeros(nbval,Float)
746             elif ndim==3 :  cozi=array(tabinfi['COOR_Z'].values()['COOR_Z'][:nbval],Float)
747
748           dxi=array(tabinfi['DX'].values()['DX'][:nbval],Float)
749           dyi=array(tabinfi['DY'].values()['DY'][:nbval],Float)
750           if ndim==2 : dzi=Numeric.zeros(nbval,Float)
751           elif ndim==3 : dzi=array(tabinfi['DZ'].values()['DZ'][:nbval],Float)
752           
753           if FOND_FISS and INFO==2 and iord==0 and not TYPE_MAILLAGE =='LIBRE':
754             for ki in range(0,nbval) :
755                texte="NOEUD RETENU POUR LA LEVRE INF: %s  %f"%(Lnoinf[ino][ki],absci[ki])
756                aster.affiche('MESSAGE',texte)
757
758 #     --- TESTS NOMBRE DE NOEUDS---
759         if nbval<3 :
760            message= 'IL FAUT AU MOINS TROIS NOEUDS DANS LE PLAN DEFINI PAR LES LEVRES ET PERPENDICULAIRE AU FOND DE FISSURE'
761            if FOND_FISS : message=message+'Noeud %.s \n'%Lnofon[ino]
762            message=message+' : augmenter ABSC_CURV_MAXI\n'
763            UTMESS('F',macro, message)
764            
765 #   ---------- CALCUL PROP. MATERIAU AVEC TEMPERATURE -----------  
766         if Tempe3D :
767            tempeno=tabtemp.NOEUD==Lnofon[ino]
768            tempeno=tempeno.INST.__eq__(VALE=inst,CRITERE='ABSOLU',PRECISION=PRECISION)
769            nompar = ('TEMP',)
770            valpar = (tempeno.TEMP.values()[0],)
771            nomres=['E','NU']
772            valres,codret = MATER.RCVALE('ELAS',nompar,valpar,nomres,'F')
773            e = valres[0]
774            nu = valres[1] 
775            coefd  = e * sqrt(2.*pi)      / ( 8.0 * (1. - nu**2))
776            coefd3 = e*sqrt(2*pi) / ( 8.0 * (1. + nu))
777            coefg  = (1. - nu**2) / e
778            coefg3 = (1. + nu)  / e
779
780 #     ------------------------------------------------------------------
781 #                           CHANGEMENT DE REPERE
782 #     ------------------------------------------------------------------
783 #
784 #       1 : VECTEUR NORMAL AU PLAN DE LA FISSURE
785 #           ORIENTE LEVRE INFERIEURE VERS LEVRE SUPERIEURE
786 #       2 : VECTEUR NORMAL AU FOND DE FISSURE EN M
787 #       3 : VECTEUR TANGENT AU FOND DE FISSURE EN M
788 #
789         if SYME_CHAR=='SANS' :
790            vo =  array([( coxs[-1]+coxi[-1] )/2.,( coys[-1]+coyi[-1] )/2.,( cozs[-1]+cozi[-1] )/2.])
791            ve =  array([( coxs[0 ]+coxi[0 ] )/2.,( coys[0 ]+coyi[0 ] )/2.,( cozs[0 ]+cozi[0 ] )/2.])
792         else :
793            vo = array([ coxs[-1], coys[-1], cozs[-1]])
794            ve = array([ coxs[0], coys[0], cozs[0]])
795         v1 =  array(VECT_K1)
796         v2 =  ve-vo
797         v2 =  v2/sqrt(v2[0]**2+v2[1]**2+v2[2]**2)
798         v1p = sum(v2*v1)
799         v1  = v1-v1p*v2
800         v1  = v1/sqrt(v1[0]**2+v1[1]**2+v1[2]**2)
801         v3  = array([v1[1]*v2[2]-v2[1]*v1[2],v1[2]*v2[0]-v2[2]*v1[0],v1[0]*v2[1]-v2[0]*v1[1]])
802         pgl  = asarray([v1,v2,v3])
803         dpls = asarray([dxs,dys,dzs])
804         dpls = matrixmultiply(pgl,dpls)
805         if SYME_CHAR=='SANS' :
806            dpli = asarray([dxi,dyi,dzi])
807            dpli = matrixmultiply(pgl,dpli)
808         else :
809            dpli = [multiply(dpls[0],-1.),dpls[1],dpls[2]]
810 #     ------------------------------------------------------------------
811 #                           CALCUL DES K1, K2, K3
812 #     ------------------------------------------------------------------
813         saut=(dpls-dpli)
814         isig=sign(transpose(resize(saut[:,-1],(nbval-1,3))))
815         isig=sign(isig+0.001)
816         saut=saut*array([[coefd]*nbval,[coefd]*nbval,[coefd3]*nbval])
817         saut=saut**2
818         ksig = isig[:,1]
819         ksig = array([ksig,ksig])
820         ksig = transpose(ksig)
821         kgsig=resize(ksig,(1,6))[0]
822         if INFO==2 :
823           mcfact=[]
824           mcfact.append(_F(PARA='ABSC_CURV'  ,LISTE_R=abscs.tolist() ))
825           mcfact.append(_F(PARA='DEPL_SUP_DX',LISTE_R=dpls[0].tolist() ))
826           mcfact.append(_F(PARA='DEPL_INF_DX',LISTE_R=dpli[0].tolist() ))
827           mcfact.append(_F(PARA='SAUT_DX'    ,LISTE_R=saut[0].tolist() ))
828           mcfact.append(_F(PARA='DEPL_SUP_DY',LISTE_R=dpls[1].tolist() ))
829           mcfact.append(_F(PARA='DEPL_INF_DY',LISTE_R=dpli[1].tolist() ))
830           mcfact.append(_F(PARA='SAUT_DY'    ,LISTE_R=saut[1].tolist() ))
831           if ndim==3 :
832             mcfact.append(_F(PARA='DEPL_SUP_DZ',LISTE_R=dpls[2].tolist() ))
833             mcfact.append(_F(PARA='DEPL_INF_DZ',LISTE_R=dpli[2].tolist() ))
834             mcfact.append(_F(PARA='SAUT_DZ'    ,LISTE_R=saut[2].tolist() ))
835           __resu0=CREA_TABLE(LISTE=mcfact,TITRE='--> SAUTS')
836           aster.affiche('MESSAGE',__resu0.EXTR_TABLE().__repr__())
837           DETRUIRE(CONCEPT=_F(NOM=__resu0),INFO=1)
838 #     ------------------------------------------------------------------
839 #                           --- METHODE 1 ---
840 #     ------------------------------------------------------------------
841         x1 = abscs[1:-1]
842         x2 = abscs[2:  ]
843         y1 = saut[:,1:-1]/x1
844         y2 = saut[:,2:  ]/x2
845         k  = abs(y1-x1*(y2-y1)/(x2-x1))
846         g  = coefg*(k[0]+k[1])+coefg3*k[2]
847         kg1 = [max(k[0]),min(k[0]),max(k[1]),min(k[1]),max(k[2]),min(k[2])]
848         kg1 = sqrt(kg1)*kgsig
849         kg1=Numeric.concatenate([kg1,[max(g),min(g)]])
850         vk  = sqrt(k)*isig[:,:-1]
851         if INFO==2 :
852           mcfact=[]
853           mcfact.append(_F(PARA='ABSC_CURV_1' ,LISTE_R=x1.tolist() ))
854           mcfact.append(_F(PARA='ABSC_CURV_2' ,LISTE_R=x2.tolist() ))
855           mcfact.append(_F(PARA='K1'          ,LISTE_R=vk[0].tolist() ))
856           mcfact.append(_F(PARA='K2'          ,LISTE_R=vk[1].tolist() ))
857           if ndim==3 :
858             mcfact.append(_F(PARA='K3'        ,LISTE_R=vk[2].tolist() ))
859           mcfact.append(_F(PARA='G'           ,LISTE_R=g.tolist() ))
860           __resu1=CREA_TABLE(LISTE=mcfact,TITRE='--> METHODE 1')
861           aster.affiche('MESSAGE',__resu1.EXTR_TABLE().__repr__())
862           DETRUIRE(CONCEPT=_F(NOM=__resu1),INFO=1)
863 #     ------------------------------------------------------------------
864 #                           --- METHODE 2 ---
865 #     ------------------------------------------------------------------
866         x1 = abscs[1: ]
867         y1 = saut[:,1:]
868         k  = abs(y1/x1)
869         g  = coefg*(k[0]+k[1])+coefg3*k[2]
870         kg2= [max(k[0]),min(k[0]),max(k[1]),min(k[1]),max(k[2]),min(k[2])]
871         kg2 = sqrt(kg2)*kgsig
872         kg2=Numeric.concatenate([kg2,[max(g),min(g)]])
873         vk = sqrt(k)*isig
874         if INFO==2 :
875           mcfact=[]
876           mcfact.append(_F(PARA='ABSC_CURV' ,LISTE_R=x1.tolist() ))
877           mcfact.append(_F(PARA='K1'        ,LISTE_R=vk[0].tolist() ))
878           mcfact.append(_F(PARA='K2'        ,LISTE_R=vk[1].tolist() ))
879           if ndim==3 :
880             mcfact.append(_F(PARA='K3'      ,LISTE_R=vk[2].tolist() ))
881           mcfact.append(_F(PARA='G'         ,LISTE_R=g.tolist() ))
882           __resu2=CREA_TABLE(LISTE=mcfact,TITRE='--> METHODE 2')
883           aster.affiche('MESSAGE',__resu2.EXTR_TABLE().__repr__())
884           DETRUIRE(CONCEPT=_F(NOM=__resu2),INFO=1)
885 #     ------------------------------------------------------------------
886 #                           --- METHODE 3 ---
887 #     ------------------------------------------------------------------
888         x1 = abscs[:-1]
889         x2 = abscs[1: ]
890         y1 = saut[:,:-1]
891         y2 = saut[:,1: ]
892         k  = (sqrt(y2)*sqrt(x2)+sqrt(y1)*sqrt(x1))*(x2-x1)
893         k  = Numeric.sum(transpose(k))
894         de = abscs[-1]
895         vk = (k/de**2)*isig[:,0]
896         g  = coefg*(vk[0]**2+vk[1]**2)+coefg3*vk[2]**2
897         kg3=Numeric.concatenate([[vk[0]]*2,[vk[1]]*2,[vk[2]]*2,[g]*2])
898         if INFO==2 :
899           mcfact=[]
900           mcfact.append(_F(PARA='K1'        ,LISTE_R=vk[0] ))
901           mcfact.append(_F(PARA='K2'        ,LISTE_R=vk[1] ))
902           if ndim==3 :
903             mcfact.append(_F(PARA='K3'      ,LISTE_R=vk[2] ))
904           mcfact.append(_F(PARA='G'         ,LISTE_R=g ))
905           __resu3=CREA_TABLE(LISTE=mcfact,TITRE='--> METHODE 3')
906           aster.affiche('MESSAGE',__resu3.EXTR_TABLE().__repr__())
907           DETRUIRE(CONCEPT=_F(NOM=__resu3),INFO=1)
908 #     ------------------------------------------------------------------
909 #                           CREATION DE LA TABLE 
910 #     ------------------------------------------------------------------
911         kg=array([kg1,kg2,kg3])
912         kg=transpose(kg)
913         mcfact=[]
914         if TITRE != None :
915           titre = TITRE
916         else :
917           v = aster.__version__
918           titre = 'ASTER %s - CONCEPT CALCULE PAR POST_K1_K2_K3 LE &DATE A &HEURE \n'%v
919         if FOND_FISS and MODELISATION=='3D': 
920           mcfact.append(_F(PARA='NOEUD_FOND',LISTE_K=[Lnofon[ino],]*3))
921           mcfact.append(_F(PARA='ABSC_CURV',LISTE_R=[dicoF[Lnofon[ino]]]*3))
922         mcfact.append(_F(PARA='METHODE',LISTE_I=(1,2,3)))
923         mcfact.append(_F(PARA='K1_MAX' ,LISTE_R=kg[0].tolist() ))
924         mcfact.append(_F(PARA='K1_MIN' ,LISTE_R=kg[1].tolist() ))
925         mcfact.append(_F(PARA='K2_MAX' ,LISTE_R=kg[2].tolist() ))
926         mcfact.append(_F(PARA='K2_MIN' ,LISTE_R=kg[3].tolist() ))
927         if ndim==3 :
928           mcfact.append(_F(PARA='K3_MAX' ,LISTE_R=kg[4].tolist() ))
929           mcfact.append(_F(PARA='K3_MIN' ,LISTE_R=kg[5].tolist() ))
930         mcfact.append(_F(PARA='G_MAX'  ,LISTE_R=kg[6].tolist() ))
931         mcfact.append(_F(PARA='G_MIN'  ,LISTE_R=kg[7].tolist() ))
932         if  (ino==0 and iord==0) and inst==None :
933            tabout=CREA_TABLE(LISTE=mcfact,TITRE = titre)
934         elif iord==0 and ino==0 and inst!=None :
935            mcfact=[_F(PARA='INST'  ,LISTE_R=[inst,]*3      )]+mcfact
936            tabout=CREA_TABLE(LISTE=mcfact,TITRE = titre)
937         else :
938            if inst!=None : mcfact=[_F(PARA='INST'  ,LISTE_R=[inst,]*3     )]+mcfact
939            __tabi=CREA_TABLE(LISTE=mcfact,)
940            npara = ['K1_MAX','METHODE']
941            if inst!=None : npara.append('INST')
942            if FOND_FISS and MODELISATION=='3D' : npara.append('NOEUD_FOND')
943            tabout=CALC_TABLE(reuse=tabout,TABLE=tabout,TITRE = titre,
944                               ACTION=_F(OPERATION = 'COMB',NOM_PARA=npara,TABLE=__tabi,))
945
946 # Tri pour conserver le meme ordre que operateur initial en fortran
947    if len(l_inst)!=1 and FOND_FISS and MODELISATION=='3D':
948       tabout=CALC_TABLE(reuse=tabout,TABLE=tabout,
949                       ACTION=_F(OPERATION = 'TRI',NOM_PARA=('INST','ABSC_CURV','METHODE'),ORDRE='CROISSANT'))
950
951    return ier
952