Télécharger la plaquette Cast3M


fr en

Télécharger rot3r.eso

Retour à la liste

Numérotation des lignes : 

rot3r

  1. C ROT3R     SOURCE    OF166741  24/10/23    21:15:07     12046          
  2.  
  3. ************************************************************************
  4. *
  5. *                             R O T 3 R
  6. *                             ---------
  7. *
  8. * FONCTION:
  9. * ---------
  10. *     CALCUL DE LA MATRICE DE RESISTANCE POUR L'ELEMENT ROT3
  11. *
  12. * PARAMETRES:   (E)=ENTREE   (S)=SORTIE   (+ = CONTENU DANS UN COMMUN)
  13. * -----------
  14. *     NEF     (E)  NUMERO DE L'ELEMENT-FINI DANS NOMTP (VOIR CCHAMP)
  15. *     IPMAIL  (E)  MAILLAGE ELEMENTAIRE CONSIDERE
  16. *     IPMODE  (E)  POINTEUR SUR UN SEGMENT IMODEL
  17. *     IPCHEM  (E)  POINTEUR SUR LE CHAMELEM DE CARACTERISTIQUE
  18. *     IPMATR (E/S)  MATRICE DE RIGIDITE ELEMENTAIRE XMATRI
  19. ************************************************************************
  20.  
  21.       SUBROUTINE ROT3R(NEF,IPMAIL,IPMODE,IPCHEM,IPMATR)
  22.  
  23.       IMPLICIT INTEGER(I-N)
  24.       IMPLICIT REAL*8 (A-H,O-Z)
  25.  
  26. -INC PPARAM
  27. -INC CCOPTIO
  28. -INC CCREEL
  29. -INC CCHAMP
  30.  
  31. -INC SMCOORD
  32. -INC SMINTE
  33. -INC SMMODEL
  34. -INC SMRIGID
  35. -INC SMELEME
  36. -INC SMCHAML
  37.  
  38.       SEGMENT NOTYPE
  39.         CHARACTER*16 TYPE(NBTYPE)
  40.       ENDSEGMENT
  41.  
  42.       SEGMENT MPTVAL
  43.         INTEGER IPOS(NS)    ,NSOF(NS)
  44.         INTEGER      IVAL(NCOSOU)
  45.         CHARACTER*16 TYVAL(NCOSOU)
  46.       ENDSEGMENT
  47.  
  48.       SEGMENT,MMAT1
  49.         REAL*8 VALMAT(NMATR)
  50.         REAL*8 XE(3,NBNN),CEL1(NBNN,NBNN),XE1(3,NBNN)
  51.         REAL*8 SHP(6,NBNN),GRAD(NDIM,NBNN)
  52.       ENDSEGMENT
  53.  
  54.       REAL*8 COSD1(3),COSD2(3),COSD3(3),YK(2,2)
  55.  
  56.       CHARACTER*8 CNM
  57.       CHARACTER*(NCONCH) CONM
  58.       PARAMETER (NINF=3)
  59.       INTEGER INFOS(NINF)
  60.  
  61.       IMODEL = IPMODE
  62.       CONM   = imodel.CONMOD
  63.  
  64.       MELEME = IPMAIL
  65. c*      meleme = imodel.IMAMOD
  66.       NBNN   = meleme.NUM(/1)
  67.       NBELEM = meleme.NUM(/2)
  68.  
  69. *     RECUPERATION DES CARACTERISTIQUES D'INTEGRATION DE L'ELEMENT
  70. *     FINI LIE A NOTRE MAILLAGE
  71.       if (infmod(/1).lt.4) then
  72.         write(ioimp,*) 'rot3r infmod(/1)'
  73.         call erreur(5)
  74.       endif
  75.  
  76. *  INFORMATION SUR L'ELEMENT
  77.       MINTE = imodel.INFMOD(4)
  78.       NBPGAU = minte.POIGAU(/1)
  79.  
  80.       xMATRI = IPMATR
  81. c*      SEGACT,xMATRI*MOD
  82.       NLIGRP = NBNN
  83.       NLIGRD = NBNN
  84.  
  85. *  REMLIR LE TABLEAU INFOS (INFORMATIONS SUR ELEMENT)
  86.       INFOS(1)=0
  87.       INFOS(2)=0
  88.       INFOS(3)=NIFOUR
  89.  
  90. *  RECHERCHE LES POINTEURS DES SEGMENTS MELVAL QUI CORRESPONDENT
  91. *  AUX COMPOSANTES DE LA CONDUCTIVITE ET L'EPAISSEUR DES ELEMENT
  92.       CNM = imodel.CMATEE
  93. c*      INM = imodel.IMATEE
  94. c*      INT = imodel.INATUU
  95.  
  96.       nbrobl = 0
  97.       nbrfac = 0
  98.       nomid  = 0
  99.       MOMATR = nomid
  100.  
  101.       NBTYPE = 1
  102.       SEGINI,notype
  103.       TYPE(1) = 'REAL*8'
  104.       MOTYR8 = notype
  105.  
  106.       IF (CNM.EQ.'ISOTROPE') THEN
  107.         NBROBL=2
  108.         SEGINI,NOMID
  109.         LESOBL(1)='ETA '
  110.         LESOBL(2)='EPAI'
  111.       ELSE IF (CNM.EQ.'ORTHOTRO') THEN
  112.         NBROBL=5
  113.         SEGINI,NOMID
  114.         LESOBL(1)='ETA1'
  115.         LESOBL(2)='ETA2'
  116.         LESOBL(3)='EPAI'
  117.         LESOBL(4)='V1X '
  118.         LESOBL(5)='V1Y '
  119.       ELSE
  120.         CALL ERREUR(251)
  121.         RETURN
  122.       ENDIF
  123.       NMATO = nbrobl
  124.       NMATF = nbrfac
  125.       NMATR = NMATO + NMATF
  126.       MOMATR = nomid
  127.  
  128.       IVAMAT = 0
  129.       CALL KOMCHA(IPCHEM,IPMAIL,CONM,MOMATR,MOTYR8,1,INFOS,NINF,IVAMAT)
  130.       IF (IERR.NE.0) GOTO 990
  131.  
  132.       MPTVAL = IVAMAT
  133.       IF (CNM.EQ.'ISOTROPE')THEN
  134.         IPMELV = IVAL(2)
  135.       ELSE IF(CNM.EQ.'ORTHOTRO') THEN
  136.         IPMELV = IVAL(3)
  137.       ENDIF
  138.       CALL QUELCH(IPMELV,ICONS)
  139.       IF (ICONS.NE.0) THEN
  140.         CALL ERREUR(566)
  141.         GOTO 990
  142.       ENDIF
  143.  
  144.       NDIM = IDIM-1
  145.       NFIN = IDIM
  146.  
  147.       SEGINI,MMAT1
  148.  
  149. *   BOUCLE (10) SUR LES ELEMENTS DU MAILLAGE ELEMENTAIRE IMAIL
  150.  
  151.       DO IEL = 1, NBELEM
  152.  
  153. *     ON CHERCHE LES COORDONNEES DES NOEUDS DE L'ELEMENT IEL,
  154. *     DANS LE REPERE GLOBAL
  155.         CALL DOXE(XCOOR,IDIM,NBNN,NUM,IEL,XE)
  156.  
  157. *     CALCUL DES COORDONNEES DES NOEUDS DANS LE REPERE LOCAL DE L'
  158. *     ELEMENT COQUE
  159.         COSD1(1) = XE(1,2)-XE(1,1)
  160.         COSD1(2) = XE(2,2)-XE(2,1)
  161.         COSD1(3) = XE(3,2)-XE(3,1)
  162.  
  163.         COSD2(1) = XE(1,3)-XE(1,1)
  164.         COSD2(2) = XE(2,3)-XE(2,1)
  165.         COSD2(3) = XE(3,3)-XE(3,1)
  166.  
  167.         COSD3(1) = COSD1(2)*COSD2(3)-COSD1(3)*COSD2(2)
  168.         COSD3(2) = COSD1(3)*COSD2(1)-COSD1(1)*COSD2(3)
  169.         COSD3(3) = COSD1(1)*COSD2(2)-COSD1(2)*COSD2(1)
  170.  
  171.         COSD1L = SQRT(COSD1(1)*COSD1(1)+COSD1(2)*COSD1(2)+
  172.      &                COSD1(3)*COSD1(3))
  173.  
  174.         COSD1(1)=COSD1(1)/COSD1L
  175.         COSD1(2)=COSD1(2)/COSD1L
  176.         COSD1(3)=COSD1(3)/COSD1L
  177.  
  178.         COSD3L = SQRT(COSD3(1)*COSD3(1)+COSD3(2)*COSD3(2)+
  179.      &                COSD3(3)*COSD3(3))
  180.  
  181.         COSD3(1)=COSD3(1)/COSD3L
  182.         COSD3(2)=COSD3(2)/COSD3L
  183.         COSD3(3)=COSD3(3)/COSD3L
  184.  
  185.         COSD2(1) = COSD3(2)*COSD1(3)-COSD3(3)*COSD1(2)
  186.         COSD2(2) = COSD3(3)*COSD1(1)-COSD3(1)*COSD1(3)
  187.         COSD2(3) = COSD3(1)*COSD1(2)-COSD3(2)*COSD1(1)
  188.  
  189.         DO NOE = 1, NBNN
  190.           XE1(1,NOE) =  XE(1,NOE)*COSD1(1) + XE(2,NOE)*COSD1(2)
  191.      &                + XE(3,NOE)*COSD1(3)
  192.           XE1(2,NOE) =  XE(1,NOE)*COSD2(1) + XE(2,NOE)*COSD2(2)
  193.      &                + XE(3,NOE)*COSD2(3)
  194.           XE1(3,NOE) = 0.D0
  195.         ENDDO
  196.  
  197. *     MISE A ZERO DU TABLEAU CEL1
  198.         CALL ZERO(CEL1,NBNN,NBNN)
  199.  
  200. *     BOUCLE (20) SUR LES POINTS DE GAUSS
  201.         IFOIS = 0
  202.         IFOI2 = 0
  203.  
  204.         DO IGAU=1,NBPGAU
  205.  
  206. *       CALCUL DE LA MATRCIE GRADIENT DES FONCTIONS DE FORME ET
  207. *       DU JACOBIEN(DANS LE PLAN), EN UN POINT DE GAUSS
  208.           DO NOE = 1, NBNN
  209.             DO I = 1, NFIN
  210.               SHP(I,NOE) = SHPTOT(I,NOE,IGAU)
  211.             ENDDO
  212.           ENDDO
  213.  
  214. *       DERIVES DES FONCTIONS DE FORME DANS LA GEOMETRIE REELLE
  215. *       ET LE JACOBIEN
  216.           CALL JACOBI(XE1,SHP,NDIM,NBNN,DJAC)
  217.  
  218.           IF (DJAC.LT.XZERO) IFOIS=IFOIS+1
  219.           IF (ABS(DJAC).LT.XPETIT) IFOI2=IFOI2+1
  220.           DJAC=ABS(DJAC)*POIGAU(IGAU)
  221.  
  222.           DO NOE = 1, NBNN
  223.             DO I = 1, NDIM
  224.               GRAD(I,NOE) = SHP(I+1,NOE)
  225.             ENDDO
  226.           ENDDO
  227.  
  228. *       ON CHERCHE LES VALEURS DE COMPOSANTES DE LA RESISTIVITE
  229. *       ET L'EPAISSEUR DE LA COQUE
  230.           MPTVAL=IVAMAT
  231.           DO IM=1,NMATR
  232.             MELVAL=IVAL(IM)
  233.             IF (MELVAL.NE.0)THEN
  234.               IBMN=MIN(IEL,VELCHE(/2))
  235.               IGMN=MIN(IGAU,VELCHE(/1))
  236.               VALMAT(IM)=VELCHE(IGMN,IBMN)
  237.             ELSE
  238.              VALMAT(IM)=0.D0
  239.             ENDIF
  240.           ENDDO
  241.  
  242. *       MATERIAU ISOTROPE
  243.           IF (CNM.EQ.'ISOTROPE') THEN
  244.             EP = VALMAT(2)
  245. *   L'ELEMENT (IEL) AU POINT DE GAUSS (IGAU) A UNE EPAISSEUR NULLE
  246.             IF (EP.LE.XPETIT) THEN
  247.               INTERR(1)=IEL
  248.               INTERR(2)=IGAU
  249.               MOTERR(1:4)=NOMTP(NEF)
  250.               CALL ERREUR(355)
  251.               GO TO 999
  252.             ENDIF
  253.             XK = VALMAT(1)*DJAC/EP
  254. *        ON AJOUTE LE PRODUIT K*DJAC*TRANSPOSEE(GRAD)*GRAD
  255. *        POUR LE POINT DE GAUSS CONSIDERE,A LA MATRICE CEL1
  256.             CALL NTNST(GRAD,XK,NBNN,NDIM,CEL1)
  257.  
  258. *        CAS ORTHOTROPE
  259.           ELSE
  260. c*            IF (CNM.EQ.'ORTHOTRO')THEN
  261.  
  262.             EP = VALMAT(3)
  263. *   L'ELEMENT (IEL) AU POINT DE GAUSS (IGAU) A UNE EPAISSEUR NULLE
  264.             IF (EP.LE.XPETIT) THEN
  265.               INTERR(1)=IEL
  266.               INTERR(2)=IGAU
  267.               MOTERR(1:4)=NOMTP(NEF)
  268.               CALL ERREUR(355)
  269.               GO TO 999
  270.             ENDIF
  271.             XK1 = VALMAT(1) / EP
  272.             XK2 = VALMAT(2) / EP
  273.  
  274.             COSA =  VALMAT(4)
  275.             SINA = -VALMAT(5)
  276.  
  277. *   CALCUL DE LA MATRICE DES COEFFICIENTS DE RESISTIVITE DANS LE
  278. *   PLAN,PAR RAPPORT AU REPERE  LOCAL DE L'ELEMENT
  279.             COS2 = COSA*COSA
  280.             SIN2 = SINA*SINA
  281.             SICO = SINA*COSA
  282.             YK(1,1) = COS2*XK1 + SIN2*XK2
  283.             YK(1,2) = SICO*(XK1-XK2)
  284.             YK(2,1) = YK(1,2)
  285.             YK(2,2) = SIN2*XK1 + COS2*XK2
  286.  
  287. *    ON AJOUTE LE PRODUIT DJAC*TRANSPOSEE(GRAD)*YK*GRAD
  288. *    POUR LE POINT DE GAUSS CONSIDERE,A LA MATRICE CEL1
  289.             CALL BDBST (GRAD,DJAC,YK,NBNN,NDIM,CEL1)
  290.  
  291.           ENDIF
  292.  
  293.         ENDDO
  294. *     FIN BOUCLE (20) SUR LES POINTS DE GAUSS
  295.  
  296. *     LE JACOBIEN EST NEGATIF ,MAILLAGE INCORRECT
  297.         IF (IFOIS.NE.0.AND.IFOIS.NE.NBPGAU) THEN
  298.           INTERR(1) = IEL
  299.           CALL ERREUR(195)
  300. *     CAS OU LE JACOBIEN EST TRES PETIT
  301.         ELSE IF (IFOI2.EQ.NBPGAU) THEN
  302.           INTERR(1) = IEL
  303.           CALL ERREUR (259)
  304.         ENDIF
  305.         IF (IERR.NE.0) GO TO 999
  306.  
  307. *     REMPLISSAGE DE XMATRI
  308.  
  309.         CALL REMPMT(CEL1,NBNN,RE(1,1,IEL))
  310.  
  311.       ENDDO
  312. *   FIN BOUCLE (10) SUR LES ELEMENTS
  313.  
  314. *  DESACTIVATION DES SEGMENTS
  315.  
  316.  999  CONTINUE
  317.       SEGSUP,MMAT1
  318.  990  CONTINUE
  319.       mptval = IVAMAT
  320.       IF (mptval.NE.0) SEGSUP,mptval
  321.  
  322.       nomid = MOMATR
  323.       SEGSUP,nomid
  324.       notype = MOTYR8
  325.       SEGSUP,notype
  326.  
  327. c      return
  328.       END
  329.  
  330.  

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