Télécharger la plaquette Cast3M


fr en

Télécharger tgrad.eso

Retour à la liste

Numérotation des lignes : 

tgrad

  1. C TGRAD     SOURCE    FANDEUR   22/01/03    21:15:52     11237          
  2.       SUBROUTINE TGRAD (IPTAB1)
  3. C
  4. C***********************************************************************
  5. C
  6. C FONCTION:
  7. C ---------
  8. C
  9. C     Calcul de grad en "0D/1D" sur des éléments de type POINT.
  10. C
  11. C
  12. C ENTREE :
  13. C --------
  14. C
  15. C   IPTAB1 : TABLE de soustype 'OPER_0D' contenant les indices suivants
  16. C            (pointeur, type ENTIER)
  17. C
  18. C     'GEOINF'     : TABLE de soustype 'GEOINF', info. géométriques
  19. C     'INCO'       : TABLE de soustype 'INCO', champs instanciés à
  20. C                    l'itéré précédent
  21. C     ('POTENTIA') : Champ potentiel (par exemple g*(z-zref))
  22. C                    (type CHPOINT de support PRIMAL(z) et DUAL(zref)).
  23. C     ('MULT1')    : Champ multiplicateur
  24. C                    (type FLOTTANT, CHPOINT de support DUAL ou MOT).
  25. C     'DUAL'       : Nom de l'inconnue duale
  26. C                    (type MOT, indice de la table 'INCO').
  27. C     'PRIMAL'     : Nom de l'inconnue primale
  28. C                    (type MOT, indice de la table 'INCO').
  29. C
  30. C
  31. C RESULTAT :
  32. C ----------
  33. C
  34. C     'LHS'        : Matrice élémentaire associée à l'opération
  35. C                    (type RIGIDITE).
  36. C
  37. C
  38. C AUTEUR, DATE DE CREATION:
  39. C -------------------------
  40. C
  41. C     1996/12 Laurent DADA     : Création
  42. C     2014/10 Frédéric DABBENE : Ajout du champ multiplicateur
  43. C
  44. C
  45. C LANGAGE:
  46. C --------
  47. C
  48. C     ESOPE + FORTRAN77
  49. C
  50. C***********************************************************************
  51. C
  52.       IMPLICIT INTEGER(I-N)
  53.       IMPLICIT REAL*8 (A-H,O-Z)
  54. C
  55.  
  56. -INC PPARAM
  57. -INC CCOPTIO
  58. -INC CCGEOME
  59. -INC SMCOORD
  60. -INC SMTABLE
  61.       POINTEUR MTABG.MTABLE,MTABS.MTABLE,IPTABI.MTABLE
  62. -INC SMCHPOI
  63. -INC SMELEME
  64.       POINTEUR IPTD1.MELEME,IPTP1.MELEME,IPJUCE.MELEME
  65. -INC SMRIGID
  66. C
  67.       SEGMENT REDI
  68.         INTEGER ORDR1(NNGOT)
  69.         INTEGER ORDR2(NNGOT)
  70.         INTEGER ORDR3(NNGOT)
  71.         INTEGER ORDR4(NNGOT)
  72.       ENDSEGMENT
  73. C
  74.       CHARACTER*8 TYPE,MOTI,MOT1,NOMPR1,NOMDU1,NOSUP1,NOSUD1
  75.       CHARACTER*7 NAMT1
  76.       LOGICAL     LPOT,LMULT,LCOEF
  77. C
  78. C Lecture de la table GEOINF de la table OPER_0D
  79. C
  80.       TYPE = 'TABLE   '
  81.       CALL ACMO (IPTAB1,'GEOINF',TYPE,IPTABG)
  82.       IF (IERR.NE.0) RETURN
  83. C
  84.       MOTI = 'SOUSTYPE'
  85.       CALL ACMM (IPTABG,MOTI,MOT1)
  86.       IF (IERR.NE.0) RETURN
  87.       IF (MOT1(1:6).NE.'GEOINF') THEN
  88.         MOTERR(1:8)  = 'GEOINF  '
  89.         MOTERR(9:16) = 'GEOINF  '
  90.         CALL ERREUR (790)
  91.         RETURN
  92.       ENDIF
  93. C
  94. C Lecture de la table INCO dans la table OPER_0D
  95. C
  96.       TYPE = 'TABLE   '
  97.       CALL ACMO (IPTAB1,'INCO',TYPE,IPTABI)
  98.       IF (IERR.NE.0) RETURN
  99. C
  100.       MOTI = 'SOUSTYPE'
  101.       CALL ACMM (IPTABI,MOTI,MOT1)
  102.       IF (IERR.NE.0) RETURN
  103.       IF (MOT1(1:4).NE.'INCO') THEN
  104.         MOTERR(1:8)  = 'INCO    '
  105.         MOTERR(9:16) = 'INCO    '
  106.         CALL ERREUR (790)
  107.         RETURN
  108.       ENDIF
  109. C
  110. C Lecture de la table SUPPORT dans la table INCO
  111. C
  112.       TYPE = 'TABLE   '
  113.       CALL ACMO (IPTABI,'SUPPORT',TYPE,IPTABS)
  114.       IF (IERR.NE.0) RETURN
  115. C
  116. C Lecture du MAILLAGE des connectivités 'JUNCEL' de la table GEOINF
  117. C Arrêt si les éléments ne sont pas des SEG3
  118. C
  119.       TYPE = 'MAILLAGE'
  120.       CALL ACMO (IPTABG,'JUNCEL',TYPE,IPJUCE)
  121.       IF (IERR.NE.0) RETURN
  122.       SEGACT IPJUCE
  123.       IF ((IPJUCE.ITYPEL).NE.3) THEN
  124.         MOTERR(1:8)  = 'JUNCEL  '
  125.         MOTERR(9:16) = 'MAILLAGE'
  126.         CALL ERREUR (787)
  127.         SEGDES IPJUCE
  128.         RETURN
  129.       ENDIF
  130.       NBEJC1 = IPJUCE.NUM(/2)
  131. C
  132. C Lecture du nom de l'inconnue PRIMAL
  133. C
  134.       TYPE = '        '
  135.       CALL ACMO (IPTAB1,'PRIMAL',TYPE,IPR1)
  136.       IF (IERR.NE.0) RETURN
  137.       IF (TYPE.EQ.'MOT     ') THEN
  138.         CALL ACMM (IPTAB1,'PRIMAL',NOMPR1)
  139.         IF (IERR.NE.0) RETURN
  140.       ELSE
  141.         MOTERR(1:8)  = 'PRIMAL  '
  142.         MOTERR(9:16) = TYPE
  143.         CALL ERREUR (787)
  144.         RETURN
  145.       ENDIF
  146. C
  147. C Lecture du nom de l'inconnue DUAL
  148. C
  149.       TYPE = '        '
  150.       CALL ACMO (IPTAB1,'DUAL',TYPE,IDU1)
  151.       IF (IERR.NE.0) RETURN
  152.       IF (TYPE.EQ.'MOT     ') THEN
  153.         CALL ACMM (IPTAB1,'DUAL',NOMDU1)
  154.         IF (IERR.NE.0) RETURN
  155.       ELSE
  156.         MOTERR(1:8)  = 'DUAL    '
  157.         MOTERR(9:16) = TYPE
  158.         CALL ERREUR (787)
  159.         RETURN
  160.       ENDIF
  161. C
  162. C Lecture éventuelle du champ POTENTIA
  163. C (CHPOINT à une composante)
  164. C
  165.       LPOT = .FALSE.
  166.       TYPE = '        '
  167.       CALL ACMO (IPTAB1,'POTENTIA',TYPE,IPCHP1)
  168.       IF (IERR.NE.0) RETURN
  169.       IF (TYPE.EQ.'CHPOINT ') THEN
  170.         LPOT = .TRUE.
  171.         MCHPOI = IPCHP1
  172.         SEGACT MCHPOI
  173.         MSOUPO = IPCHP(1)
  174.         SEGDES MCHPOI
  175.         SEGACT MSOUPO
  176.         NC = NOHARM(/1)
  177.         IF (NC.NE.1) THEN
  178.           MOTERR(1:8)  = 'POTENTIA'
  179.           MOTERR(9:16) = 'CHPOINT '
  180.           CALL ERREUR (784)
  181.           SEGDES MSOUPO
  182.           RETURN
  183.         ENDIF
  184.         IPT1   = IGEOC
  185.         MPOVA1 = IPOVAL
  186.         SEGDES MSOUPO
  187.       ENDIF
  188. C
  189. C Lecture éventuelle du champ MULT1 (CHPO ou FLOTTANT ou MOT)
  190. C (si MOT, c'est l'indice de la table INCO ou le CHPO est stocké)
  191. C
  192.       LMULT = .FALSE.
  193.       LCOEF = .FALSE.
  194.       VAL1  = 1.D0
  195.       TYPE = '        '
  196.       CALL ACMO (IPTAB1,'MULT1',TYPE,IPCHP1)
  197.       IF (IERR.NE.0) RETURN
  198.       IF ((TYPE.EQ.'CHPOINT ') .OR. (TYPE.EQ.'MOT     ')) THEN
  199.         IF (TYPE.EQ.'MOT     ') THEN
  200.           CALL ACMM (IPTAB1,'MULT1',MOT1)
  201.           IF (IERR.NE.0) RETURN
  202.           TYPE = 'CHPOINT '
  203.           CALL ACMO (IPTABI,MOT1,TYPE,IPCHP1)
  204.           IF (IERR.NE.0) RETURN
  205.         ENDIF
  206.         LMULT = .TRUE.
  207.         MCHPOI = IPCHP1
  208.         SEGACT MCHPOI
  209.         MSOUPO = IPCHP(1)
  210.         SEGDES MCHPOI
  211.         SEGACT MSOUPO
  212.         NC = NOHARM(/1)
  213.         IF (NC.NE.1) THEN
  214.           MOTERR(1:8)  = 'MULT    '
  215.           MOTERR(9:16) = 'CHPOINT '
  216.           CALL ERREUR (784)
  217.           SEGDES MSOUPO
  218.           RETURN
  219.         ENDIF
  220.         IPT2   = IGEOC
  221.         MPOVA2 = IPOVAL
  222.         SEGDES MSOUPO
  223.       ENDIF
  224.       IF (TYPE.EQ.'FLOTTANT') THEN
  225.         LCOEF = .TRUE.
  226.         CALL ACMF(IPTAB1,'MULT1',VAL1)
  227.         IF (IERR.NE.0) RETURN
  228.       ENDIF
  229. C
  230. C Lecture du nom du support de l'inconnue PRIMAL
  231. C Lecture du MAILLAGE de l'inconnue PRIMAL
  232. C
  233.       TYPE = '        '
  234.       CALL ACMO (IPTABS,NOMPR1,TYPE,ISUP1)
  235.       IF (IERR.NE.0) RETURN
  236.       IF (TYPE.EQ.'MOT     ') THEN
  237.         CALL ACMM (IPTABS,NOMPR1,NOSUP1)
  238.         IF (IERR.NE.0) RETURN
  239.       ENDIF
  240.       TYPE = 'MAILLAGE'
  241.       CALL ACMO (IPTABG,NOSUP1,TYPE,IPTP1)
  242.       IF (IERR.NE.0) RETURN
  243. C
  244. C Lecture du nom du support de l'inconnue DUAL
  245. C Arrêt si différent de 'JUNCTION'
  246. C Lecture du MAILLAGE de l'inconnue DUAL
  247. C
  248.       TYPE = '        '
  249.       CALL ACMO (IPTABS,NOMDU1,TYPE,ISUD1)
  250.       IF (TYPE.EQ.'MOT     ') THEN
  251.         CALL ACMM (IPTABS,NOMDU1,NOSUD1)
  252.         IF (IERR.NE.0) RETURN
  253.         IF (NOSUD1.NE.'JUNCTION') THEN
  254.           MOTERR(1:8)  = 'DUAL    '
  255.           MOTERR(9:16) = 'CHPOINT '
  256.           CALL ERREUR (788)
  257.           RETURN
  258.         ENDIF
  259.       ENDIF
  260.       TYPE = 'MAILLAGE'
  261.       CALL ACMO (IPTABG,NOSUD1,TYPE,IPTD1)
  262.       IF (IERR.NE.0) RETURN
  263. C
  264. C Création du support géométrique pour la RIGIDITE
  265. C (maillage de type SUPER-ELEMENT).
  266. C
  267.       NNGOT  = nbpts
  268.       SEGINI REDI
  269. C
  270. C On fusionne les maillages de POI1 des supports des inconnues PRIMAL
  271. C et DUAL en un maillage de type SUPER-ELEMENT.
  272. C
  273.       SEGACT IPTD1
  274.       IF (IPTD1.ITYPEL.NE.1)   CALL CHANGE (IPTD1,1)
  275.       SEGACT IPTD1
  276.       SEGACT IPTP1
  277.       IF (IPTP1.ITYPEL.NE.1)   CALL CHANGE (IPTP1,1)
  278.       SEGACT IPTP1
  279.       NBNNP1 = IPTP1.NUM(/2)
  280.       NBNND1 = IPTD1.NUM(/2)
  281.       NBNN   = NBNNP1 + NBNND1
  282.       NBSOUS = 0
  283.       NBREF  = 0
  284.       NBELEM = 1
  285.       SEGINI MELEME
  286.       ICOLOR(1) = IDCOUL
  287.       ITYPEL    = 28
  288.       DO 50 I50=1,NBNNP1
  289.         NUM(I50,1) = IPTP1.NUM(1,I50)
  290.         ORDR2(NUM(I50,1)) = I50
  291.  50   CONTINUE
  292.       DO 60 I60=1,NBNND1
  293.         NUM(I60+NBNNP1,1) = IPTD1.NUM(1,I60)
  294.         ORDR1(NUM(I60+NBNNP1,1)) = I60
  295.  60   CONTINUE
  296. C
  297. C On crée le tableau de redirection des CHPO POTENTIA et MULT
  298. C On vérifie qu'on a les informations pour tous les noeuds
  299. C
  300.       IF (LPOT) THEN
  301.          SEGACT IPT1
  302.          NBIPT1 = IPT1.NUM(/2)
  303.          DO 65 I65=1,NBIPT1
  304.             ORDR3(IPT1.NUM(1,I65)) = I65
  305.  65      CONTINUE
  306.          SEGDES IPT1
  307.          DO 655 I655=1,NBNN
  308.             IF (ORDR3(NUM(I655,1)).EQ.0) THEN
  309.                MOTERR(1:8)  = 'POTENTIA'
  310.                MOTERR(9:16) = 'CHPOINT '
  311.                CALL ERREUR (788)
  312.                SEGDES IPJUCE
  313.                SEGSUP MELEME
  314.                SEGDES IPTP1
  315.                SEGDES IPTD1
  316.                SEGSUP REDI
  317.                RETURN
  318.             ENDIF
  319.  655     CONTINUE
  320.       ENDIF
  321.       IF (LMULT) THEN
  322.          SEGACT IPT2
  323.          NBIPT2 = IPT2.NUM(/2)
  324.          DO 66 I66=1,NBIPT2
  325.             ORDR4(IPT2.NUM(1,I66)) = I66
  326.  66      CONTINUE
  327.          SEGDES IPT2
  328.          DO 675 I675=1,NBNND1
  329.             IF (ORDR4(IPTD1.NUM(I675,1)).EQ.0) THEN
  330.                MOTERR(1:8)  = 'MULT1'
  331.                MOTERR(9:16) = 'CHPOINT '
  332.                CALL ERREUR (788)
  333.                SEGDES IPJUCE
  334.                SEGSUP MELEME
  335.                SEGDES IPTP1
  336.                SEGDES IPTD1
  337.                SEGSUP REDI
  338.                RETURN
  339.             ENDIF
  340.  675     CONTINUE
  341.       ENDIF
  342.       SEGDES IPTP1
  343.       SEGDES IPTD1
  344. C
  345. C Création de la RIGIDITE
  346. C
  347.       NRIGE  = 7
  348.       NRIGEL = 1
  349.       SEGINI MRIGID
  350. C
  351.       MTYMAT = 'RIGIDITE'
  352.       IFORIG = IFOUR
  353.       ICHOLE = 0
  354.       IMGEO1 = 0
  355.       IMGEO2 = 0
  356.       ISUPEQ = 0
  357.       COERIG(1)   = 1.D0
  358.       IRIGEL(1,1) = MELEME
  359.       IRIGEL(2,1) = 0
  360.       IRIGEL(5,1) = NIFOUR
  361.       IRIGEL(6,1) = 0
  362.       IRIGEL(7,1) = 2
  363. C
  364.       SEGDES MELEME
  365. C
  366. C Remplissage du descripteur de l'objet RIGIDITE
  367. C
  368.       NLIGRP = NBNNP1
  369.       NLIGRD = NBNND1
  370.       SEGINI DESCR
  371.       IRIGEL(3,1) = DESCR
  372.       DO 10 I10=1,NBNNP1
  373.         NOELEP(I10) = I10
  374.         LISINC(I10) = NOMPR1
  375. 10    CONTINUE
  376.       DO 11 I11=1,NBNND1
  377.         NOELED(I11) = I11+NBNNP1
  378.         LISDUA(I11) = NOMDU1
  379. 11    CONTINUE
  380. C
  381.       SEGDES DESCR
  382. C
  383.       NELRIG = 1
  384.       SEGINI XMATRI
  385. C
  386.       IRIGEL(4,1) = XMATRI
  387.       xmatri.symre=2
  388. C
  389. C Calcul de la matrice élémentaire
  390. C
  391. c      DO 70 I70=1,NBNND1
  392. c        DO 701 I701=1,NBNNP1
  393. c          RE(I70,I701,1) = 0.D0
  394. c 701    CONTINUE
  395. c 70   CONTINUE
  396. C
  397.       IF (LMULT) THEN
  398.          SEGACT MPOVA2
  399.       ENDIF
  400. C
  401.       IF (LPOT) THEN
  402.         SEGACT MPOVA1
  403.         DO 72 I72=1,NBEJC1
  404.           NPT1  = IPJUCE.NUM(1,I72)
  405.           NPTF1 = IPJUCE.NUM(2,I72)
  406.           NPT2  = IPJUCE.NUM(3,I72)
  407.           XGZP1 = MPOVA1.VPOCHA(ordr3(npt1),1)
  408.           XGZP2 = MPOVA1.VPOCHA(ordr3(npt2),1)
  409.           XGZD1 = MPOVA1.VPOCHA(ordr3(nptf1),1)
  410.           IF (LMULT) THEN
  411.              VAL1 = MPOVA2.VPOCHA(ordr4(nptf1),1)
  412.           ENDIF
  413.           IPRI1 = ORDR2(NPT1)
  414.           IDUA1 = ORDR1(NPTF1)
  415.           IPRI2 = ORDR2(NPT2)
  416.           RE(IDUA1,IPRI1,1) = -1.D0 * (XGZP1 - XGZD1) * VAL1
  417.           RE(IDUA1,IPRI2,1) =         (XGZP2 - XGZD1) * VAL1
  418.  72     CONTINUE
  419.         SEGDES MPOVA1
  420.       ELSE
  421.         DO 82 I82=1,NBEJC1
  422.           NPT1  = IPJUCE.NUM(1,I82)
  423.           NPTF1 = IPJUCE.NUM(2,I82)
  424.           NPT2  = IPJUCE.NUM(3,I82)
  425.           IF (LMULT) THEN
  426.              VAL1 = MPOVA2.VPOCHA(ordr4(nptf1),1)
  427.           ENDIF
  428.           IPRI1 = ORDR2(NPT1)
  429.           IDUA1 = ORDR1(NPTF1)
  430.           IPRI2 = ORDR2(NPT2)
  431.           RE(IDUA1,IPRI1,1) = -1.D0 * VAL1
  432.           RE(IDUA1,IPRI2,1) =  1.D0 * VAL1
  433.  82     CONTINUE
  434.       ENDIF
  435. C
  436.       IF (LMULT) THEN
  437.          SEGDES MPOVA2
  438.       ENDIF
  439.       SEGDES IPJUCE
  440.       SEGDES XMATRI
  441.       SEGSUP REDI
  442. C
  443. C Ecriture du résultat
  444. C
  445.       TYPE = 'RIGIDITE'
  446.       CALL ECMO (IPTAB1,'LHS',TYPE,MRIGID)
  447.       IF (IERR.NE.0) RETURN
  448. C
  449.       SEGDES MRIGID
  450. C
  451.       END
  452.  
  453.  
  454.  
  455.  
  456.  
  457.  
  458.  
  459.  
  460.  
  461.  
  462.  

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