Télécharger la plaquette Cast3M


fr en

Télécharger jacono.eso

Retour à la liste

Numérotation des lignes : 

jacono

  1. C JACONO    SOURCE    MB234859  25/09/08    21:15:41     12358          
  2.  
  3. C=======================================================================
  4. C    ENTREES :
  5. C    ---------
  6. C      IPMODL= pointeur sur un MMODEL
  7. C      INORM = 1 si les vecteurs doivent etre normes 0 sinon
  8. C
  9. C    SORTIES :
  10. C    --------
  11. C
  12. C      IPCHE = CHAMELEM contenant les JACOBIENS
  13. C          ( = normale aux faces des elements dans le cas des coques)
  14. C          ( = tangente a la fibre neutre dans le cas des poutres)
  15. C      IRET  = 1 si succes 0 sinon
  16. C
  17. C    Passage au nouveau Chamelem PAR S.RAMAHANDRY le 11/09/90
  18. C
  19. C    2013-01-02 (BP) : ajout zones cohesives (ZCO2,3 et 4 => coque mince)
  20. C                      + calcul de la tangente pour les poutres
  21. C
  22. C
  23. C=====================================================================
  24.       SUBROUTINE JACONO(IPMODL,INORM,IPCHE,IRET)
  25.  
  26.       IMPLICIT INTEGER(I-N)
  27.       IMPLICIT REAL*8(A-H,O-Z)
  28.  
  29. -INC PPARAM
  30. -INC CCOPTIO
  31. -INC CCHAMP
  32.  
  33. -INC SMCHAML
  34. -INC SMMODEL
  35. -INC SMELEME
  36. -INC SMCOORD
  37. -INC SMINTE
  38.  
  39. -INC TMPTVAL
  40.  
  41.       SEGMENT TRA
  42.          REAL*8 XEL(3,NBNN) ,SHP(6,NBNN) ,XE(3,NBNN)
  43.       ENDSEGMENT
  44. C
  45.       SEGMENT TR1
  46.          REAL*8 TH(NBN1) ,TXR(3,3,NBN1) ,XJ(3,3)
  47.       ENDSEGMENT
  48. C
  49.       PARAMETER(UN=1.D0,XZER=0.D0)
  50.       DIMENSION BPSS(3,3)
  51.  
  52.       DIMENSION XU(3), XV(3), XW(3)
  53.  
  54.       IDIMP1 = IDIM+1
  55.       NHRM=NIFOUR
  56.       IRET=1
  57. C
  58. C     ACTIVATION DU MODELE
  59. C
  60.       MMODEL= IPMODL
  61.       SEGACT MMODEL
  62.       NSOUS=KMODEL(/1)
  63. C
  64. C     CREATION DU MCHELM
  65. C
  66.       N1=NSOUS
  67.       N3=6
  68.       IF (INORM .EQ. 1) THEN
  69.         L1=8
  70.       ELSE
  71.         L1=16
  72.       ENDIF
  73.       SEGINI MCHELM
  74.       IF (INORM .EQ. 1) THEN
  75.          TITCHE='NORMALES'
  76.       ELSE
  77.          TITCHE='VECTEURS SURFACE'
  78.       ENDIF
  79.       IFOCHE=IFOUR
  80.       IPCHE=MCHELM
  81. C____________________________________________________________________
  82. C
  83. C     DEBUT DE LA BOUCLE SUR LES DIFFERENTES ZONES
  84. C____________________________________________________________________
  85. C
  86.       DO 500 ISOUS=1,NSOUS
  87. C
  88. C     ON RECUPERE L INFORMATION GENERALE
  89. C
  90.          IMODEL=KMODEL(ISOUS)
  91.          SEGACT IMODEL
  92.          IPMAIL=IMAMOD
  93.          IMACHE(ISOUS)=IPMAIL
  94.          CONCHE(ISOUS)=CONMOD
  95. C
  96. C     TRAITEMENT DU MODELE
  97. C
  98.          MELE=NEFMOD
  99.          MELEME=IMAMOD
  100.          NFOR=FORMOD(/2)
  101.          NMAT=MATMOD(/2)
  102. C____________________________________________________________________
  103. C
  104. C     INFORMATION SUR L'ELEMENT FINI
  105. C____________________________________________________________________
  106. C
  107.          MELE =INFELE(1)
  108.          MFR  =INFELE(13)
  109.          MINTE=INFMOD(7)
  110.          MINTE1=INFMOD(3)
  111. C
  112.          INFCHE(ISOUS,1)=0
  113.          INFCHE(ISOUS,2)=0
  114.          INFCHE(ISOUS,3)=NHRM
  115.          INFCHE(ISOUS,4)=MINTE
  116.          INFCHE(ISOUS,5)=0
  117.          INFCHE(ISOUS,6)=5
  118. C
  119. C     INITIALISATION DE MINTE
  120. C
  121.          SEGACT MINTE
  122.          NBPGAU=POIGAU(/1)
  123. C
  124. C     ACTIVATION DU MELEME
  125. C
  126.          SEGACT MELEME
  127.          NBNN  =NUM(/1)
  128.          NBELEM=NUM(/2)
  129. C
  130. C     RECHERCHE DE LA TAILLE DES MELVAL A ALLOUER
  131. C
  132.          IF (MFR.EQ.3.OR.MFR.EQ.5.OR.MFR.EQ.9 .OR. MFR.EQ.77) THEN
  133.             N1PTEL=NBPGAU
  134.             N1EL=NBELEM
  135.          ELSEIF(MFR.EQ.7 .OR. MFR.EQ.13) THEN
  136.             N1PTEL=NBPGAU
  137.             N1EL=NBELEM
  138.          ELSE
  139.             N1PTEL = 0
  140.             N1EL = 0
  141.          ENDIF
  142.          N2PTEL=0
  143.          N2EL  =0
  144. C
  145. C     CREATION DU MCHAML DE LA SOUS ZONE
  146. C
  147.          NJAC=IDIM
  148.          N2 = NJAC
  149.          SEGINI MCHAML
  150.          ICHAML(ISOUS)=MCHAML
  151.          NSR   =1
  152.          NCOSOR=NJAC
  153.          SEGINI MPTVAL
  154.          IVAJAC=MPTVAL
  155. C
  156. C     2 OU 3 COMPOSANTES
  157. C
  158.          ICOMP=1
  159.          IF (IFOUR.EQ.0.OR.IFOUR.EQ.1) THEN
  160.             NOMCHE(ICOMP)='VR  '
  161.          ELSE
  162.             NOMCHE(ICOMP)='VX  '
  163.          ENDIF
  164.          TYPCHE(ICOMP)='REAL*8'
  165.          SEGINI MELVA1
  166.          IELVAL(ICOMP)=MELVA1
  167.          IVAL(ICOMP)=MELVA1
  168. C
  169.          ICOMP=2
  170.          IF (IFOUR.EQ.0.OR.IFOUR.EQ.1) THEN
  171.             NOMCHE(ICOMP)='VZ  '
  172.          ELSE
  173.             NOMCHE(ICOMP)='VY  '
  174.          ENDIF
  175.          TYPCHE(ICOMP)='REAL*8'
  176.          SEGINI MELVA2
  177.          IELVAL(ICOMP)=MELVA2
  178.          IVAL(ICOMP)=MELVA2
  179. C
  180.          MELVA3 = 0
  181.          IF (IDIM .EQ. 3) THEN
  182.             ICOMP=3
  183.             NOMCHE(ICOMP)='VZ  '
  184.             TYPCHE(ICOMP)='REAL*8'
  185.             SEGINI MELVA3
  186.             IELVAL(ICOMP)=MELVA3
  187.             IVAL(ICOMP)=MELVA3
  188.          ENDIF
  189. C
  190.          SEGINI TRA
  191. C
  192. C ================ FORMULATION MASSIVE =======================
  193. C
  194.          IF(MFR.EQ.1.OR.MFR.EQ.33) THEN
  195.             GOTO 520
  196. C
  197. C ================ FORMULATION COQUE  MINCE =====================
  198. C
  199. C
  200.          ELSE IF(MFR.EQ.3.OR.MFR.EQ.9 .OR. MFR.EQ.77) THEN
  201.           IDI2=IDIM-1
  202.           DO 3000 IB=1,NBELEM
  203. *--------------Calcul de la normale a l'élément
  204.             IREF1 = IDIMP1*(MELEME.NUM(1, IB) - 1)
  205.             IREF2 = IDIMP1*(MELEME.NUM(2, IB) - 1)
  206.             XNORU = 0.D0
  207.             DO IC = 1, IDIM
  208.                r_z = XCOOR(IREF2+IC)-XCOOR(IREF1+IC)
  209.                XU(IC) = r_z
  210.                XNORU = XNORU + (r_z * r_z)
  211.             ENDDO
  212.             XNORU = SQRT(XNORU)
  213.             DO IC = 1, IDIM
  214.                XU(IC) = XU(IC)/XNORU
  215.             ENDDO
  216.             IF (IDIM .EQ. 2) THEN
  217.                XW(1) = -XU(2)
  218.                XW(2) = XU(1)
  219.             ELSE
  220.                IN = 3
  221.    33          IREF3 = IDIMP1*(MELEME.NUM(IN, IB) - 1)
  222.                XNORV = 0.D0
  223.                DO IC = 1, IDIM
  224.                   r_z = XCOOR(IREF3+IC)-XCOOR(IREF1+IC)
  225.                   XV(IC) = r_z
  226.                   XNORV = XNORV + (r_z * r_z)
  227.                ENDDO
  228.                XNORV = SQRT(XNORV)
  229.                DO IC = 1, IDIM
  230.                  XV(IC) = XV(IC)/XNORV
  231.                ENDDO
  232.                XW(1) = XU(2)*XV(3) - XU(3)*XV(2)
  233.                XW(2) = XU(3)*XV(1) - XU(1)*XV(3)
  234.                XW(3) = XU(1)*XV(2) - XU(2)*XV(1)
  235.                XNORW = 0.
  236.                DO IC = 1, IDIM
  237.                  XNORW = XNORW + XW(IC)*XW(IC)
  238.                ENDDO
  239.                IF (XNORW .LT. 1.E-4) THEN
  240.                   IN = IN + 1
  241.                   if(IN.le.NBNN) GOTO 33
  242.                ENDIF
  243.                XNORW = SQRT(XNORW)
  244.                IF(XNORW .LT.1.E-4) call erreur(345)
  245.                DO IC = 1, IDIM
  246.                   XW(IC) = XW(IC)/XNORW
  247.                ENDDO
  248.             ENDIF
  249. *--------------Fin du calcul de la normale a l'élément
  250.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  251. C
  252.             IF(IDIM.EQ.2)THEN
  253.                CALL VPAST2(XE,BPSS)
  254.             ELSE IF(IDIM.EQ.3) THEN
  255.                CALL VPAST(XE,BPSS)
  256.             ENDIF
  257.             CALL VCORL1(XE,XEL,BPSS,NBNN)
  258.             MPTVAL=IVAJAC
  259.             DO 3002 IC=1,NBPGAU
  260.                IF (INORM .EQ. 1) THEN
  261.                   DJAC = 1.D0
  262.                ELSE
  263.                   DO IE=1,NBNN
  264.                      DO ID=1,6
  265.                         SHP(ID,IE)=SHPTOT(ID,IE,IC)
  266.                      ENDDO
  267.                   ENDDO
  268.                   CALL JACOBI(XEL,SHP,IDI2,NBNN,DJAC)
  269.                ENDIF
  270.                MELVAL = IVAL(1)
  271.                IBMN=MIN(IB, MELVAL.VELCHE(/2))
  272.                IGMN=MIN(IC,MELVAL.VELCHE(/1))
  273.                MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(1)
  274.                MELVAL = IVAL(2)
  275.                IBMN=MIN(IB, MELVAL.VELCHE(/2))
  276.                IGMN=MIN(IC,MELVAL.VELCHE(/1))
  277.                MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(2)
  278.                IF (IDIM .EQ. 3) THEN
  279.                   MELVAL = IVAL(3)
  280.                   IBMN=MIN(IB, MELVAL.VELCHE(/2))
  281.                   IGMN=MIN(IC,MELVAL.VELCHE(/1))
  282.                   MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(3)
  283.                ENDIF
  284.  3002       CONTINUE
  285.  3000    CONTINUE
  286.          GOTO 520
  287. C
  288. C ================ FORMULATION POUTRE ET TUYAU ====================
  289. C
  290.       ELSE IF(MFR.EQ.7.OR.MFR.EQ.13) THEN
  291.          IDI2=IDIM-1
  292.          DO 4000 IB=1,NBELEM
  293. *-----------Calcul de la tangente a l'élément
  294.             IREF1 = IDIMP1*(MELEME.NUM(1, IB) - 1)
  295.             IREF2 = IDIMP1*(MELEME.NUM(2, IB) - 1)
  296.             XNORU = 0.D0
  297.             DO IC = 1, IDIM
  298.                r_z = XCOOR(IREF2+IC)-XCOOR(IREF1+IC)
  299.                XU(IC) = r_z
  300.                XNORU = XNORU + (r_z * r_z)
  301.             ENDDO
  302.             XNORU = SQRT(XNORU)
  303.             DO IC = 1, IDIM
  304.                XU(IC) = XU(IC)/XNORU
  305.             ENDDO
  306. *-----------Fin du calcul de la tangente a l'élément
  307. *-----------Calcul de la tangente a l'élément en chaque point de Gauss
  308. c  BP : On suppose le jacobien constant dans l'element (idem POUJAC.eso)
  309. C  => on sort le calcul du jacobien de la boucle sur les points de Gauss.
  310. C  Cela implique que la POUTre de Bernoulli n'est pas isoparamétrique...
  311.             IF (INORM .EQ. 1) THEN
  312.                DJAC = 1.D0
  313.             ELSE
  314.                DJAC = 1.D0/DBLE(NBPGAU)
  315.             ENDIF
  316.             MPTVAL=IVAJAC
  317.             DO 4002 IC=1,NBPGAU
  318.                MELVAL = IVAL(1)
  319.                IBMN=MIN(IB, MELVAL.VELCHE(/2))
  320.                IGMN=MIN(IC,MELVAL.VELCHE(/1))
  321.                MELVAL.VELCHE(IGMN,IBMN)=DJAC*XU(1)
  322.                MELVAL = IVAL(2)
  323.                IBMN=MIN(IB, MELVAL.VELCHE(/2))
  324.                IGMN=MIN(IC,MELVAL.VELCHE(/1))
  325.                MELVAL.VELCHE(IGMN,IBMN)=DJAC*XU(2)
  326.                IF (IDIM .EQ. 3) THEN
  327.                   MELVAL = IVAL(3)
  328.                   IBMN=MIN(IB, MELVAL.VELCHE(/2))
  329.                   IGMN=MIN(IC,MELVAL.VELCHE(/1))
  330.                   MELVAL.VELCHE(IGMN,IBMN)=DJAC*XU(3)
  331.                ENDIF
  332.  4002       CONTINUE
  333.  4000    CONTINUE
  334.  
  335.          GOTO 520
  336. C
  337. C ================ FORMULATION COQUE  EPAISSE ====================
  338. C
  339.          ELSE IF(MFR.EQ.5) THEN
  340.             SEGACT MINTE1
  341.             NBPGA1=MINTE1.POIGAU(/1)
  342.             NBN1  =MINTE1.SHPTOT(/2)
  343.             SEGINI TR1
  344. C
  345. C     UNE PETITE HORREUR ON CONSIDERE LES EPAISSEURS CONSTANTES
  346. C
  347.             DO 5010 IC=1,NBNN
  348.                TH(IC)=UN
  349.  5010       CONTINUE
  350.             DO 5000 IB=1,NBELEM
  351. *--------------Calcul de la normale a l'élément
  352.                IREF1 = IDIMP1*(MELEME.NUM(1, IB) - 1)
  353. *               IREF2 = IDIMP1*(MELEME.NUM(2, IB) - 1)
  354. * bp : les EF de coque epaisse etant quadratiques (coq6 et coq8), on
  355. * prend les noeuds "coins" pour eviter pb avec les noeuds 1,2,3 si courbures
  356.                IREF2 = IDIMP1*(MELEME.NUM(3, IB) - 1)
  357.                XNORU = 0.
  358.                DO IC = 1, IDIM
  359.                   r_z = XCOOR(IREF2+IC)-XCOOR(IREF1+IC)
  360.                   XU(IC) = r_z
  361.                   XNORU = XNORU + (r_z * r_z)
  362.                ENDDO
  363.                XNORU = SQRT(XNORU)
  364.                DO IC = 1, IDIM
  365.                   XU(IC) = XU(IC)/XNORU
  366.                ENDDO
  367.                IF (IDIM .EQ. 2) THEN
  368.                   XW(1) = -XU(2)
  369.                   XW(2) = XU(1)
  370.                ELSE
  371. *                  IN = 3
  372.                   IN = 5
  373.    13             IREF3 = IDIMP1*(MELEME.NUM(IN, IB) - 1)
  374.                   XNORV = 0.D0
  375.                   DO IC = 1, IDIM
  376.                     r_z = XCOOR(IREF3 + IC) - XCOOR(IREF1 + IC)
  377.                     XV(IC) = r_z
  378.                     XNORV = XNORV + (r_z * r_z)
  379.                   ENDDO
  380.                   XNORV = SQRT(XNORV)
  381.                   DO IC = 1, IDIM
  382.                     XV(IC) = XV(IC)/XNORV
  383.                   ENDDO
  384.                   XW(1) = XU(2)*XV(3) - XU(3)*XV(2)
  385.                   XW(2) = XU(3)*XV(1) - XU(1)*XV(3)
  386.                   XW(3) = XU(1)*XV(2) - XU(2)*XV(1)
  387.                   XNORW = 0.
  388.                   DO IC = 1, IDIM
  389.                      XNORW = XNORW + XW(IC)*XW(IC)
  390.                   ENDDO
  391.                   XNORW = SQRT(XNORW)
  392.                   IF (XNORW .LT. 1.E-4) THEN
  393.                    if(IN.LT.NBNN) then
  394.                      IN = IN + 1
  395.                      GOTO 13
  396.                    else
  397.                     write(6,*) 'Difficultes pour etablir la normale de'
  398.                     write(6,*) 'l element',IB
  399.                     write(6,*) 'Verifiez votre maillage'
  400.                     GOTO 9990
  401.                    endif
  402.                   ENDIF
  403.                   DO IC = 1, IDIM
  404.                      XW(IC) = XW(IC)/XNORW
  405.                   ENDDO
  406.                ENDIF
  407. *--------------Fin du calcul de la normale a l'élément
  408.                IF (INORM .EQ. 0) THEN
  409.                   CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  410. C
  411.                   CALL CQ8LOC(XE,NBN1,MINTE1.SHPTOT,TXR,IRR)
  412.                ENDIF
  413. C
  414.                MPTVAL=IVAJAC
  415.                DO 5002 IC=1,NBPGAU
  416.                   IF (INORM .EQ. 1) THEN
  417.                      DJAC = 1.D0
  418.                   ELSE
  419.                      E=DZEGAU(IC)
  420.          CALL COQ8JC(IC,NBN1,E,XE,TH,TXR,SHPTOT,XJ,DJAC,IRR)
  421.                   ENDIF
  422.                   MELVAL = IVAL(1)
  423.                   IBMN=MIN(IB, MELVAL.VELCHE(/2))
  424.                   IGMN=MIN(IC, MELVAL.VELCHE(/1))
  425.                   MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(1)
  426.                   MELVAL = IVAL(2)
  427.                   IBMN=MIN(IB, MELVAL.VELCHE(/2))
  428.                   IGMN=MIN(IC, MELVAL.VELCHE(/1))
  429.                   MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(2)
  430.                   IF (IDIM .EQ. 3) THEN
  431.                      MELVAL = IVAL(3)
  432.                      IBMN=MIN(IB, MELVAL.VELCHE(/2))
  433.                      IGMN=MIN(IC, MELVAL.VELCHE(/1))
  434.                      MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(3)
  435.                   ENDIF
  436.  5002          CONTINUE
  437.  5000       CONTINUE
  438.             SEGSUP TR1
  439.             GOTO 520
  440.          ENDIF
  441. C---------------------------------------------------------------------
  442. C  DESACTIVATION DES SEGMENTS PROPRES A LA ZONE GEOMETRIQUE ISOUS
  443. C---------------------------------------------------------------------
  444.  520     CONTINUE
  445.          MPTVAL=IVAJAC
  446.          SEGSUP MPTVAL
  447.          SEGSUP TRA
  448.  
  449.   500 CONTINUE
  450.       RETURN
  451. C
  452. C-------------------------------------------------------------------
  453. C         ERREUR DANS UNE ZONE , DESACTIVATION ET RETOUR
  454. C-------------------------------------------------------------------
  455.  9990 CONTINUE
  456.       IRET = 0
  457.       MPTVAL=IVAJAC
  458.       SEGSUP MPTVAL
  459.  
  460. *      RETURN
  461.       END
  462.  
  463.  
  464.  
  465.  

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