Télécharger la plaquette Cast3M


fr en

Télécharger epsi2.eso

Retour à la liste

Numérotation des lignes : 

epsi2

  1. C EPSI2     SOURCE    CB215821  24/04/12    21:15:46     11897          
  2.  
  3.       SUBROUTINE EPSI2(IPMAIL,IPMINT,MELE,IELE,
  4.      &  IVADEP,NBPTEL,LRE,NSTRS,LHOOK,
  5.      &  MFR,NDEP,IPORE,IREPS2,NBPGAU,IVAEPS,UZDPG,RYDPG,RXDPG,IIPDPG,
  6.      &  IDERI,ivamat,ivade2,mate,nmatt,cmate,ngra,noer,kerr)
  7.  
  8. C---------------------------------------------------------------------*
  9. C
  10. C            calcul des deformations
  11. C
  12. C             massif, poreux, joints poreux, incompressibles
  13. C---------------------------------------------------------------------*
  14. C                                                                     *
  15. C   entrees :                                                         *
  16. C   ________                                                          *
  17. C                                                                     *
  18. C        ipmail   pointeur sur un segment  meleme                     *
  19. C        ipmint   pointeur sur un segment minte                       *
  20. C        mele     numero de l'element fini                            *
  21. C        iele     numero geometrique de l'element                     *
  22. C        nbpgau   nombre de point d'integration pour la rigidite      *
  23. C        ivadep   pointeur sur le chamelem de deplacements            *
  24. C        nbptel   nombre de points par element                        *
  25. C        lre      nombre de ddl dans la matrice de rigidite           *
  26. C        nstrs    nombre de composante de contraintes/deformations    *
  27. C                 pour une matrice de hooke                           *
  28. C        lhook    dimension de la matrice de hooke                    *
  29. C        mfr      numero de la formulation de l'element fini          *
  30. C        ndep     nombre de composantes de deplacements               *
  31. C        ipore    nombre de fonctions de forme                        *
  32. C        iresp2   flag pour indiquer si on veut les contraintes       *
  33. C                  de piola-kirchhoff                                 *
  34. C        uzdpg    = deformation au point nsdpge support de la         *
  35. C        rydpf    = deformation plane generalisee                     *
  36. C        rxdpg    =                                                   *
  37. C                                                                     *
  38. C   sorties :                                                         *
  39. C   ________                                                          *
  40. C                                                                     *
  41. C        ivaeps   pointeur sur un segment mptval contenant les        *
  42. C                 les melvals de deformations                         *
  43. C---------------------------------------------------------------------*
  44. C Pour MEMOIRE : Si MELE element incompressible alors MFR = 31
  45. C---------------------------------------------------------------------*
  46.  
  47.       IMPLICIT INTEGER(I-N)
  48.       IMPLICIT REAL*8(A-H,O-Z)
  49.  
  50. -INC PPARAM
  51. -INC CCOPTIO
  52. -INC CCREEL
  53. -INC CCHAMP
  54. -INC CCGEOME
  55.  
  56. -INC SMCHAML
  57. -INC SMCHPOI
  58. -INC SMELEME
  59. -INC SMCOORD
  60. -INC SMMODEL
  61. -INC SMINTE
  62. -INC SMLREEL
  63.  
  64.       SEGMENT WRK1
  65.        REAL*8 DDHOOK(NSTRS,NSTRS),XDDL(LRE),XSTRS(NSTRS)
  66.        REAL*8 XE(3,NBBB),DDHOMU(NSTRS,NSTRS)
  67.        REAL*8 SHPWRK(6,NBNO),BGENE(LHOOK,LRE)
  68.        REAL*8 XE1(3,NBBB),XE2(3,NBBB),xstrs2(NSTRS)
  69.        REAL*8 xjac(3,3),valmat(20)
  70.       ENDSEGMENT
  71.       SEGMENT WRK2
  72.         REAL*8 BGR(NGRA,LRE),BB(2,NGRA),gradi(ngra),R(ngra),u(ngra)
  73.         REAL*8 TENS(9),tentra(9),xddls2(lre)
  74.       ENDSEGMENT
  75. C
  76.       SEGMENT WRK3
  77.        REAL*8 BPSS(3,3),XEL(3,NBBB)
  78.        REAL*8 XNTH(LPP,LPP),XNTB(LPP,LPP),XNTT(LPP)
  79.       ENDSEGMENT
  80. C
  81.       SEGMENT WRK5
  82.        REAL*8 XGENE(NSTN,LRN)
  83.       ENDSEGMENT
  84. C
  85.       SEGMENT NOTYPE
  86.         CHARACTER*16 TYPE(NBTYPE)
  87.       ENDSEGMENT
  88. C
  89.       SEGMENT MPTVAL
  90.        INTEGER IPOS(NS),NSOF(NS)
  91.        INTEGER IVAL(NCOSOU)
  92.        CHARACTER*16 TYVAL(NCOSOU)
  93.       ENDSEGMENT
  94. C
  95.       SEGMENT MTRACE
  96.         REAL*8 TRACE(NBPTEL)
  97.       ENDSEGMENT
  98.  
  99.       CHARACTER*(NCONCH) CONM
  100.       PARAMETER (NINF=3)
  101.       INTEGER INFOS(NINF)
  102.       CHARACTER*8 CMATE
  103.  
  104.       DIMENSION A(4,60),BBX(3,60),UDPGE(3)
  105.       DIMENSION IN(6),JN(6),ITAB(3,3),PP(4,4)
  106.       real*8 valcar(12),var(3)
  107.       real*8 cobma(lhook)
  108. C
  109.       DATA IN/1,2,3,1,1,2/
  110.       DATA JN/1,2,3,2,3,3/
  111. C
  112.       DATA ITAB(1,1),ITAB(1,2),ITAB(1,3)/1,4,5/
  113.       DATA ITAB(2,1),ITAB(2,2),ITAB(2,3)/4,2,6/
  114.       DATA ITAB(3,1),ITAB(3,2),ITAB(3,3)/5,6,3/
  115. C
  116.       real*8 s(2)
  117.  
  118.       s(1)=0.d0
  119.       s(2)=0.d0
  120.       kerr=0
  121. C  Introduction du point autour duquel se fait le mouvement
  122. C  de la section en defo plane generalisee
  123. C  IIPDPG = numero du noeud/point support si defini pour le modele
  124. C  NDPGE > 0 si prise en compte du point support
  125.       IF (IIPDPG.GT.0) THEN
  126.         IF (IFOUR.EQ.-3) THEN
  127.           NDPGE=3
  128.           UDPGE(1)=UZDPG
  129.           UDPGE(2)=RYDPG
  130.           UDPGE(3)=RXDPG
  131. C          SEGACT,MCOORD
  132.           IREF=(IIPDPG-1)*(IDIM+1)
  133.           XDPGE=XCOOR(IREF+1)
  134.           YDPGE=XCOOR(IREF+2)
  135.         ELSE IF (IFOUR.EQ.11) THEN
  136.           NDPGE=2
  137.           UDPGE(1)=UZDPG
  138.           UDPGE(2)=RXDPG
  139.           UDPGE(3)=XZero
  140.           XDPGE=XZero
  141.           YDPGE=XZero
  142.         ELSE IF (IFOUR.EQ. 7 .OR. IFOUR.EQ. 8 .OR. IFOUR.EQ. 9 .OR.
  143.      &           IFOUR.EQ.10 .OR. IFOUR.EQ.14) THEN
  144.           NDPGE=1
  145.           UDPGE(1)=UZDPG
  146.           UDPGE(2)=XZero
  147.           UDPGE(3)=XZero
  148.           XDPGE=XZero
  149.           YDPGE=XZero
  150.         else
  151.           write(ioimp,*) 'EPSI2 : ERREUR NDPGE'
  152.           call erreur(5)
  153.           return
  154.         ENDIF
  155.       ELSE
  156.         NDPGE=0
  157.         UDPGE(1)=UZDPG
  158.         UDPGE(2)=XZero
  159.         UDPGE(3)=XZero
  160.         XDPGE=XZero
  161.         YDPGE=XZero
  162.       ENDIF
  163. C
  164.  
  165.       MELEME=IPMAIL
  166.       NBNN=NUM(/1)
  167.       NBELEM=NUM(/2)
  168. C
  169.       NHRM=NIFOUR
  170.       MINTE=IPMINT
  171.       NBBB=NBNN
  172. C
  173. C_______________________________________________________________________
  174. C
  175. C     numero des etiquettes      :
  176. C     etiquettes de 1 a 98 pour traitement specifique a l element
  177. C     dans la zone specifique a chaque element commencant par :
  178. C     5  continue
  179. C     element 5   etiquettes 1005 2005 3005 4005   ...
  180. C     44 continue
  181. C     element 44  etiquettes 1044 2044 3044 4044   ...
  182. C_______________________________________________________________________
  183. C
  184.       IF(MELE.GE.1.AND.MELE.LE.100) THEN
  185. C            CABL SEG2 SEG3 TRI3 TRI4 TRI6 TRI7 QUA4 QUA5 QUA8
  186.       GOTO (   99,  99,  99,   4,  99,   4,  99,   4,  99,   4
  187. C            QUA9 RAC2 RAC3 CUB8 CU20 PRI6 PR15 LIA3 LIA4 LIA6
  188.      1      ,  99,  99,  99,   4,   4,   4,   4,  99,  99,  99
  189. C            LIA8 MULT TET4 TE10 PYR5 PY13 COQ3  DKT POUT LISP
  190.      2      ,  99,  99,   4,   4,   4,   4,  99,  99,  99,  99
  191. C            FAC3 FAC4 FAC6 FAC8 LTR3 LQU4 LCU8 LPR6 LTE4 LPY5
  192.      3      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  193. C            COQ8 TUYA TUFI COQ2 POI1 BARR RACO LSU2 COQ4 LISM
  194.      4      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  195. C            COF3 RES2 LSU3 LSU4 LICO COQ6 CVS2 CVS3 CVT3 CVT6
  196.      5      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  197. C            CVQ4 CVQ8 THP5 TH13 THP6 TH15 THC8 TH20 ICT3 ICQ4
  198.      6      ,  99,  99,  99,  99,  99,  99,  99,  99,   4,   4
  199. C            ICT6 ICQ8 ICC8 ICT4 ICP6 IC20 IC10 IC15 TRIP QUAP
  200.      7      ,   4,   4,   4,   4,   4,   4,   4,   4,  79,  79
  201. C            CUBP TETP PRIP TIMO JOI2 JOI3 JOT3 JOI4 JOI6 JOI8
  202.      8      ,  79,  79,  79,  99,  99,  99,  99,  99,  99,  99
  203. C            LISC TRIH  DST LIC4 CERC TUYO LSE2 LITU HYT3 HYQ4
  204.      9      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99)
  205. c cccccc
  206.      .      ,MELE
  207.       ELSEIF(MELE.GE.101.AND.MELE.LE.200) THEN
  208. C            HYT4 HYP6 HYC8 TRIS QUAS POIS FOR3 JOP3 JOP6 JOP8
  209.       GOTO (   99,  99,  99,  99,  99,  99,  99,  80,  80,  80
  210. C            POL3 POL4 POL5 POL6 POL7 POL8 POL9 PO10 PO11 PO12
  211.      1      ,   4,   4,   4,   4,   4,   4,   4,   4,   4,   4
  212. C            PO13 PO14 BAR3 BAEX LIA2 QUAH CUBH ROT3 SEF2 TRF3
  213.      2      ,   4,   4,  99,  99,  99,  99,  99,  99,  99,  99
  214. C            QUF4 CUF8 PRF6 TEF4 PYF5 MSE3 MTR6 MQU9 MC27 MP18
  215.      3      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  216. C            MT10 MP14 SEF3 TRF7 QUF9 CF27 PF21 TF15 PF19 SEG6
  217.      4      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  218. C            TR21 QU36 C216 P126 TE56 PY91 TRH6 BSE2 BTR4 BQU5
  219.      5      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  220. C            BCU9 BPR7 BTE5 BPY6 FRO4 SEGS POJS JCT3 JCI4 JGI2
  221.      6      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  222. C            JGT3 JGI4 TRIQ QUAQ CUBQ TETQ PRIQ TRIR QUAR CUBR
  223.      7      ,  99,  99, 173, 173, 173, 173, 173, 173, 173, 173
  224. C            TETR PRIR Q4RI Q8RI JOQ3 JOQ6 JOQ8 JOR3 JOR6 JOR8
  225.      8      , 173, 173,   4,   4, 185, 185, 185, 185, 185, 185
  226. C            T1D2 T1D3 M1D2 M1D3 LC03 LC07 LC09 LC27 LC21 LC15
  227.      9      ,  99,  99,   4,   4,  99,  99,  99,  99,  99,  99)
  228. c cccccc
  229.      .      ,MELE-100
  230.       ELSEIF(MELE.GE.201.AND.MELE.LE.300) THEN
  231. C            LC19 LS03 LS07 LS09 LS27 LS21 LS15 LS19 BS03 BS07
  232.       GOTO (   99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  233. C            BS09 BS27 BS21 BS15 BS19 MC03 MC07 MC09 MC27 MC21
  234.      1      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  235. C            MC15 MC19 M103 M107 M109 M127 M121 M115 M119 MS03
  236.      2      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  237. C            MS07 MS09 MS27 MS21 MS15 MS19 QC03 QC07 QC09 QC27
  238.      3      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  239. C            QC21 QC15 QC19 Q103 Q107 Q109 Q127 Q121 Q115 Q119
  240.      4      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  241. C            QS03 QS07 QS09 QS27 QS21 QS15 QS19 CIFL SURE SHB8
  242.      5      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  243. C            CAF2 CAF3 XQ4R XC8R JOI1 ZCO2 ZCO3 ZCO4 TUY2 TUY3
  244.      6      ,  99,  99,  99,  99,  99,  99,  99,  99,  99,  99
  245. C            COS2 COA2 ICY5 IC13 CU27 PR21 TE15 PY19 C20R P15R
  246.      7      ,  99,  99,   4,   4,   4,   4,   4,   4,   4,   4)
  247. c cccccc
  248.      .      ,MELE-200
  249.       ENDIF
  250.       GOTO 99
  251. C
  252. C_______________________________________________________________________
  253. C
  254. C  elements massifs et elements incompressibles MECANIQUE
  255. C_______________________________________________________________________
  256. C
  257.   4   CONTINUE
  258.       IF (MFR.EQ.71 .OR. MFR.EQ.73) GOTO 97173
  259.       IF( ideri.le.2.or.ideri.eq.5 ) then
  260. C ideri le 2 est pour lineaire et quadratique  et =5 pour utlisateur
  261.  
  262. C Elements massifs en FORMULATION 'MECANIQUE'
  263. C -------------------------------------------
  264.       NBNO=NBNN
  265.       NDDD=NDEP-NDPGE
  266. C
  267. C     Donnees liees a l'element de reference
  268. C
  269.       SEGINI,WRK1
  270.       IF (Ideri.eq.2) SEGINI,MTRACE
  271. C
  272. C  boucle sur les elements
  273. C
  274.       DO 3004  IB=1,NBELEM
  275. C
  276. C     on cherche les deplacements
  277. C
  278.         MPTVAL=IVADEP
  279.         IE=1
  280.         DO IGAU=1,NBNN
  281.           DO ICOMP=1,NDDD
  282.             MELVAL=IVAL(ICOMP)
  283.             IGMN=MIN(IGAU,VELCHE(/1))
  284.             IBMN=MIN(IB  ,VELCHE(/2))
  285.             XDDL(IE)=VELCHE(IGMN,IBMN)
  286.             IE=IE+1
  287.           ENDDO
  288.         ENDDO
  289.         IF (NDPGE.GT.0) THEN
  290.           DO i=1,NDPGE
  291.             XDDL(IE)=UDPGE(i)
  292.             IE=IE+1
  293.           ENDDO
  294.         ENDIF
  295. C
  296. C     on cherche  les coordonnees des noeuds de l element ib
  297. C
  298.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  299.         if( ideri.eq.5) then
  300. C on se met a mi-pas
  301.           do iyu=1,nbnn
  302.             i_z = (iyu-1)*nddd
  303.             do iou=1,idim
  304.               XE(iou,iyu)= xe(iou,iyu) + xddl( iou + i_z )*0.5D0
  305.             enddo
  306.           enddo
  307.         endif
  308. C
  309. C     boucle sur les points de gauss
  310. C
  311.         ISDJC=0
  312. C
  313. C     Calcul des coeff de modification de b-barre (elts incompres)
  314. C= NOM  :   ICT3, ICQ4, ICT6, ICQ8, ICC8, ICT4, ICP6, IC20, IC10, IC15
  315. C= MELE :    69 ,  70 ,  71 ,  72 ,  73 ,  74 ,  75 ,  76 ,  77 ,  78
  316.         IF (MFR.EQ.31) THEN
  317.           CALL BBCALC(MELE,NBNN,MFR,IDIM,XE,
  318.      1                NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  319.      2                NSTRS,LRE,IFOUR,NHRM,A,BBX,SHPTOT,SHPWRK,
  320.      3                BGENE,XDPGE,YDPGE,PP)
  321.         ENDIF
  322.  
  323.         DO IGAU=1,NBPTEL
  324. C
  325.           CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  326.      1                MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,1.D0,XE,
  327.      2                SHPTOT,SHPWRK,BGENE,DJAC,XDPGE,YDPGE)
  328.  
  329.           IF (DJAC.EQ.0.D0) THEN
  330.             INTERR(1)=IB
  331.             if (noer.eq.0) CALL ERREUR(259)
  332.             kerr=259
  333.             GOTO 9904
  334.           ENDIF
  335.           IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  336.  
  337. C En cas d'elements incompressibles : BGENE selon la methode B-BARRE
  338.           IF (MFR.EQ.31) THEN
  339.             CALL BBAR(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  340.      &                MELE,NBNN,LRE,IFOUR,NSTRS,XE,DJAC,A,BBX,BGENE)
  341.           ENDIF
  342. C
  343.           CALL BST(BGENE,XDDL,LRE,NSTRS,XSTRS)
  344. C
  345. C     calcul des eps 2
  346. C
  347.           IF (Ideri.eq.2)
  348.      &      CALL BST2(SHPWRK,XDDL,XE,NBNO,IFOUR,XSTRS,TRACE,
  349.      &                IGAU,XDPGE,YDPGE,UDPGE,NHRM)
  350. C
  351. C     remplissage du segment contenant les deformations
  352. C
  353.           MPTVAL=IVAEPS
  354.           DO ICOMP=1,NSTRS
  355.             MELVAL=IVAL(ICOMP)
  356.             IBMN=MIN(IB,VELCHE(/2))
  357.             VELCHE(IGAU,IBMN)=XSTRS(ICOMP)
  358.           ENDDO
  359. C
  360.         ENDDO
  361. C
  362. C     fin de la boucle sur les points de gauss
  363. C
  364.         IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  365.          if (noer.eq.1) then
  366.           kerr=195
  367.          else
  368.           INTERR(1)=IB
  369.           CALL ERREUR(195)
  370.          endif
  371.           GOTO 9904
  372.         ENDIF
  373. C
  374. C   correction sur la partie quadratique de la deformation dans le cas
  375. C   des elements incompressibles
  376. C
  377.         IF (Ideri.eq.2) THEN
  378.           IF (MFR.EQ.31) THEN
  379.             CALL BBST2(TRACE,NBPTEL,IFOUR,MELE,POIGAU,QSIGAU,
  380.      &                 ETAGAU,DZEGAU,SHPTOT,NBNO,SHPWRK,XE,PP)
  381.             L=2
  382.             IF (IDIM.EQ.3 .OR. IFOUR.EQ.0) L=3
  383.             DO ICOMP=1,L
  384.               MELVAL=IVAL(ICOMP)
  385.               IBMN=MIN(IB  ,VELCHE(/2))
  386.               DO IGAU=1,NBPTEL
  387.                 VELCHE(IGAU,IBMN)=VELCHE(IGAU,IBMN)+TRACE(IGAU)
  388.               ENDDO
  389.             ENDDO
  390.           ENDIF
  391.         ENDIF
  392.  
  393.  3004 CONTINUE
  394. C
  395. C    fin de la boucle sur les elements
  396. C
  397.  9904 CONTINUE
  398.       SEGSUP WRK1
  399.       IF (IREPS2.EQ.1) SEGSUP MTRACE
  400.  
  401. C cas de la dérivée de Truesdell
  402.       elseif( ideri.eq.3) then
  403.       NBNO=NBNN
  404.       NDDD=NDEP-NDPGE
  405.       SEGINI,WRK1
  406. C      IF (IREPS2.EQ.1) SEGINI,MTRACE
  407.  
  408. C       on cherche le max des variations des champs pour savoir s'il faut
  409. C       appeler hookis plusieurs fois
  410.       mptval=ivamat
  411.       segact mptval
  412.       nbgmat=0
  413.       nelmat=0
  414.       DO IM=1,NMATT
  415.         VALMAT(IM) = 0.D0
  416.         IF (IVAL(IM).NE.0) THEN
  417.           MELVAL=IVAL(IM)
  418.           nelmat=Max(nelmat  ,VELCHE(/2))
  419.           nbgmat=Max(nbgmat,VELCHE(/1))
  420.         ENDIF
  421.       ENDDO
  422. C
  423. C  boucle sur les elements
  424. C
  425.       DO 34  IB=1,NBELEM
  426. C
  427. C     on cherche les deplacements
  428. C
  429.         MPTVAL=IVADEP
  430.         IE=1
  431.         DO IGAU=1,NBNN
  432.           DO ICOMP=1,NDDD
  433.             MELVAL=IVAL(ICOMP)
  434.             IGMN=MIN(IGAU,VELCHE(/1))
  435.             IBMN=MIN(IB  ,VELCHE(/2))
  436.             XDDL(IE)=VELCHE(IGMN,IBMN)
  437.             IE=IE+1
  438.           ENDDO
  439.         ENDDO
  440.         IF (NDPGE.GT.0) THEN
  441.           DO i=1,NDPGE
  442.             XDDL(IE)=UDPGE(i)
  443.             IE=IE+1
  444.           ENDDO
  445.         ENDIF
  446. C
  447. C     on cherche  les coordonnees des noeuds de l element ib
  448. C
  449.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  450. C on ajoute aux coordonnées la moitié du deplacements pour faire
  451. C la configuration à mi-pas
  452.         do iyu=1,idim
  453.           i_z = (iyu-1)*nddd
  454.           do iou=1,nbnn
  455.              XE(iou,iyu)= XE(iou,iyu) + xddl(iou+i_z)*0.5D0
  456.           enddo
  457.         enddo
  458. C
  459. C     boucle sur les points de gauss
  460. C
  461.         ISDJC=0
  462. C
  463.         DO 54  IGAU=1,NBPTEL
  464.           CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  465.      1                MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,1.D0,XE,
  466.      2                SHPTOT,SHPWRK,BGENE,DJAC,XDPGE,YDPGE)
  467. C
  468.           IF (DJAC.EQ.0.D0) THEN
  469.             INTERR(1)=IB
  470.             if (noer.eq.0) CALL ERREUR(259)
  471.             kerr=259
  472.             GOTO 994
  473.           ELSE IF (DJAC.LT.0.D0) THEN
  474.             ISDJC=ISDJC+1
  475.           ENDIF
  476. C
  477. C       on cherche les matrices de Hooke
  478. C
  479.           if(nbgmat+nelmat.gt.2 . or . ib+igau.eq.2) then
  480.             mptval=ivamat
  481.             DO IM=1,NMATT
  482.               IF (IVAL(IM).NE.0) THEN
  483.                 MELVAL=IVAL(IM)
  484.                 IBMN=MIN(IB  ,VELCHE(/2))
  485.                 IGMN=MIN(IGAU,VELCHE(/1))
  486.                 VALMAT(IM)=VELCHE(IGMN,IBMN)
  487.               ELSE
  488.                 VALMAT(IM)=0.D0
  489.               ENDIF
  490.             ENDDO
  491.             kcas=2
  492.             CALL HOOKIS(VALMAT,VALCAR,VAR,MFR,IB,IGAU,EXCEN,EPAIST,
  493.      +                  INAT,MELE,NPINT,IFOUR,KCAS,NBGMAT,Nelmat,
  494.      +                  S,SECT,LHOOK,DDHOMU,DDHOOK,
  495.      +                  COBMA,XMOB,IRETOU)
  496.           endif
  497.           do iou=1,nstrs
  498.             do iyu=1,nstrs
  499.              ddhomu(iyu,iou)=ddhook(iyu,iou)
  500.             enddo
  501.           enddo
  502.  
  503.           CALL DBST(BGENE,DDHomu,XDDL,LRE,NSTRS,XSTRS)
  504. C xstrs contient la contrainte on va faire pica xstrs zdep05
  505.           DO INO = 1, NBNN
  506.             i_z = (ino-1)*nddd
  507.             DO ID=1,IDIM
  508.               XE1(ID,INO)=XE(ID,INO)
  509.               XE2(ID,INO)=XE(ID,INO)-xddl( id + i_z )*0.5D0
  510.             ENDDO
  511.           ENDDO
  512.           DO IYU=1,3
  513.             DO iou=1,3
  514.               XJAC(iou,iyu)=0.D0
  515.             enddo
  516.           enddo
  517.           CALL HPRIME(XE1,NBNN,IDIM,SHPtot,IGAU,SHpwrk,DJAC)
  518. C
  519. C       CALCUL DE LA MATRICE  F
  520. C
  521.           DO IF=1,IDIM
  522.             DO IE=1,IDIM
  523.               R1 = 0.D0
  524.               DO ID=1,NBNN
  525.                 R1 = R1 + SHpwrk(IF+1,ID)*XE2(IE,ID)
  526.               ENDDO
  527.               XJAC(IE,IF) = R1
  528.             ENDDO
  529.           ENDDO
  530.         IF(IDIM.EQ.2) THEN
  531.            XJAC(3,3)=1.D0
  532.            IF(IFOUR.EQ.0) THEN
  533. C
  534. CC CAS AXISYMETRIQUE
  535. C
  536.              R1=0.D0
  537.              R2=0.D0
  538.              DO 150 ID=1,NBNN
  539.                R1=R1+SHpwrk(1,ID)*XE1(1,ID)
  540.                R2=R2+SHpwrk(1,ID)*XE2(1,ID)
  541.  150         CONTINUE
  542.              XJAC(3,3)=R2/(R1+1.D-20)
  543.            ENDIF
  544.         ENDIF
  545. CC  CALCUL DE DETERMINANT DE F
  546. C
  547.         IF(IDIM.EQ.3) THEN
  548.           DETF=XJAC(1,1)*(XJAC(2,2)*XJAC(3,3)-XJAC(3,2)*XJAC(2,3))
  549.           DETF=DETF-XJAC(2,1)*(XJAC(1,2)*XJAC(3,3)-XJAC(3,2)*XJAC(1,3))
  550.           DETF=DETF+XJAC(3,1)*(XJAC(1,2)*XJAC(2,3)-XJAC(1,3)*XJAC(2,2))
  551.         ELSE IF(IDIM.EQ.2) THEN
  552.           DETF=XJAC(1,1)*XJAC(2,2)-XJAC(1,2)*XJAC(2,1)
  553.           DETF = DETF * XJAC (3,3)
  554.         ENDIF
  555.         DETF=1.D0/(DETF+1.D-20)
  556. C
  557. C     CALCUL DES CONTRAINTES DE CAUCHY
  558. C
  559.         DO ID=1,NSTRS
  560.           IND=IN(ID)
  561.           JND=JN(ID)
  562.           R1=0.D0
  563.           DO IE=1,IDIM
  564.             DO IF=1,IDIM
  565.               ICO=ITAB(IE,IF)
  566.               R1 = R1 + XsTRS(ICO)*XJAC(IND,IE)*XJAC(JND,IF)
  567.             ENDDO
  568.           ENDDO
  569.           xstrs2(ID)= R1 * DETF
  570.         ENDDO
  571. C
  572. C PEGON  :  ON NE FAIT PAS LA TRANSFORMATION SUR LA 3-EME COMPOSANTE
  573. C
  574.        IF(IDIM.EQ.2) THEN
  575.          xstrs2(3)=xstrs(3)*XJAC(3,3)*XJAC(3,3)*DETF
  576.        ENDIF
  577. C fin du calcul de capi dans dans xstrs2 la contrainte de kirchoff
  578. C on continu en calculant epsi sur config initiale
  579.        DO INO=1,NBNN
  580.          i_z = (ino-1) * nddd
  581.          DO ID=1,IDIM
  582.            XE(ID,INO)=XE2(ID,INO)+xddl( id + i_z )*0.5D0
  583.          ENDDO
  584.        ENDDO
  585. C      inversion loi de hooke
  586.        CALL INVALM(DDHOMU,LHOOK,LHOOK,KERRE,0.D0)
  587.        DO I=1,LHOOK
  588.           R1 = 0.D0
  589.           DO J=1,LHOOK
  590.             R1 = R1 + DDHOMU(I,J)*xstrs2(J)
  591.           ENDDO
  592.           xstrs(I)= R1
  593.        ENDDO
  594. C
  595. C     remplissage du segment contenant les deformations
  596. C
  597.           MPTVAL=IVAEPS
  598.           DO ICOMP=1,NSTRS
  599.             MELVAL=IVAL(ICOMP)
  600.             IBMN=MIN(IB,VELCHE(/2))
  601.             VELCHE(IGAU,IBMN)=XSTRS(ICOMP)
  602.           ENDDO
  603.    54   continue
  604.  
  605.         IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  606.          if (noer.eq.1) then
  607.           kerr=195
  608.          else
  609.           INTERR(1)=IB
  610.           CALL ERREUR(195)
  611.           GOTO 994
  612.          endif
  613.         ENDIF
  614.    34 CONTINUE
  615.   994 CONTINUE
  616.       SEGSUP WRK1
  617. C fin du truesdell
  618. C debut du jaumann
  619.       elseif(ideri.eq.4) then
  620. C==========================
  621.       NBNO=NBNN
  622. C*    NDDD=NDEP
  623. C*    IF (IFOUR.EQ.-3) NDDD=NDEP-3
  624.       NDDD=NDEP-NDPGE
  625. C
  626.       SEGINI,WRK1,MTRACE,wrk2
  627.  
  628. C  boucle sur les elements
  629. C
  630.       DO 394  IB=1,NBELEM
  631. C
  632. C     on cherche les deplacements
  633. C
  634.         MPTVAL=IVADEP
  635.         IE=1
  636.         DO IGAU=1,NBNN
  637.           DO ICOMP=1,NDDD
  638.             MELVAL=IVAL(ICOMP)
  639.             IGMN=MIN(IGAU,VELCHE(/1))
  640.             IBMN=MIN(IB  ,VELCHE(/2))
  641.             XDDL(IE)=VELCHE(IGMN,IBMN)
  642.             IE=IE+1
  643.           ENDDO
  644.         ENDDO
  645.         IF (NDPGE.GT.0) THEN
  646.           DO i=1,NDPGE
  647.             XDDL(IE)=UDPGE(i)
  648.             IE=IE+1
  649.           ENDDO
  650.         ENDIF
  651. C
  652. C     on cherche  les coordonnees des noeuds de l element ib
  653. C
  654.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  655. C
  656. C   on se met sur la config à mi pas (XE) xe1 est la config initiale
  657. C
  658.         do iyu=1,nbnn
  659.           i_z = (iyu-1)*nddd
  660.           do iou=1,idim
  661.              XE1(iou,iyu)= xe(iou,iyu)
  662.              XE(iou,iyu)= xe(iou,iyu) + xddl( iou+ i_z )*0.5d0
  663.           enddo
  664.         enddo
  665. C
  666. C     boucle sur les points de gauss
  667. C
  668.         ISDJC=0
  669. C
  670. C     Calcul des coeff de modification de b-barre (elts incompres)
  671. C= NOM  :   ICT3, ICQ4, ICT6, ICQ8, ICC8, ICT4, ICP6, IC20, IC10, IC15
  672. C= MELE :    69 ,  70 ,  71 ,  72 ,  73 ,  74 ,  75 ,  76 ,  77 ,  78
  673.         IF (MFR.EQ.31) THEN
  674.           CALL BBCALC(MELE,NBNN,MFR,IDIM,XE,
  675.      1                NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  676.      2                NSTRS,LRE,IFOUR,NHRM,A,BBX,SHPTOT,SHPWRK,
  677.      3                BGENE,XDPGE,YDPGE,PP)
  678.         ENDIF
  679. C
  680.         DO 594  IGAU=1,NBPTEL
  681.           CALL BMATST(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  682.      1                MELE,MFR,NBNN,LRE,IFOUR,NSTRS,NHRM,1.D0,XE,
  683.      2                SHPTOT,SHPWRK,BGENE,DJAC,XDPGE,YDPGE)
  684.  
  685.           IF (DJAC.EQ.0.D0) THEN
  686.             INTERR(1)=IB
  687.             if (noer.eq.0) CALL ERREUR(259)
  688.             kerr=259
  689.             GOTO 9964
  690.           ENDIF
  691.           IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  692.  
  693. C En cas d'elements incompressibles : BGENE selon la methode B-BARRE
  694.           IF (MFR.EQ.31) THEN
  695.             CALL BBAR(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU,
  696.      &                MELE,NBNN,LRE,IFOUR,NSTRS,XE,DJAC,A,BBX,BGENE)
  697.           ENDIF
  698. C
  699.           CALL BST(BGENE,XDDL,LRE,NSTRS,XSTRS)
  700. C dans xstrs on a les deformations II sur config mi pas
  701. C on va calculer grad u/2 puis decomposition polaire puis rtens
  702. C on retravaille sur config initiale
  703.          r_z=XZero
  704.          iipdpg=0
  705.          CALL BGRMAS(iGau,NOELE,NBNO,LRE,IFOUR,NGRA,NIFOUR,XE1,
  706.      .                r_z,SHPTOT,SHPWRK,BB,BGR,DJAC,IIPDPG)
  707.          do iou=1,lre
  708.            xddls2(iou)= 0.5D0 * xddl(iou)
  709.          enddo
  710.          CALL BGRDEP(BGR,NGRA,XDDLs2,LRE,GRADI)
  711. C on ajoute l'identité au gradient
  712.          if(idim.eQ.2) then
  713.            gradi(1)=gradi(1)+1.D0
  714.            gradi(4)=gradi(4)+1.D0
  715.          ELSE IF(IDIM.EQ.3) THEN
  716.            gradi(1)=gradi(1)+1.D0
  717.            gradi(5)=gradi(5)+1.D0
  718.            gradi(9)=gradi(9)+1.D0
  719.          ENDIF
  720.  
  721.          CALL POLA2(gradi,R,U,IDIM)
  722. C fait le rtens Rt.A.R on utilise u pour mettre Rt
  723. C et on met le tenseur dans le tableau tens
  724. C attention vu le stockage R est en fait Rt
  725.          if(idim.eq.2) then
  726.            U(1)=r(1)
  727.            u(2)=r(3)
  728.            U(3)=R(2)
  729.            u(4)=R(4)
  730.            tens(1)=xstrs(1)
  731.            tens(2)=xstrs(4)*0.5d0
  732.            tens(3)=xstrs(4)*0.5d0
  733.            tens(4)=xstrs(2)
  734.  
  735.          elseif(idim.eq.3) then
  736.            U(1)=r(1)
  737.            u(2)=r(4)
  738.            U(3)=R(7)
  739.            u(4)=R(2)
  740.            u(5)=r(5)
  741.            u(6)=r(8)
  742.            u(7)=r(3)
  743.            u(8)=r(6)
  744.            u(9)=r(9)
  745.            tens(1)=xstrs(1)
  746.            tens(2)=xstrs(4)*0.5D0
  747.            tens(3)=xstrs(5)*0.5D0
  748.            tens(4)=tens(2)
  749.            tens(5)=xstrs(2)
  750.            tens(6)=xstrs(6)*0.5D0
  751.            tens(7)=tens(3)
  752.            tens(8)=tens(6)
  753.            tens(9)=xstrs(3)
  754.          else
  755.            write(6,*)' idim est ni 2 ni 3 stop'
  756.            stop
  757.          endif
  758.  
  759.          CALL MULMAT(tentra,tens,U,IDIM,IDIM,IDIM)
  760.          CALL MULMAT(tens,R,Tentra,IDIM,IDIM,IDIM)
  761. C  tens contient le nouveau tenseur  on  va remplir xstrs
  762. C en 2 D epzz ne change pas
  763.          if(idim.eq.2) then
  764.            xstrs(1)=tens(1)
  765.            xstrs(2)=tens(4)
  766.            xstrs(4)=tens(2)*2.D0
  767.          else
  768.            xstrs(1)=tens(1)
  769.            xstrs(2)=tens(5)
  770.            xstrs(3)=tens(9)
  771.            xstrs(4)=tens(2)*2.D0
  772.            xstrs(5)=tens(3)*2.D0
  773.            xstrs(6)=tens(6)*2.D0
  774.         endif
  775. C
  776. C     remplissage du segment contenant les deformations
  777. C
  778.           MPTVAL=IVAEPS
  779.           DO ICOMP=1,NSTRS
  780.             MELVAL=IVAL(ICOMP)
  781.             IBMN=MIN(IB,VELCHE(/2))
  782.             VELCHE(IGAU,IBMN)=XSTRS(ICOMP)
  783.          ENDDO
  784. C
  785.  594   CONTINUE
  786. C
  787. C     fin de la boucle sur les points de gauss
  788. C
  789.         IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  790.           INTERR(1)=IB
  791.          if (noer.eq.1) then
  792.           kerr=195
  793.          else
  794.           CALL ERREUR(195)
  795.           GOTO 9964
  796.          endif
  797.         ENDIF
  798.  
  799.  394   CONTINUE
  800. C
  801. C    fin de la boucle sur les elements
  802. C
  803.  9964  CONTINUE
  804.        SEGSUP WRK1,wrk2,MTRACE
  805.       endif
  806. C
  807.       GOTO 510
  808.  
  809. C Elements massifs en FORMULATIONs 'ELECTROSTATIQUE' et 'DIFFUSION'
  810. C -----------------------------------------------------------------
  811. 97173 CONTINUE
  812.       NBNO = NBNN
  813.       NDDD = NDEP
  814.       SEGINI,WRK1
  815. C-------------------------
  816. C Boucle sur les elements
  817. C-------------------------
  818.       DO IEL = 1, NBELEM
  819. C - Recuperation des coordonnees des noeuds de l element IEL
  820.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IEL,XE)
  821. C - Recuperation des inconnues primales aux noeuds de l element IEL
  822.         MPTVAL = IVADEP
  823.         IE = 1
  824.         DO IGAU = 1, NBNN
  825.           DO ICOMP = 1, NDDD
  826.             MELVAL = IVAL(ICOMP)
  827.             IGMN = MIN(IGAU,VELCHE(/1))
  828.             IEMN = MIN(IEL ,VELCHE(/2))
  829.             XDDL(IE) = VELCHE(IGMN,IEMN)
  830.             IE = IE+1
  831.           ENDDO
  832.         ENDDO
  833. C--  --  --  --  --  --  --  --  --
  834. C - Boucle sur les points de Gauss
  835. C--  --  --  --  --  --  --  --  --
  836.         ISDJC=0
  837.         DO IGAU = 1, NBPTEL
  838. C -- Calcul de la matrice B et du jacobien au point de Gauss IGAU
  839.           IF (MFR.EQ.71) THEN
  840.             CALL BELEC(XE,SHPTOT(1,1,IGAU),NBNN,LHOOK,-1,
  841.      &                 SHPWRK,BGENE,DJAC)
  842.           ELSE IF (MFR.EQ.73) THEN
  843.             CALL BDIFF(XE,SHPTOT(1,1,IGAU),NBNN,LHOOK,-1,
  844.      &                 SHPWRK,BGENE,DJAC)
  845.           ENDIF
  846.           IF (DJAC.EQ.0.D0) THEN
  847.             INTERR(1) = IEL
  848.             if (noer.eq.0) CALL ERREUR(259)
  849.             kerr=259
  850.             GOTO 98173
  851.           ENDIF
  852.           IF (DJAC.LT.0.D0) ISDJC = ISDJC+1
  853.           CALL BST(BGENE,XDDL,LRE,NSTRS,XSTRS)
  854. C -- Remplissage du segment contenant les "deformations" = -grad(.)
  855.           MPTVAL = IVAEPS
  856.           DO ICOMP = 1, NSTRS
  857.             MELVAL = IVAL(ICOMP)
  858.             IEMN = MIN(IEL,VELCHE(/2))
  859.             VELCHE(IGAU,IEMN) = XSTRS(ICOMP)
  860.           ENDDO
  861. C--  --  --  --  --  --  --  --  --
  862.         ENDDO
  863. C--  --  --  --  --  --  --  --  --
  864.         IF (ISDJC.NE.0 .AND. ISDJC.NE.NBPGAU) THEN
  865.           INTERR(1) = IEL
  866.          if (noer.eq.1) then
  867.           kerr=195
  868.          else
  869.           CALL ERREUR(195)
  870.           GOTO 98173
  871.          endif
  872.         ENDIF
  873. C-------------------------
  874.       ENDDO
  875. C-------------------------
  876. 98173 CONTINUE
  877.       SEGSUP,WRK1
  878.       GOTO 510
  879.  
  880. C_______________________________________________________________________
  881. C
  882. C     milieux poreux
  883. C_______________________________________________________________________
  884. C
  885.   79  CONTINUE
  886. C
  887. C    pour ces elements  nbbb = nombre de noeuds
  888. C                       nbno = nombre de fonctions de forme
  889. C
  890.       NBNO=IPORE
  891.       NSTN=1
  892.       LRN=NBNO-NBBB
  893.       LRB=LRE-LRN
  894. C
  895.       SEGINI WRK1,WRK5
  896. C  Initialisation de MTRACE necessaire mais inutilise pour ces elements
  897.       IF (IREPS2.EQ.1) SEGINI MTRACE
  898. C
  899.       DO 3079  IB=1,NBELEM
  900. C
  901. C     on cherche  les coordonnees des noeuds de l element ib
  902. C
  903.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  904. C
  905. C     on cherche les deplacements
  906. C
  907.         MPTVAL=IVADEP
  908.         IE=1
  909.         DO IGAU=1,NBNN
  910.           DO ICOMP=1,NDEP-1
  911.             MELVAL=IVAL(ICOMP)
  912.             IGMN=MIN(IGAU,VELCHE(/1))
  913.             IBMN=MIN(IB  ,VELCHE(/2))
  914.             XDDL(IE)=VELCHE(IGMN,IBMN)
  915.             IE=IE+1
  916.           ENDDO
  917.         ENDDO
  918. C
  919. C     puis  les pressions
  920. C
  921.         MELVAL=IVAL(NDEP)
  922.         IBMN=MIN(IB,VELCHE(/2))
  923.         DO IGAU=1,LRN
  924.           IGAUSO=IBSOM(NSPOS(IELE)+IGAU-1)
  925.           IGMN=MIN(IGAUSO,VELCHE(/1))
  926.           XDDL(IE)=VELCHE(IGMN,IBMN)
  927.           IE=IE+1
  928.         ENDDO
  929. C
  930. C     boucle sur les points de gauss
  931. C
  932.         ISDJC=0
  933. C
  934.         DO 5079 IGAU=1,NBPTEL
  935. C
  936.           CALL BNPORE(IGAU,NBNO,NBBB,LRE,IFOUR,LHOOK,NSTN,NHRM,
  937.      &                1.D0,XE,SHPTOT,SHPWRK,BGENE,XGENE,DJAC,1)
  938. C
  939.           IF (DJAC.EQ.0.D0) THEN
  940.             INTERR(1)=IB
  941.             if (noer.eq.0) CALL ERREUR(259)
  942.             kerr=259
  943.             GOTO 9979
  944.           ENDIF
  945.           IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  946. C
  947.           CALL BST(BGENE,XDDL,LRE,LHOOK,XSTRS)
  948. C
  949. C     calcul des eps 2
  950. C
  951.           IF (IREPS2.EQ.1)
  952.      &      CALL BST2(SHPWRK,XDDL,XE,NBNN,IFOUR,XSTRS,TRACE,
  953.      &                IGAU,XDPGE,YDPGE,UDPGE,NHRM)
  954. C
  955. C     calcul de la pression
  956. C
  957.           XP=0.D0
  958.           DO ID=1,LRN
  959.             XP=XP+XGENE(1,ID)*XDDL(LRB+ID)
  960.           ENDDO
  961.           XSTRS(NSTRS)=XP
  962. C
  963. C     remplissage du segment contenant les deformations
  964. C
  965.           MPTVAL=IVAEPS
  966.           DO 7079 ICOMP=1,NSTRS
  967.             MELVAL=IVAL(ICOMP)
  968.             IGMN=MIN(IGAU,VELCHE(/1))
  969.             IBMN=MIN(IB  ,VELCHE(/2))
  970.             VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  971.  7079     CONTINUE
  972. C
  973.  5079   CONTINUE
  974. C
  975. C     fin de la boucle sur les points de gauss
  976. C
  977.         IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  978.           INTERR(1)=IB
  979.          if (noer.eq.1) then
  980.           kerr=195
  981.          else
  982.           CALL ERREUR(195)
  983.           GOTO 9979
  984.          endif
  985.         ENDIF
  986. C
  987.  3079 CONTINUE
  988. C
  989. C    fin de la boucle sur les elements
  990. C
  991.  9979 CONTINUE
  992.       SEGSUP WRK1,WRK5
  993.       IF (IREPS2.EQ.1) SEGSUP MTRACE
  994. C
  995.       GOTO 510
  996. C_______________________________________________________________________
  997. C
  998. C     milieux poreux - SUITE
  999. C_______________________________________________________________________
  1000. C
  1001.  173  CONTINUE
  1002. C
  1003. C    pour ces elements  nbbb = nombre de noeuds
  1004. C                       nbno = nombre de fonctions de forme
  1005. C
  1006.       IF (MELE.GE.173.AND.MELE.LE.177) THEN
  1007.         IDECAP = 2
  1008.       ELSE IF (MELE.GE.178.AND.MELE.LE.182) THEN
  1009.         IDECAP = 3
  1010.       ENDIF
  1011. C
  1012.       NBNO=IPORE
  1013.       NSTN=IDECAP
  1014.       NSTB=4
  1015.       IF(IFOUR.EQ.1.OR.IFOUR.EQ.2) NSTB=6
  1016. C
  1017.       LPP=NBNO-NBBB
  1018.       LRN=IDECAP*LPP
  1019.       LRB=LRE-LRN
  1020. C
  1021.       SEGINI WRK1,WRK5
  1022. C  Initialise de MTRACE necessaire mais inutilise pour cet element
  1023.       IF (IREPS2.EQ.1) SEGINI MTRACE
  1024. C
  1025.       DO 3173  IB=1,NBELEM
  1026. C
  1027. C     on cherche  les coordonnees des noeuds de l element ib
  1028. C
  1029.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1030. C
  1031. C     on cherche les deplacements
  1032. C
  1033.         MPTVAL=IVADEP
  1034.         IE=1
  1035.         DO IGAU=1,NBNN
  1036.           DO ICOMP=1,NDEP-IDECAP
  1037.             MELVAL=IVAL(ICOMP)
  1038.             IGMN=MIN(IGAU,VELCHE(/1))
  1039.             IBMN=MIN(IB  ,VELCHE(/2))
  1040.             XDDL(IE)=VELCHE(IGMN,IBMN)
  1041.             IE=IE+1
  1042.           ENDDO
  1043.         ENDDO
  1044. C
  1045. C     puis  les pressions
  1046. C
  1047.         DO IPR = 1,IDECAP
  1048.           MELVAL=IVAL(NDEP-IDECAP+IPR)
  1049.           DO IGAU=1,LPP
  1050.             IGAUSO=IBSOM(NSPOS(IELE)+IGAU-1)
  1051.             IGMN=MIN(IGAUSO,VELCHE(/1))
  1052.             IBMN=MIN(IB  ,VELCHE(/2))
  1053.             XDDL(IE)=VELCHE(IGMN,IBMN)
  1054.             IE=IE+1
  1055.           ENDDO
  1056.         ENDDO
  1057. C
  1058. C     boucle sur les points de gauss
  1059. C
  1060.         ISDJC=0
  1061. C
  1062.         DO 5173  IGAU=1,NBPTEL
  1063. C
  1064.           CALL BNQORE(IGAU,NBNO,NBBB,LRE,IFOUR,NSTB,NSTN,NHRM,
  1065.      &        1.D0,XE,SHPTOT,SHPWRK,BGENE,XGENE,DJAC,IDECAP,LHOOK,1)
  1066. C
  1067.           IF (DJAC.EQ.0.D0) THEN
  1068.             INTERR(1)=IB
  1069.             if (noer.eq.0) CALL ERREUR(259)
  1070.             kerr=259
  1071.             GOTO 9173
  1072.           ENDIF
  1073.           IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  1074. C
  1075.           CALL BST(BGENE,XDDL,LRE,LHOOK,XSTRS)
  1076. C
  1077. C     calcul des eps 2
  1078. C
  1079.           IF (IREPS2.EQ.1)
  1080.      &      CALL BST2(SHPWRK,XDDL,XE,NBNN,IFOUR,XSTRS,TRACE,
  1081.      &                IGAU,XDPGE,YDPGE,UDPGE,NHRM)
  1082. C
  1083. C     calcul des  pressions
  1084. C
  1085.           IE=LRB
  1086.           DO IPR=1,IDECAP
  1087.             XP=0.D0
  1088.             IPR1=(IPR-1)*LPP
  1089.             DO ID=1,LPP
  1090.               IE=IE+1
  1091.               XP=XP+XGENE(IPR,ID+IPR1)*XDDL(IE)
  1092.             ENDDO
  1093.             XSTRS(NSTRS-IDECAP+IPR)=XP
  1094.           ENDDO
  1095. C
  1096. C     remplissage du segment contenant les deformations
  1097. C
  1098.           MPTVAL=IVAEPS
  1099.           DO 7173 ICOMP=1,NSTRS
  1100.             MELVAL=IVAL(ICOMP)
  1101.             IGMN=MIN(IGAU,VELCHE(/1))
  1102.             IBMN=MIN(IB  ,VELCHE(/2))
  1103.             VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  1104.  7173     CONTINUE
  1105. C
  1106.  5173    CONTINUE
  1107. C
  1108. C     fin de la boucle sur les points de gauss
  1109. C
  1110.          IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  1111.            INTERR(1)=IB
  1112.          if (noer.eq.1) then
  1113.           kerr=195
  1114.          else
  1115.            CALL ERREUR(195)
  1116.            GOTO 9173
  1117.          endif
  1118.          ENDIF
  1119. C
  1120.  3173 CONTINUE
  1121. C
  1122. C    fin de la boucle sur les elements
  1123. C
  1124.  9173 CONTINUE
  1125.       SEGSUP WRK1,WRK5
  1126.       IF (IREPS2.EQ.1) SEGSUP MTRACE
  1127. C
  1128.       GOTO 510
  1129.  
  1130. C_______________________________________________________________________
  1131. C
  1132. C     joints poreux
  1133. C_______________________________________________________________________
  1134. C
  1135.   80  CONTINUE
  1136. C
  1137. C    pour ces elements  nbbb = nombre de noeuds
  1138. C                       nbno = nombre de fonctions de forme
  1139. C
  1140.       NBNO=IPORE
  1141.       NSTN=1
  1142.       LRN=(NBNO-NBBB)*3/2
  1143.       LPP = LRN
  1144.       LRB=LRE-LRN
  1145.       NFAC=(3*NBBB-NBNO)/2
  1146. C
  1147.       SEGINI WRK1,WRK3,WRK5
  1148. C
  1149.       DO 3080  IB=1,NBELEM
  1150. C
  1151. C     on cherche d'abord les deplacements
  1152. C
  1153.         MPTVAL=IVADEP
  1154.         IE=1
  1155.         DO 4180 IGAU=1,NFAC
  1156.           DO 4280 ICOMP=1,NDEP-1
  1157.             MELVAL=IVAL(ICOMP)
  1158.             IGMN=MIN(IGAU,VELCHE(/1))
  1159.             IBMN=MIN(IB  ,VELCHE(/2))
  1160.             XDDL(IE)=VELCHE(IGMN,IBMN)
  1161.             IE=IE+1
  1162.  4280     CONTINUE
  1163.  4180   CONTINUE
  1164. C
  1165. C     puis les pressions
  1166. C
  1167.         MELVAL=IVAL(NDEP)
  1168.         DO 4080 IGAU=1,NBNN
  1169.           DO 4190 INSOM=1,NBSOM(IELE)
  1170.             IF (IGAU.EQ.IBSOM(NSPOS(IELE)+INSOM-1)) GO TO 4191
  1171. 4190      CONTINUE
  1172.           IF (IGAU.GT.NFAC) GO TO 4191
  1173.           GO TO 4080
  1174. 4191      CONTINUE
  1175.           IBMN=MIN(IB  ,VELCHE(/2))
  1176.           IGMN=MIN(IGAU,VELCHE(/1))
  1177.           XDDL(IE)=VELCHE(IGMN,IBMN)
  1178.           IE=IE+1
  1179.  4080   CONTINUE
  1180. C
  1181. C     on cherche  les coordonnees des noeuds de l element ib
  1182. C
  1183.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1184. C
  1185. C     calcul des exes locaux et des coordonnees locales
  1186. C
  1187.         CALL JOPLOC(XE,SHPTOT,NBBB,NBNO,IFOUR,XEL,BPSS)
  1188. C
  1189.         CALL INTDEL(XNTH,XNTB,XNTT,LRN,MELE)
  1190. C
  1191. C     boucle sur les points de gauss
  1192. C
  1193.         ISDJC=0
  1194. C
  1195.         DO 5080  IGAU=1,NBPTEL
  1196. C
  1197.           CALL BNPORJ(IGAU,NBNO,NBBB,LRE,IFOUR,LHOOK,NSTN,XE,XEL,
  1198.      &                SHPTOT,SHPWRK,BPSS,BGENE,XGENE,DJAC,1)
  1199. C
  1200.           IF (DJAC.EQ.0.D0) THEN
  1201.             INTERR(1)=IB
  1202.             if (noer.eq.0) CALL ERREUR(259)
  1203.             kerr=259
  1204.             GOTO 9980
  1205.           ENDIF
  1206.           IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  1207. C
  1208.           CALL BST(BGENE,XDDL,LRB,LHOOK,XSTRS)
  1209.  
  1210. C
  1211. C     calcul de la pression
  1212. C
  1213.           XP=0.D0
  1214.           DO 4480  ID=1,LRN
  1215.             XP=XP+XNTT(ID)*XGENE(1,ID)*XDDL(LRB+ID)
  1216.  4480     CONTINUE
  1217.           XSTRS(NSTRS)=XP
  1218. C
  1219. C     remplissage du segment contenant les deformations
  1220. C
  1221.           MPTVAL=IVAEPS
  1222.           DO 7080 ICOMP=1,NSTRS
  1223.             MELVAL=IVAL(ICOMP)
  1224.             IGMN=MIN(IGAU,VELCHE(/1))
  1225.             IBMN=MIN(IB  ,VELCHE(/2))
  1226.             VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  1227.  7080     CONTINUE
  1228. C
  1229.  5080   CONTINUE
  1230. C
  1231.         IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  1232.           INTERR(1)=IB
  1233.          if (noer.eq.1) then
  1234.           kerr=195
  1235.          else
  1236.           CALL ERREUR(195)
  1237.           GOTO 9980
  1238.          endif
  1239.         ENDIF
  1240. C
  1241.  3080 CONTINUE
  1242. C
  1243.  9980 CONTINUE
  1244.       SEGSUP WRK1,WRK3,WRK5
  1245. C
  1246.       GOTO 510
  1247.  
  1248.  
  1249. C_______________________________________________________________________
  1250. C
  1251. C     joints poreux  - SUITE
  1252. C_______________________________________________________________________
  1253. C
  1254.  185  CONTINUE
  1255. C
  1256. C    pour ces elements  nbbb = nombre de noeuds
  1257. C                       nbno = nombre de fonctions de forme
  1258. C
  1259.       IF (MELE.GE.185.AND.MELE.LE.187) THEN
  1260.         IDECAP = 2
  1261.       ELSE IF (MELE.GE.188.AND.MELE.LE.190) THEN
  1262.         IDECAP = 3
  1263.       ENDIF
  1264.  
  1265. C
  1266.       NBNO=IPORE
  1267.       NSTN=IDECAP
  1268.       NSTB=2
  1269.       IF(IFOUR.EQ.1.OR.IFOUR.EQ.2) NSTB=3
  1270.  
  1271. C
  1272.       LPP=(NBNO-NBBB)*3/2
  1273.       LRN=IDECAP*LPP
  1274.       LRB=LRE-LRN
  1275.  
  1276.       NFAC=(3*NBBB-NBNO)/2
  1277. C
  1278.       SEGINI WRK1,WRK3,WRK5
  1279. C
  1280.       DO 3185  IB=1,NBELEM
  1281. C
  1282. C     on cherche d'abord les deplacements
  1283. C
  1284.         MPTVAL=IVADEP
  1285.         IE=1
  1286.         DO 4185 IGAU=1,NFAC
  1287.           DO 4285 ICOMP=1,NDEP-IDECAP
  1288.             MELVAL=IVAL(ICOMP)
  1289.             IGMN=MIN(IGAU,VELCHE(/1))
  1290.             IBMN=MIN(IB  ,VELCHE(/2))
  1291.             XDDL(IE)=VELCHE(IGMN,IBMN)
  1292.             IE=IE+1
  1293.  4285     CONTINUE
  1294.  4185   CONTINUE
  1295. C
  1296. C     puis les pressions
  1297. C
  1298.         DO 4785 IPR=1,IDECAP
  1299.           MELVAL=IVAL(NDEP-IDECAP+IPR)
  1300.           DO 4085 IGAU=1,NBNN
  1301.             DO 4195 INSOM=1,NBSOM(IELE)
  1302.               IF (IGAU.EQ.IBSOM(NSPOS(IELE)+INSOM-1)) GO TO 4891
  1303. 4195        CONTINUE
  1304.             IF (IGAU.GT.NFAC) GO TO 4891
  1305.             GO TO 4085
  1306. 4891        CONTINUE
  1307.             IBMN=MIN(IB  ,VELCHE(/2))
  1308.             IGMN=MIN(IGAU,VELCHE(/1))
  1309.             XDDL(IE)=VELCHE(IGMN,IBMN)
  1310.             IE=IE+1
  1311.  4085     CONTINUE
  1312. 4785    CONTINUE
  1313. C
  1314. C     on cherche  les coordonnees des noeuds de l element ib
  1315. C
  1316.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1317. C
  1318. C     calcul des exes locaux et des coordonnees locales
  1319. C
  1320.         CALL JOPLOC(XE,SHPTOT,NBBB,NBNO,IFOUR,XEL,BPSS)
  1321. C
  1322.         CALL INTDEL(XNTH,XNTB,XNTT,LPP,MELE)
  1323. C
  1324. C     boucle sur les points de gauss
  1325. C
  1326.         ISDJC=0
  1327. C
  1328.         DO 5185  IGAU=1,NBPTEL
  1329. C
  1330.           CALL BNPQRJ(IGAU,NBNO,NBBB,LRE,IFOUR,LHOOK,NSTN,XE,XEL,
  1331.      &           SHPTOT,SHPWRK,BPSS,BGENE,XGENE,DJAC,IDECAP,NSTB,1)
  1332. C
  1333.           IF (DJAC.EQ.0.D0) THEN
  1334.             INTERR(1)=IB
  1335.             if (noer.eq.0) CALL ERREUR(259)
  1336.             kerr=259
  1337.             GOTO 9985
  1338.           ENDIF
  1339.           IF (DJAC.LT.0.D0) ISDJC=ISDJC+1
  1340. C
  1341.           CALL BST(BGENE,XDDL,LRB,LHOOK,XSTRS)
  1342. C
  1343. C     calcul de la pression
  1344. C
  1345.           IE=LRB
  1346.           DO 4985 IPR=1,IDECAP
  1347.             XP=0.D0
  1348.             IPR1=(IPR-1)*LPP
  1349.             DO 4485  ID=1,LPP
  1350.               IE=IE+1
  1351.               XP=XP+XNTT(ID)*XGENE(IPR,ID+IPR1)*XDDL(IE)
  1352.  4485       CONTINUE
  1353.             XSTRS(NSTRS-IDECAP+IPR)=XP
  1354.  4985     CONTINUE
  1355. C
  1356. C     remplissage du segment contenant les deformations
  1357. C
  1358.           MPTVAL=IVAEPS
  1359.           DO 7185 ICOMP=1,NSTRS
  1360.             MELVAL=IVAL(ICOMP)
  1361.             IGMN=MIN(IGAU,VELCHE(/1))
  1362.             IBMN=MIN(IB  ,VELCHE(/2))
  1363.             VELCHE(IGMN,IBMN)=XSTRS(ICOMP)
  1364.  7185     CONTINUE
  1365. C
  1366.  5185   CONTINUE
  1367. C
  1368.         IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  1369.           INTERR(1)=IB
  1370.          if (noer.eq.1) then
  1371.           kerr=195
  1372.          else
  1373.           CALL ERREUR(195)
  1374.           GOTO 9985
  1375.          endif
  1376.         ENDIF
  1377. C
  1378.  3185 CONTINUE
  1379. C
  1380.  9985 CONTINUE
  1381.       SEGSUP WRK1,WRK3,WRK5
  1382. C
  1383.       GOTO 510
  1384. C____________________________________________________________________
  1385.  99   CONTINUE
  1386.       MOTERR(1:4)=NOMTP(MELE)
  1387.       MOTERR(9:12)='EPSI'
  1388.       CALL ERREUR(86)
  1389.  
  1390.  510  CONTINUE
  1391. c      RETURN
  1392.       END
  1393.  
  1394.  
  1395.  
  1396.  

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