Télécharger la plaquette Cast3M


fr en

Télécharger quahm2.eso

Retour à la liste

Numérotation des lignes : 

quahm2

  1. C QUAHM2    SOURCE    CHAT      05/01/13    02:40:56     5004
  2. C QUAHM2    SOURCE    AM1       95/11/24    22:54:38     1918
  3.       SUBROUTINE QUAHM2(IGAU,ITEL,MFR,NBNO,XEL,SHPTOT,SHP,IFOU,
  4.      #        NHARM,VKL12,VKL23,VKL33,POIGAU,ISDJC,LRE,REL,IRET)
  5. C=======================================================================
  6. C
  7. C    CALCULE LES TERMES  EN  P * PI ,PI * (UR,RT) ,(UR,RT) *(UR,RT)
  8. C        (UT,RR) * (UT,RR) , PI * (UT,RR)     DE   LA     MATRICE
  9. C           MASSE  DANS  LE  CAS  AXISYMETRIQUE OU FOURIER  POUR
  10. C                  LA FORMULATION (37) HOMOGENE
  11. C=======================================================================
  12. C  INPUT
  13. C     IGAU=NUMERO DU POINT DE GAUSS
  14. C     ITEL=NUMERO DE L ELEMENT DANS NOMTP
  15. C     MFR =NUMERO DE LA FORMULATION
  16. C     NBNO=NOMBRE DE NOEUDS
  17. C     XEL =COORDONNEES  DE L ELEMENT
  18. C     IFOU=IFOUR DE CCOPTIO
  19. C     NHARM=NUMERO DU MODE DE FOURIER
  20. C     VKL12=-((COEFPI*COEFPR)/(RHOF*C**2))*SFLU/SCEL
  21. C     VKL23=(BET11+BET22)*COEFPI/(2.*SCEL)
  22. C     VKL33=(RHOS*2.+RHOF*(BET11+BET22))/SCEL
  23. C     POIGAU=MINTE.POIGAU(IGAU)
  24. C     LRE =NOMBRE DE D.D.L DE LA MATRICE  DE  RIGIDITE
  25. C     SHPTOT(6,NBNO,NBGAU)=FONCTIONS DE FORMES ET DERIVEES
  26. C     ISDJC = INDICATEUR SUR LE SIGNE DU JACOBIEN
  27. C  ZONE DE TRAVAIL
  28. C     SHP(5,NBNO)=TABLEAU DE TRAVAIL
  29. C OUTPUT
  30. C     ISDJC = INDICATEUR SUR LE SIGNE DU JACOBIEN
  31. C     REL=MATRICE DE MASSE
  32. C     IRET:INDICATEUR = 1 : SUCCES
  33. C                       0 : ECHEC  (ELEMENT MELE INCOMPATIBLE )
  34. C                       2 : ECHEC  (JACOBIEN  NUL )
  35. C                       3 :ECHEC   (ROUTINE N EST VALABLE QU
  36. C                                          EN AXISYMETRIQUE OU FOURIER )
  37. C=======================================================================
  38.       IMPLICIT INTEGER(I-N)
  39.       IMPLICIT REAL*8(A-H,O-Z)
  40.       DIMENSION XEL(3,*),SHP(6,*),SHPTOT(6,NBNO,*),REL(LRE,*)
  41.       IF (ITEL.EQ.126)              GOTO 10
  42. C
  43. C     ERREUR : TYPE D' ELEMENT  INCOMPATIBLE AVEC LA FORMULATION
  44. C
  45.       IRET = 0
  46.       GOTO 666
  47.   10  CONTINUE
  48.       IF (IFOU.EQ.0.OR.IFOU.EQ.1)  GOTO 11
  49. C
  50. C     MESSAGE D ERREUR :  ROUTINE  N  EST  VALABLE  QU  EN FOURIER
  51. C                     OU   EN   AXISYMETRIQUE
  52. C
  53.       IRET = 3
  54.       GOTO 666
  55.   11  CONTINUE
  56. C
  57. C     ELEMENTS HOMOGENEISES QUAH   EN AXISYMETRIE OU EN FOURIER
  58. C          NBDL = LRE/NBNO   NOMBRE DE D.D.L PAR NOEUD
  59. C
  60.       NBDL = LRE /NBNO
  61. C
  62. C     SHP(1,I) : FONCTION DE FORME
  63. C     SHP(2,I) : DERIVEE % R DE LA FONCTION DE FORME
  64. C     SHP(3,I) : DERIVEE % Z DE LA FONCTION DE FORME
  65. C
  66.       DO 101 NP=1,NBNO
  67.       SHP(1,NP)=SHPTOT(1,NP,IGAU)
  68.       SHP(2,NP)=SHPTOT(2,NP,IGAU)
  69.       SHP(3,NP)=SHPTOT(3,NP,IGAU)
  70.   101 CONTINUE
  71. C
  72.       CALL DEVOLU(XEL,SHP,MFR,NBNO,IFOU,NHARM,2,1.D0,RR,DJAC)
  73.       IF (DJAC.EQ.0.) GOTO 667
  74.       IF (DJAC.LT.0.) ISDJC = ISDJC + 1
  75. C
  76. C     FONCTIONS DE FORME
  77. C
  78.       SHP(4,1)=SHP(1,1) + SHP(1,4)
  79.       SHP(4,2)=SHP(1,2) + SHP(1,3)
  80.       SHP(4,3)=0.D0
  81.       SHP(4,4)=0.D0
  82. C
  83.       XZH=SHP(1,3) + SHP(1,4)
  84. C
  85.       A1=MIN(XEL(2,1),XEL(2,2))
  86.       A2=MIN(XEL(2,1),XEL(2,3))
  87.       A3=MIN(XEL(2,1),XEL(2,4))
  88.       A4=MIN(A1,A2)
  89.       A5=MIN(A1,A3)
  90.       Z1=MIN(A4,A5)
  91. C
  92.       B1=MAX(XEL(2,1),XEL(2,2))
  93.       B2=MAX(XEL(2,1),XEL(2,3))
  94.       B3=MAX(XEL(2,1),XEL(2,4))
  95.       B4=MAX(B1,B2)
  96.       B5=MAX(B1,B3)
  97.       Z2=MAX(B4,B5)
  98. C
  99.       RH=Z2-Z1
  100. C
  101.       S1=SHP(1,1) + SHP(1,4)
  102.       S2=SHP(1,2) + SHP(1,3)
  103. C
  104. C     FONCTIONS DE FORME EN Z
  105. C
  106.       XZH= SHP(1,3) + SHP(1,4)
  107. C
  108.       SHP(5,1)=1.D0-(3.D0*XZH*XZH)+(2.D0*XZH*XZH*XZH)
  109.       SHP(5,2)=3.D0*(XZH**2)-(2.D0*XZH*XZH*XZH)
  110.       SHP(5,3)=(RH*XZH)*(1.D0-2.D0*XZH+XZH*XZH)
  111.       SHP(5,4)=(RH*XZH)*(XZH*XZH-XZH)
  112. C
  113. C     FONCTIONS DE FORME POUR LA FLEXION
  114. C
  115.       SHP(6,1)=S1*SHP(5,1)
  116.       SHP(6,2)=S2*SHP(5,1)
  117.       SHP(6,3)=S2*SHP(5,2)
  118.       SHP(6,4)=S1*SHP(5,2)
  119. C
  120.       SHP(6,5)=S1*SHP(5,3)
  121.       SHP(6,6)=S2*SHP(5,3)
  122.       SHP(6,7)=S2*SHP(5,4)
  123.       SHP(6,8)=S1*SHP(5,4)
  124. C
  125. C
  126.       C1=MIN(XEL(1,1),XEL(1,2))
  127.       C2=MIN(XEL(1,1),XEL(1,3))
  128.       C3=MIN(XEL(1,1),XEL(1,4))
  129.       C4=MIN(C1,C2)
  130.       C5=MIN(C1,C3)
  131.       R1=MIN(C4,C5)
  132. C
  133.       D1=MAX(XEL(1,1),XEL(1,2))
  134.       D2=MAX(XEL(1,1),XEL(1,3))
  135.       D3=MAX(XEL(1,1),XEL(1,4))
  136.       D4=MAX(D1,D2)
  137.       D5=MAX(D1,D3)
  138.       R2=MAX(D4,D5)
  139. C
  140.       DELTAR = R2 - R1
  141.       DV = (DELTAR*RH)/4.d0
  142. C
  143.       DJAC = ABS(DJAC)
  144. C
  145. C
  146. C     TERMES EN P * PI
  147. C
  148.       DJAC1 = ABS(DJAC)*POIGAU
  149.       IX1=0
  150.       IY1=0
  151.       DO   102 IX=2,LRE ,NBDL
  152.       IX1=IX1 + 1
  153.       DO   103 IY=1,IX  ,NBDL
  154.       IY1=IY1 + 1
  155.       REL(IY,IX) = REL(IY,IX) + VKL12*DJAC1*SHP(1,IX1)*SHP(1,IY1)
  156.       REL(IX,IY) = REL(IY,IX)
  157.   103 CONTINUE
  158.       IY1=0
  159.   102 CONTINUE
  160.       DO   104 IX=2+NBDL,LRE ,NBDL
  161.       IX2=IX - NBDL
  162.       DO   105 IY=1,IX2 ,NBDL
  163.       REL(IY+1,IX-1) = REL(IY,IX)
  164.       REL(IX-1,IY+1) = REL(IY+1,IX-1)
  165.   105 CONTINUE
  166.   104 CONTINUE
  167. C
  168. C     TERMES EN PI * (UR , RT )
  169. C
  170. C
  171.       IX1=0
  172.       IY1=0
  173.       DO   106 IX=3,LRE ,NBDL
  174.       IX1=IX1 + 1
  175.       DO   107 IY=2,IX  ,NBDL
  176.       IY1=IY1 + 1
  177.       REL(IY,IX) = REL(IY,IX) + VKL23*DJAC1*SHP(2,IY1)
  178.      #*SHP(6,IX1)
  179.       REL(IY,IX+1) = REL(IY,IX+1) + VKL23*DJAC1*SHP(2,IY1)
  180.      #*SHP(6,IX1+4)
  181.       REL(IX,IY) = REL(IY,IX)
  182.       REL(IX+1,IY) = REL(IY,IX+1)
  183.   107 CONTINUE
  184.       IY1=0
  185.   106 CONTINUE
  186.       IX1=1
  187.       IY1=0
  188.       DO   108 IX=2+NBDL,LRE ,NBDL
  189.       IX1=IX1 + 1
  190.       DO   109 IY=3,IX  ,NBDL
  191.       IY1=IY1 + 1
  192.       REL(IY,IX) = REL(IY,IX) + VKL23*DJAC1*SHP(2,IX1)
  193.      #*SHP(6,IY1)
  194.       REL(IY+1,IX) = REL(IY+1,IX) + VKL23*DJAC1*SHP(2,IX1)
  195.      #*SHP(6,IY1+4)
  196.       REL(IX,IY) = REL(IY,IX)
  197.       REL(IX,IY+1) = REL(IY+1,IX)
  198.   109 CONTINUE
  199.       IY1=0
  200.   108 CONTINUE
  201. C
  202.       IF ( IFOU.EQ.1) THEN
  203. C
  204. C     TERMES EN PI * (UT , RR )
  205. C     NON NULS QU EN FOURIER
  206. C
  207.       DJAC2 = ABS(DJAC)*POIGAU
  208. C
  209.       VKL25 = -1.D0* VKL23*NHARM / RR
  210.       IX1=0
  211.       IY1=0
  212.       DO   110 IX=5,LRE ,NBDL
  213.       IX1=IX1 + 1
  214.       DO   111 IY=2,IX  ,NBDL
  215.       IY1=IY1 + 1
  216.       REL(IY,IX) = REL(IY,IX) + VKL25*DJAC2*SHP(1,IY1)*SHP(6,IX1)
  217.       REL(IY,IX+1) = REL(IY,IX+1) + VKL25*DJAC2*SHP(1,IY1)
  218.      #*SHP(6,IX1+4)*(-1.D0)
  219.       REL(IX,IY) = REL(IY,IX)
  220.       REL(IX+1,IY) = REL(IY,IX+1)
  221.   111 CONTINUE
  222.       IY1=0
  223.   110 CONTINUE
  224.       IX1=1
  225.       IY1=0
  226.       DO   112 IX=2+NBDL,LRE ,NBDL
  227.       IX1=IX1 + 1
  228.       DO   113 IY=5,IX  ,NBDL
  229.       IY1=IY1 + 1
  230.       REL(IY,IX) = REL(IY,IX) + VKL25*DJAC2*SHP(1,IX1)*SHP(6,IY1)
  231.       REL(IY+1,IX) = REL(IY+1,IX) + VKL25*DJAC2*SHP(1,IX1)
  232.      #*SHP(6,IY1+4)*(-1.D0)
  233.       REL(IX,IY) = REL(IY,IX)
  234.       REL(IX,IY+1) = REL(IY+1,IX)
  235.   113 CONTINUE
  236.       IY1=0
  237.   112 CONTINUE
  238.       ENDIF
  239. C
  240. C     TERMES EN (UR,RT ) * (UR , RT )
  241. C
  242. C
  243.       IX1=0
  244.       IY1=0
  245.       DO   114 IX=3,LRE ,NBDL
  246.       IX1=IX1 + 1
  247.       DO   115 IY=3,IX  ,NBDL
  248.       IY1=IY1 + 1
  249. C
  250.       REL(IY,IX) = REL(IY,IX) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1)
  251.       REL(IY,IX+1) = REL(IY,IX+1) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1+4)
  252.       REL(IY+1,IX) = REL(IY+1,IX) + VKL33*DJAC1*SHP(6,IY1+4)*SHP(6,IX1)
  253.       REL(IY+1,IX+1) = REL(IY+1,IX+1)+VKL33*DJAC1*SHP(6,IY1+4)
  254.      #*SHP(6,IX1+4)
  255.       REL(IX,IY) = REL(IY,IX)
  256.       REL(IX+1,IY) = REL(IY,IX+1)
  257.       REL(IX,IY+1) = REL(IY+1,IX)
  258.       REL(IX+1,IY+1) = REL(IY+1,IX+1)
  259.   115 CONTINUE
  260.       IY1=0
  261.   114 CONTINUE
  262. C
  263.       IF ( IFOU.EQ.1) THEN
  264. C
  265. C     TERMES EN (UT,RR ) * (UT , RR )
  266. C     NON NULS QU EN FOURIER
  267. C
  268.       IX1=0
  269.       IY1=0
  270.       DO   303 IX=5,LRE ,NBDL
  271.       IX1=IX1 + 1
  272.       DO   403 IY=5,IX  ,NBDL
  273.       IY1=IY1 + 1
  274.       REL(IY,IX) = REL(IY,IX) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1)
  275.       REL(IY,IX+1) = REL(IY,IX+1) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1+4)
  276.      #*(-1.D0)
  277.       REL(IY+1,IX) = REL(IY+1,IX) + VKL33*DJAC1*SHP(6,IY1+4)*SHP(6,IX1)
  278.      #*(-1.D0)
  279.       REL(IY+1,IX+1) = REL(IY+1,IX+1)+VKL33*DJAC1*SHP(6,IY1+4)
  280.      #*SHP(6,IX1+4)
  281.       REL(IX,IY) = REL(IY,IX)
  282.       REL(IX+1,IY) = REL(IY,IX+1)
  283.       REL(IX,IY+1) = REL(IY+1,IX)
  284.       REL(IX+1,IY+1) = REL(IY+1,IX+1)
  285.   403 CONTINUE
  286.       IY1=0
  287.   303 CONTINUE
  288. C
  289.       ENDIF
  290.       IRET = 1
  291.       GOTO  666
  292. C
  293. C     MESSAGE D ERREUR : ELEMENT A  SURFACE  NULLE
  294. C
  295.  667  CONTINUE
  296.       IRET = 2
  297.       GOTO  666
  298. C
  299.  666  CONTINUE
  300.       RETURN
  301.       END
  302.  
  303.  
  304.  
  305.  
  306.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales