Télécharger la plaquette Cast3M


fr en

Télécharger hhosig.eso

Retour à la liste

Numérotation des lignes : 

hhosig

  1. C HHOSIG    SOURCE    OF166741  24/05/06    21:15:15     11082          
  2. C HHOSIG    SOURCE    FANDEUR                                           
  3.  
  4. *----------------------------------------------------------------------*
  5. *   CALCUL DES CONTRAINTES POUR LA FORMULATION 'HHO'                   *
  6. *----------------------------------------------------------------------*
  7. *   ENTREES :                                                          *
  8. *        MATE     Numero du materiau                                   *
  9. *        IVAMAT   Pointeur sur un segment MPTVAL pour le materiau      *
  10. *        NVMAT    Nombre de composantes obligatoires materiau          *
  11. *                                                                      *
  12. *   SORTIES :                                                          *
  13. *----------------------------------------------------------------------*
  14.  
  15.       SUBROUTINE HHOSIG(imoHHO, ichDEP,nmoDEP,
  16.      &                  IIPDPG,UZDPG,RYDPG,RXDPG,
  17.      &                  MATE,IVAMAT,NVMAT, IPMINT, NBPTEL,
  18.      &                  IVASTR,NCSTR, iret)
  19.  
  20.       IMPLICIT INTEGER(I-N)
  21.       IMPLICIT REAL*8(A-H,O-Z)
  22.  
  23. -INC PPARAM
  24. -INC CCOPTIO
  25. -INC CCREEL
  26.  
  27. -INC CCHHOPA
  28. -INC CCHHOPR
  29.  
  30. -INC SMCOORD
  31. -INC SMCHAML
  32. -INC SMELEME
  33. -INC SMINTE
  34. -INC SMLENTI
  35.       POINTEUR mlent4.mlenti
  36. -INC SMLMOTS
  37. -INC SMLREEL
  38.       POINTEUR mlrdef.mlreel, mlrdec.mlreel, mlrmbh.mlreel
  39. -INC SMMODEL
  40.  
  41.       SEGMENT MPTVAL
  42.         INTEGER IPOS(NS),NSOF(NS)
  43.         INTEGER IVAL(NCOSOU)
  44.         CHARACTER*16 TYVAL(NCOSOU)
  45.       ENDSEGMENT
  46.  
  47.       SEGMENT MWKMAT
  48.         REAL*8 VALMAT(NTMAT), HOELAS(lhook,lhook)
  49. c**        REAL*8 XLOC(3,3),XGLOB(3,3),TXR(IDIM,IDIM)
  50.       ENDSEGMENT
  51.  
  52.       SEGMENT MWKHHO
  53.         INTEGER TABINT(NBINT)
  54.         REAL*8  TABFLO(NBFLO)
  55.       ENDSEGMENT
  56.  
  57.       DIMENSION UDPGE(3)
  58.  
  59.       iret = 0
  60.  
  61.       imodel = imoHHO
  62. c*      segact,imodel <- actif en entree/sortie
  63.  
  64. C- Premieres verifications :
  65.       CALL HHONOB(imoHHO, nobHHO, iret)
  66.       IF (nobHHO.LE.0) THEN
  67.         write(ioimp,*) 'HHOSIG: IMODEL incorrect (not HHO)'
  68.         iret = 5
  69.         return
  70.       END IF
  71.  
  72. C- Introduction du point autour duquel se fait le mouvement
  73. C  de la section en defo plane generalisee
  74. C  IIPDPG = numero du noeud/point support si defini pour le modele
  75. C  NDPGE > 0 si prise en compte du point support
  76.       IF (IIPDPG.GT.0) THEN
  77.         write(ioimp,*) 'HHOSIG: 2D PLAN GENE mode not implemented!'
  78.         iret = 5
  79.         return
  80.         IF (IFOUR.EQ.-3) THEN
  81.           NDPGE = 3
  82.           UDPGE(1) = UZDPG
  83.           UDPGE(2) = RYDPG
  84.           UDPGE(3) = RXDPG
  85. C          SEGACT,MCOORD
  86.           IREF = (IIPDPG-1)*(IDIM+1)
  87.           XDPGE = mcoord.XCOOR(IREF+1)
  88.           YDPGE = mcoord.XCOOR(IREF+2)
  89.         ELSE IF (IFOUR.EQ.11) THEN
  90.           NDPGE = 2
  91.           UDPGE(1) = UZDPG
  92.           UDPGE(2) = RXDPG
  93.           UDPGE(3) = XZero
  94.           XDPGE = XZero
  95.           YDPGE = XZero
  96.         ELSE IF (IFOUR.EQ. 7 .OR. IFOUR.EQ. 8 .OR. IFOUR.EQ. 9 .OR.
  97.      &           IFOUR.EQ.10 .OR. IFOUR.EQ.14) THEN
  98.           NDPGE = 1
  99.           UDPGE(1) = UZDPG
  100.           UDPGE(2) = XZero
  101.           UDPGE(3) = XZero
  102.           XDPGE = XZero
  103.           YDPGE = XZero
  104.         else
  105.           write(ioimp,*) 'HHOSIG: ERREUR NDPGE'
  106.           iret = 5
  107.           return
  108.         END IF
  109.       ELSE
  110.         NDPGE = 0
  111.         UDPGE(1) = XZero
  112.         UDPGE(2) = XZero
  113.         UDPGE(3) = XZero
  114.         XDPGE = XZero
  115.         YDPGE = XZero
  116.       END IF
  117.  
  118. C- Recuperation des donnees du modele en entree
  119. c*      MELE   = imodel.NEFMOD
  120. c*      MFR    = imodel.infele(13)
  121.       meleme = imodel.IMAMOD
  122.       NBNOE  = meleme.NUM(/1)
  123.       NBELT  = meleme.NUM(/2)
  124.  
  125.       mlenti = imodel.IVAMOD(nobHHO+1)
  126. c*      segact,mlenti
  127.       mlent3 = imodel.IVAMOD(nobHHO+3)
  128. c*      segact,mlent3
  129.       mlent4 = imodel.IVAMOD(nobHHO+4)
  130. c*      segact,mlent4
  131.  
  132. c*      n_o_face = mlenti.lect(2)
  133.       n_d_face = mlenti.lect(3)
  134. c*      n_o_cell = mlenti.lect(4)
  135.       n_d_cell = mlenti.lect(5)
  136.       nb_faces = mlenti.lect(7)
  137.       NBPGAU   = mlenti.lect(8)
  138.       idifo    = mlenti.lect(9)
  139.       NBDDL    = mlenti.lect(11)
  140.       NDEPF    = idifo * n_d_face
  141.       NDEPC    = idifo * n_d_cell
  142.       lhook = 9
  143.  
  144.       nbel3 = mlent3.lect(/1) / 2
  145.       nbfa3 = 2 * nb_faces
  146.       nbel4 = mlent4.lect(/1) / 2
  147.  
  148.       IF (mlenti.lect(6).NE.NBNOE) THEN
  149.         write(ioimp,*) 'HHOSIG: Bizarre nb_vertices'
  150.       END IF
  151.       IF (NBDDL .NE. imodel.INFELE(9)) then
  152.         write(ioimp,*) 'HHOSIG: Bizarre NBDDL'
  153.       END IF
  154. c NBPGAU =? (NBPTEL = imodel.INFELE(4))
  155.       IF (NBPGAU .NE. imodel.INFELE(4)) then
  156.         write(ioimp,*) 'HHOSRIG: Bizarre nb.p.gau(1)'
  157.       END IF
  158. c NBPGAU =? minte.POIGAU(/1)
  159.       minte = IPMINT
  160. c*      SEGACT minte <- actif en E/S
  161.       if (NBPGAU .NE. minte.POIGAU(/1)) then
  162.         write(ioimp,*) 'HHOSIG: Bizarre nb.p.gau (2)'
  163.       end if
  164. c-dbg      write(ioimp,*) 'HHOSIG:',nb_faces,NBDDL,NBPGAU
  165.       if (nbel3.NE.(NBELT*nb_faces)) then
  166.         write(ioimp,*) 'HHOSIG: Bizarre nbel3'
  167.       end if
  168.       if (nbel4.NE.NBELT) then
  169.         write(ioimp,*) 'HHOSIG: Bizarre nbel4'
  170.       end if
  171.  
  172.       ivid = 1
  173. C- Deplacements des faces et des cellules :
  174.       nomid = nmoDEP
  175. c*      segact,nomid
  176.       JGN = nomid.lesobl(/1)
  177. c-dbg      write(ioimp,*) 'HHOSIG=',NDEPC,NDEPF,lesobl(/2),lesfac(/2)
  178.       nfac = 0
  179. C Deplacements des cellules - Points supports des cellules
  180.       JGM = NDEPC + nfac
  181.       SEGINI,mlmots
  182.       DO i = 1, NDEPC
  183.         mlmots.MOTS(i)(1:JGN) = nomid.lesobl(i)(1:JGN)
  184.       END DO
  185. c*      DO i = 1, nfac
  186. c*        mlmots.MOTS(NDEP+i)(1:JGN) = nomid.lesfac(i)(1:JGN)
  187. c*      END DO
  188.       CALL EXTR23(ichDEP,mlmots,MPCHHO,mlrDEC,ivid)
  189.       SEGACT,mlrDEC
  190. c*      IF (IERR.NE.0) THEN
  191. c*        iret = 21
  192. c*        return
  193. c*      END IF
  194.  
  195. c-dbg      write(ioimp,*) 'mlmots DEPC ',(mots(i),i=1,mots(/2))
  196. c-dbg      write(ioimp,*) 'U.CELL',mlrdec.prog(/1),NCEHHO,NDEPC
  197. c-dbg      write(ioimp,*) '      ',(mlrdec.prog(i),i=1,mlrdec.prog(/1))
  198. c-dbg      write(ioimp,*)
  199.  
  200. C  Deplacements des faces - Points supports des faces
  201.       JGM = NDEPF + nfac
  202.       SEGADJ,mlmots
  203.       DO i = 1, NDEPF
  204.         mlmots.MOTS(i)(1:JGN) = nomid.lesobl(NDEPC+i)(1:JGN)
  205.       END DO
  206. c*      DO i = 1, nfac
  207. c*        mlmots.MOTS(NDEP+i)(1:JGN) = nomid.lesfac(i)(1:JGN)
  208. c*      END DO
  209.       CALL EXTR23(ichDEP,mlmots,MPFHHO,mlrDEF,ivid)
  210.       SEGACT,mlrDEF
  211. c*      IF (IERR.NE.0) THEN
  212. c*        iret = 21
  213. c*        return
  214. c*      END IF
  215.  
  216. c-dbg      write(ioimp,*) 'MLMOTS DEPF ',(mots(i),i=1,mots(/2))
  217. c-dbg      write(ioimp,*) 'U.FACE',mlrDEF.prog(/1),NFAHHO,NDEPF
  218. c-dbg      write(ioimp,*) '      ',(mlrDEF.prog(i),i=1,mlrDEF.prog(/1))
  219. c-dbg      write(ioimp,*)
  220.  
  221. C- Verification des proprietes materielles : Champ CONSTANT par ELEMENT !
  222.       mptval = IVAMAT
  223.       NTMAT = mptval.IVAL(/1)
  224.       NUMAT = NTMAT - 4
  225.  
  226. C- Tableau de stockage des proprietes materielles
  227.       MWKMAT = 0
  228.       SEGINI,MWKMAT
  229.       IWKMAT = MWKMAT
  230.  
  231.       isz = 0
  232.       iunif = 0
  233.       DO is = 1, NTMAT
  234.         IF (is.LE.(NVMAT-4) .OR. is.GT.NVMAT) THEN
  235.           melval = mptval.IVAL(is)
  236.           IGMN = melval.VELCHE(/1)
  237.           IEMN = melval.VELCHE(/2)
  238.           IF (IGMN.NE.1) isz = isz + 1
  239.           IF (IGMN.EQ.1 .AND. IEMN.EQ.1) iunif = iunif + 1
  240.           VALMAT(is) = melval.VELCHE(1,1)
  241.         END IF
  242.       END DO
  243.       IF (isz.NE.0) THEN
  244.         write(ioimp,*) 'HHOSIG: Material properties must be constant ',
  245.      &                 'by element (for the moment) !'
  246.         iret = 21
  247.         RETURN
  248.       END IF
  249. * Le champ de proprietes est completement uniforme :
  250.       IF (iunif.NE.NUMAT) iunif = 0
  251.  
  252. C- Materiau isotrope et proprietes uniformes : -> calcul de HOELAS
  253.       IF (MATE.EQ.1 .AND. iunif.NE.0) THEN
  254.         CALL DOHMAS(VALMAT,'ISOTROPE',IFOUR,lhook,1,HOELAS,iretc)
  255.         r_z = HOELAS(4,4) * 2.D0
  256.         DO i = 4, lhook
  257.           HOELAS(i,i) = r_z
  258.         END DO
  259.       END IF
  260. C - Materiaux orthotrope et anisotrope
  261.       IF (MATE.EQ.2.OR.MATE.EQ.3) THEN
  262. C** AFAIRE !!!        
  263.       END IF
  264.  
  265. C- Les donnees liees a la formulation HHO :
  266.       IVMBHO = mptval.IVAL(NVMAT-1)
  267.       melval = IVMBHO
  268.       IGMB = melval.IELCHE(/1)
  269.       IEMB = melval.IELCHE(/2)
  270. c-dbg      write(ioimp,*) 'IVMBHO',melval,igmb,iemb,tyval(nvmat-1)
  271.       mlrmbh = melval.IELCHE(1,1)
  272. c*      segact,mlrmbh
  273. c*      write(ioimp,*) 'HHOSIG MBHHO SIZE:',mlrmbh.prog(/1),
  274. c*     &                NBDDL,9*NBDDL,mlenti.lect(14)
  275.       IF ((mlrmbh.prog(/1).NE.(9*NBDDL)) .OR. 
  276.      &    (mlrmbh.prog(/1).NE.mlenti.lect(14))) THEN
  277.         write(ioimp,*) 'HHOSIG: BHHO matrix size incorrect'
  278.         iret = 21
  279.         RETURN
  280.       END IF
  281.  
  282. C- Indices et tableau de travail
  283.       ir_coo =          0
  284. c* si besoin des coordonnees      ir_eps = ir_coo + (IDIM*NBNOE)
  285.       ir_eps = ir_coo + 0
  286.       ir_sig = ir_eps + lhook
  287.       ir_uce = ir_sig + lhook
  288.       ir_ufa = ir_uce + NDEPC
  289.       ir_uge = ir_ufa + (NDEPF*nb_faces)
  290.       ir_fin = ir_uge + NDPGE
  291.  
  292.       NBINT = 1
  293.       NBFLO = ir_fin
  294.       SEGINI,MWKHHO
  295.  
  296.       IF (NDPGE.GT.0) THEN
  297.         DO ic = 1, NDPGE
  298.           TABFLO(ir_uge+ic) = UDPGE(ic)
  299.         END DO
  300.       END IF
  301. c* si besoin des coordonnees        SEGACT,MCOORD*NOMOD
  302.  
  303. C-------------------------
  304. C Boucle sur les cellules du modele
  305. C-------------------------
  306.       DO IEL = 1,NBELT
  307.  
  308. C- Recuperation des coordonnees des noeuds de l element IEL
  309. c*        CALL HHOCOO(meleme.num,NBNOE, mcoord.xcoor, IEL,
  310. c*    &               TABFLO(ir_coo+1), iret)
  311. c*        IF (iret.NE.0) RETURN
  312.  
  313. C- Valeurs des inconnues primales pour l'element IEL (cell+faces)
  314.         in1 = IEL * 2
  315.         je = mlent4.lect(in1-1)
  316.         ip = mlent4.lect(in1)
  317.       if (ip.le.0) write(ioimp,*) 'HHOSIG ICEL Bizarre...',iel,je,ip
  318.         jp = ip + NBCHHO(je-1)
  319.         DO ic = 1, NDEPC
  320.           jc = NCEHHO * (ic - 1)
  321.           TABFLO(ir_uce + ic) = mlrDEC.prog(jp + jc)
  322.         END DO
  323. c-dbg      write(ioimp,*) 'HHOSIG ICEL',iel,je,ip,jp,ir_uce
  324. c-dbg      write(ioimp,*) (TABFLO(ir_uce+ic),ic=1,ndepc)
  325.  
  326.         in1 = (IEL-1) * nbfa3
  327.         ir_kc = ir_ufa
  328.         DO j1 = 1, nb_faces
  329.           in2 = in1 + (2 * j1)
  330.           je = mlent3.lect(in2-1)
  331.           ip = ABS(mlent3.lect(in2))
  332.       if (ip.eq.0) write(ioimp,*) 'HHOSIG IFAE Bizarre...',iel,j1,je,ip
  333.           jp = ip + NBFHHO(je-1)
  334.           DO ic = 1, NDEPF
  335.             jc = NFAHHO * (ic - 1)
  336.             TABFLO(ir_kc + ic) = mlrDEF.prog(jp + jc)
  337.           END DO
  338. c-dbg      write(ioimp,*) 'HHOSIG IFAE',iel,j1,je,ip,jp,ir_kc
  339. c-dbg      write(ioimp,*) (TABFLO(ir_kc+ic),ic=1,ndepf)
  340.           ir_kc = ir_kc + NDEPF
  341.         END DO
  342.  
  343. c-dbg      write(ioimp,*) 'HHOSIG ...',ir_uce,ir_kc,NBDDL
  344.  
  345. C- Recuperation des proprietes constantes par element :
  346.         IF (iunif.EQ.0) THEN
  347.           DO is = 1, NVMAT
  348.             IF (is.LE.(NVMAT-4) .OR. is.GT.NVMAT) THEN
  349.               melval = mptval.IVAL(is)
  350.               VALMAT(is) = melval.VELCHE(1,IEL)
  351.             END IF
  352.           END DO
  353.         END IF
  354. C- Materiau isotrope
  355.         IF (MATE.EQ.1 .AND. iunif.EQ.0) THEN
  356.           CALL DOHMAS(VALMAT,'ISOTROPE',IFOUR,lhook,1,HOELAS,iretc)
  357.           r_z = HOELAS(4,4) * 2.D0
  358.           DO i = 4, lhook
  359.             HOELAS(i,i) = r_z
  360.           END DO
  361.         ENDIF
  362. C - Calcul des axes locaux dans les cas orthotrope et anisotrope
  363.         IF (MATE.EQ.2.OR.MATE.EQ.3) THEN
  364.         ENDIF
  365.  
  366.         JEMB = MIN(IEL,IEMB)
  367. C--  --  --  --  --  --  --  --  --
  368. C - Boucle sur les points de Gauss
  369. C--  --  --  --  --  --  --  --  --
  370.         DO IGAU = 1, NBPGAU
  371.  
  372. C-AFAIRE C -- Recuperation des proprietes materielles (IGAU)
  373. C-AFAIRE           DO is = 1, NVMAT
  374. C-AFAIRE             IF (is.LE.(NVMAT-4) .OR. is.GT.NVMAT) THEN
  375. C-AFAIRE               melval = mptval.IVAL(is)
  376. C-AFAIRE               JEMN = MIN(IEL ,melval.VELCHE(/2))
  377. C-AFAIRE               JGMN = MIN(IGAU,melval.VELCHE(/1))
  378. C-AFAIRE               VALMAT(is) = melval.VELCHE(JGMN,JEMN)
  379. C-AFAIRE             END IF
  380. C-AFAIRE           END DO
  381. C-AFAIREC- Materiau isotrope
  382. C-AFAIRE           IF (MATE.EQ.1) THEN
  383. C-AFAIRE             CALL DOHMAS(VALMAT,'ISOTROPE',IFOUR,lhook,1,HOELAS,iretc)
  384. C-AFAIRE             r_z = HOELAS(4,4) * 2.D0
  385. C-AFAIRE             DO i = 4, lhook
  386. C-AFAIRE               HOELAS(i,i) = r_z
  387. C-AFAIRE             END DO
  388. C-AFAIRE           END IF
  389.  
  390.           melval = IVMBHO
  391.           JGMB = MIN(IGAU,IGMB)
  392.           mlrmbh = melval.IELCHE(JGMB,JEMB)
  393. c*          segact,mlrmbh
  394. c* !! matrice BHHO stockee colonne par colonne : lhook*NBDDL
  395.           DO ic = 1, lhook
  396.             r_z = XZero
  397.             DO jc = 1, NBDDL
  398.               jnc = lhook * (jc-1)
  399.               r_z = r_z + mlrmbh.prog(ic + jnc) * TABFLO(ir_uce + jc)
  400.             END DO
  401.             TABFLO(ir_eps + ic) = r_z
  402.           END DO
  403.           DO ic = 1, lhook
  404.             r_z = XZero
  405.             DO jc = 1, lhook
  406.               r_z = r_z + HOELAS(ic,jc) * TABFLO(ir_eps + jc)
  407.             END DO
  408.             TABFLO(ir_sig + ic) = r_z
  409.           END DO
  410.  
  411. C -- Remplissage du segment contenant les "contraintes"
  412. c* A voir selon ordre - permutation :
  413.           mptval = IVASTR
  414.           DO ic = 1, NCSTR
  415.             melval = mptval.IVAL(ic)
  416.             melval.velche(IGAU,IEL) = TABFLO(ir_sig + ic)
  417.           END DO
  418.  
  419. C--  --  --  --  --  --  --  --  --
  420.         END DO
  421. C--  --  --  --  --  --  --  --  --
  422. C-------------------------
  423.       END DO
  424. C-------------------------
  425.  
  426.       IF (MATE.EQ.2.OR.MATE.EQ.3) THEN
  427.       ENDIF
  428.       SEGSUP,MWKMAT,MWKHHO
  429.       SEGSUP,mlrDEC,mlrDEF
  430.       SEGSUP,mlmots
  431.  
  432. c*      RETURN
  433.       END
  434.  
  435.  
  436.  

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