Télécharger la plaquette Cast3M


fr en

Télécharger mhybr4.eso

Retour à la liste

Numérotation des lignes : 

mhybr4

  1. C MHYBR4    SOURCE    CB215821  24/04/12    21:16:43     11897          
  2.       SUBROUTINE MHYBR4(IMAIL,IPMODE,IPRIGI,IPGEOS)
  3. C-----------------------------------------------------------------------
  4. C Calcul de la matrice masse hybride utilisé pour prendre en compte les forces
  5. C de volume.
  6. C Traitement du cas des elements finis massifs a integration numerique
  7. C pour un maillage elementaire et une formulation hybride.
  8. C-----------------------------------------------------------------------
  9. C
  10. C---------------------------
  11. C Parametres Entree/Sortie :
  12. C---------------------------
  13. C
  14. C E/  IMAIL   : Numero du maillage elementaire considere,
  15. C               dans l'objet modele.
  16. C E/  IPMODE  : Pointeur sur un segment IMODEL.
  17. C E/  IPGEOS  : Pointeur sur le maillage sommet
  18. C E/S IPRIGI  : Pointeur sur l'objet rigidite resultat.
  19. C
  20. C----------------------
  21. C Variables en COMMON :
  22. C----------------------
  23. C
  24. C E/  XCOOR   : VOIR SMCOORD
  25. C E/  IERR    : VOIR CCOPTIO
  26. C E/  IDIM    : VOIR CCOPTIO
  27. C E/  INTERR  : VOIR CCOPTIO
  28. C E/  IFOMOD  : VOIR CCOPTIO
  29. C E/  XPETIT  : VOIR CCREEL
  30. C
  31. C----------------------
  32. C Tableaux de travail :
  33. C----------------------
  34. C
  35. C     NBNN    : Nombre de noeuds dans l'element considere
  36. C     NEFHYB  : Numero de l'element fini dans NOMTP.
  37. C     NEF     : Numero de l'element fini support geometrique
  38. C               dans NOMTP (voir CCHAMP)
  39. C     NBELEM  : Nombre d'element dans le maillage elementaire
  40. C     NBPGAU  : Nombre de points de gauss pour l'element fini NEF
  41. C**    CEL     : Matrice de conductivite elementaire
  42. C    CEL     : Matrice des masses elementaires
  43. C     XE      : Coordonnees des noeuds dans le repere global
  44. C**    CMAT    : Matrice de permeabilite dans le repere global
  45. C     SHP     : Tableau de travail contenant les fonctions de forme au
  46. C               point de gauss utilise + derivees
  47. C     SHY     : Contient les fonctions de base hybride en un point
  48. C               mais pas les derivees de la fonction de base.
  49. C**     VALMAT  : Valeurs des coeff. de la matrice CMAT et des
  50. C**               cosinus directeurs des axes d'ortho. / repere local
  51. C**     XGLOB   : Cosinus directeurs des axes d'ortho. / repere global
  52. C**     XLOC    : Cosinus directeurs des axes d'ortho. / repere local
  53. C**     TXR     : Cosinus directeurs des axes locaux / repere global
  54. C
  55. C
  56. C-----------------------------------------------------------------------
  57. C
  58. C Langage : ESOPE + FORTRAN77
  59. C
  60. C           02/96 L.V.BENET : fonction propre à l'option 'MASSE'
  61. C
  62. C-----------------------------------------------------------------------
  63.       IMPLICIT INTEGER(I-N)
  64.       IMPLICIT REAL*8(A-H,O-Z)
  65. *
  66. -INC CCHAMP
  67.  
  68. -INC PPARAM
  69. -INC CCOPTIO
  70. -INC CCREEL
  71. -INC SMCOORD
  72. -INC SMINTE
  73. -INC SMMODEL
  74. -INC SMRIGID
  75. -INC SMELEME
  76. ****-INC SMCHAML
  77. *
  78. *
  79.       SEGMENT MMAT1
  80.          REAL*8  CEL(NBDDL,NBDDL),CEL1(NBDDL,NBDDL),XE(3,NBNN)
  81.          REAL*8  SHP(6,NBNN),SHY(IDIM,NBDDL),ZJAC(IDIM,IDIM)
  82.          REAL*8  CMAT(IDIM,IDIM),CMAT1(IDIM,IDIM),CMAT2(IDIM,IDIM)
  83.          INTEGER ICSTO(NBDDL)
  84.       ENDSEGMENT
  85. *
  86.       SEGMENT NOTYPE
  87.          CHARACTER*16 TYPE(NBTYPE)
  88.       ENDSEGMENT
  89. *
  90.       SEGMENT HYBSTO
  91.          REAL*8 HYBASE(NDIM,NBDDL,NBPP)
  92.       ENDSEGMENT
  93. *
  94.       CHARACTER*8  CNM
  95.       PARAMETER(NINF=3)
  96.       INTEGER INFOS(NINF)
  97. *
  98. * Recup. des caracteristiques geometriques du maillage elementaire
  99. * et du maillage hybride dual
  100. *
  101.       IMODEL = IPMODE
  102.       SEGACT IMODEL
  103.       IPMAIL = IMAMOD
  104.       MELEME = IPGEOS
  105.       SEGACT MELEME
  106.       NBNN   = NUM(/1)
  107.       NBELEM = NUM(/2)
  108.       NEFHYB = NEFMOD
  109.       NEF    = NUMGEO(NEFHYB)
  110.       MFR    = NUMMFR(NEFHYB)
  111. *
  112.       MRIGID = IPRIGI
  113.       SEGACT MRIGID
  114.       IPT1   = IRIGEL(1,IMAIL)
  115.       SEGACT IPT1
  116.       NBDDL  = IPT1.NUM(/1)
  117.       SEGDES IPT1
  118. *
  119. * Recup. des caracteristiques d'integration de l'EF support geometrique
  120. * de l'EF hybride
  121. *
  122.       CALL TSHAPE(NEF,'GAUSS',IPINTE)
  123.       IF (IERR.NE.0) THEN
  124.          SEGDES IMODEL , MELEME
  125.          RETURN
  126.       ENDIF
  127. *
  128. * Recup. des fonctions de bases hybrides
  129. *
  130.       CALL HYSHPT(NEFHYB,NBDDL,IPINTE,IPTHYB)
  131.       IF (IERR.NE.0) THEN
  132.          SEGDES IMODEL , MELEME
  133.          RETURN
  134.       ENDIF
  135. *
  136. * Activation des segments "d'integration"
  137. *
  138.       MINTE  = IPINTE
  139.       SEGACT MINTE
  140.       NBPGAU = POIGAU(/1)
  141.       HYBSTO = IPTHYB
  142.       SEGACT HYBSTO
  143. *
  144. * Recup. des caracteristiques d'integration au centre de l'EF
  145. *
  146.       CALL RESHPT(1,NBNN,NEF,NEF,0,IPT1,IRT1)
  147.       MINTE1 = IPT1
  148.       SEGACT MINTE1
  149. *
  150. * Initialisation des chapeaux de l'objet rigidité
  151. *
  152.       xMATRI = IRIGEL(4,IMAIL)
  153.       SEGACT xMATRI*MOD
  154.       NLIGRP = NBDDL
  155.       NLIGRD = NBDDL
  156. *
  157. * Remplissage du tableau INFOS (informations sur element).
  158. *
  159.       INFOS(1) = 0
  160.       INFOS(2) = 0
  161.       INFOS(3) = NIFOUR
  162. *
  163. * Initialisation des tableaux de travail
  164. *
  165.       NDIM = IDIM * (IDIM+1)
  166.       SEGINI MMAT1
  167.       CALL ZERO(CMAT2,IDIM,IDIM)
  168.       DO 5 I=1,IDIM
  169.            CMAT2(I,I)=1.D0
  170.   5   CONTINUE
  171. *
  172. *-------------------------------------------------------
  173. * BOUCLE SUR LES ELEMENTS DU MAILLAGE ELEMENTAIRE IMAIL
  174. *-------------------------------------------------------
  175. *
  176.       DO 30 IEL=1,NBELEM
  177. *
  178. * Initialisations
  179. *
  180.          IFOIS = 0
  181.          IFOI2 = 0
  182.          CALL ZERO(CEL,NBDDL,NBDDL)
  183. *
  184. * Recuperation des coordonnees globales des noeuds de l'element IEL
  185. *
  186.          CALL DOXE(XCOOR,IDIM,NBNN,NUM,IEL,XE)
  187. *--------------------------------
  188. * BOUCLE SUR LES POINTS DE GAUSS
  189. *--------------------------------
  190.          DO 20 IGAU=1,NBPGAU
  191. *
  192. * Calcul de la matrice jacobienne de la fonction de passage du
  193. * repere local au repere global, de son determinant ET recup.
  194. * des fonctions de base hybride au point de gauss IGAU.
  195. *
  196.             CALL MHYBR3(IGAU,NBNN,NBDDL,NDIM,IDIM,IDIM,XE,HYBASE,
  197.      S                  SHPTOT,SHY,SHP,ZJAC,DJAC)
  198. *
  199. * Controle du maillage
  200. *
  201.             IF (DJAC.LT.0.D0) IFOIS = IFOIS + 1
  202.             IF (ABS(DJAC).LT.XPETIT) THEN
  203.                IFOI2 = IFOI2 + 1
  204.                DJAC  = XPETIT
  205.             ENDIF
  206. *
  207. * Calcul du poids d'integration global affecte dans DJAC.
  208. *
  209.             DJAC = POIGAU(IGAU) / ABS(DJAC)
  210. *
  211. *- Calcul de la contribution du point de gauss IGAU a la matrice
  212. *- elementaire CEL de l'element IEL :
  213. *-    POIGAU/DJAC * transpose( ZJAC*SHY ) *inv(CMAT)* ( ZJAC*SHY )
  214. *- On calcule CMAT2=inv(CMAT)                    avec INVRS puis
  215. *- on calcule CMAT1=transpose(ZJAC)*CMAT2*ZJAC   avec PRODT puis
  216. *- on somme POIGAU/DJAC * transp.(SHY)*CMAT1*SHY avec BDBST.
  217. *
  218.             CALL PRODT(CMAT1,CMAT2,ZJAC,IDIM,IDIM)
  219.             CALL BDBST(SHY,DJAC,CMAT1,NBDDL,IDIM,CEL)
  220.   20     CONTINUE
  221. *
  222. * Le jacobien est negatif --> MAILLAGE INCORRECT
  223. *
  224.          IF (IFOIS.NE.0.AND.IFOIS.NE.NBPGAU) THEN
  225.             INTERR(1) = IEL
  226.             CALL ERREUR(195)
  227.             SEGSUP xMATRI , MRIGID
  228.             GOTO 40
  229.          ENDIF
  230. *
  231. * Le jacobien est tres petit --> MAILLAGE INCORRECT
  232. *
  233.          IF (IFOI2.EQ.NBPGAU) THEN
  234.             INTERR(1) = IEL
  235.             CALL ERREUR(259)
  236.             SEGSUP xMATRI , MRIGID
  237.             GOTO 40
  238.          ENDIF
  239. *
  240. * Remplissage de XMATRI
  241. *
  242. *         SEGINI XMATRI
  243. *         IMATTT(IEL) = XMATRI
  244.          CALL REMPMT(CEL,NBDDL,RE(1,1,iel))
  245. *         SEGDES XMATRI
  246.   30  CONTINUE
  247. *
  248. * Desactivation des segments
  249. *
  250.       SEGDES xMATRI , MRIGID
  251.   40  CONTINUE
  252.       SEGSUP MMAT1 , HYBSTO
  253.       SEGDES MELEME
  254.       SEGDES IMODEL
  255.       SEGDES MINTE , MINTE1
  256.       RETURN
  257.       END
  258.  
  259.  
  260.  
  261.  
  262.  
  263.  

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