Télécharger la plaquette Cast3M


fr en

Télécharger frvisq.eso

Retour à la liste

Numérotation des lignes : 

frvisq

  1. C FRVISQ    SOURCE    OF166741  25/02/21    21:16:58     12166          
  2.  
  3.       SUBROUTINE FRVISQ(IPMODL,JPMAIL,IPCHE1, IPRIG)
  4. C
  5. C***********************************************************************
  6. C                                                                      *
  7. C  Routine principale appelée par AMOR                                 *
  8. C                                                                      *
  9. C  Calcule la matrice d'amortissement associée à la frontière du       *
  10. C  maillage dans plusieurs cas :                                       *
  11. C                                                                      *
  12. C   FORMULATION MECANIQUE                                              *
  13. C   +++++++++++++++++++++                                              *
  14. C                                                                      *
  15. C     * cas des massifs, dont l'enveloppe est constituée de SEG2 ou    *
  16. C       SEG3 (cas 2D), FAC3, FAC4, FAC6, ou FAC8 (cas 3D)              *
  17. C                                                                      *
  18. C   FORMULATION LIQUIDE                                                *
  19. C   +++++++++++++++++++                                                *
  20. C                                                                      *
  21. C     * cas des éléments dont l'enveloppe est constituée d'éléments    *
  22. C       à 2 (cas 2D), 3 ou 4 noeuds (cas 3D).                          *
  23. C______________________________________________________________________*
  24. C                                                                      *
  25. C  Entrées :                                                           *
  26. C  --------                                                            *
  27. C                                                                      *
  28. C     IPMODL : pointeur sur le modèle, objet MMODEL                    *
  29. C     JPMAIL : pointeur sur le maillage de la frontière, objet MELEME  *
  30. C     IPCHE1 : pointeur sur le champ par éléments de caractéristiques  *
  31. C              matériau, objet MCHAML                                  *
  32. C  Sorties :                                                           *
  33. C  --------                                                            *
  34. C     IPRIG : pointeur sur la matrice d'amortissement construite,      *
  35. C             objet MRIGID (=0 en cas d'erreur)                        *
  36. C                                                                      *
  37. C***********************************************************************
  38. C
  39.       IMPLICIT INTEGER(I-N)
  40.       IMPLICIT REAL*8(A-H,O-Z)
  41.  
  42. -INC PPARAM
  43. -INC CCOPTIO
  44. -INC CCHAMP
  45.  
  46. -INC SMCOORD
  47. -INC SMELEME
  48. -INC SMMODEL
  49. -INC SMCHAML
  50. -INC SMRIGID
  51. -INC SMINTE
  52.  
  53. -INC TMPTVAL
  54.  
  55.       INTEGER oooval
  56.  
  57.       SEGMENT INFO
  58.         INTEGER INFELL(JG)
  59.       ENDSEGMENT
  60. C
  61.       SEGMENT NOTYPE
  62.         CHARACTER*16 TYPE(NBTYPE)
  63.       ENDSEGMENT
  64. C
  65.       CHARACTER*(NCONCH) CONM
  66.  
  67. C Support du champ de caracteristiques
  68.       PARAMETER ( IPLAZ=3 )
  69.  
  70.       PARAMETER (NINF=3)
  71.       INTEGER INFOS(NINF)
  72.  
  73.       IPRIG = 0
  74.  
  75.       IF (IFOUR.EQ.-3) THEN
  76.         CALL ERREUR(21)
  77.         RETURN
  78.       ENDIF
  79.  
  80. c_______________________________________________________________________
  81. c
  82. c     activation du modele
  83. c_______________________________________________________________________
  84. C
  85.       MMODEL=IPMODL
  86.       SEGACT MMODEL
  87.       NSOUS=KMODEL(/1)
  88.  
  89. C______________________________________________________________________C
  90. C                                                                      C
  91. C  CREATION DE L'OBJET MATRICE DE RIGIDITE                             C
  92. C______________________________________________________________________C
  93. C                                                                      C
  94.       NRIGEL=0
  95.       SEGINI,MRIGID
  96.       MTYMAT='RIGIDITE'
  97.       IFORIG=IFOUR
  98.       ICHOLE=0
  99.       IMGEO1=0
  100.       IMGEO2=0
  101.       ISUPEQ=0
  102.       JRCOND=0
  103.       JRDEPP=0
  104.       JRDEPD=0
  105.  
  106. C______________________________________________________________________C
  107. C                                                                      C
  108. C  BOUCLE SUR LES SOUS ZONES                                           C
  109. C______________________________________________________________________C
  110. C                                                                      C
  111.       DO 100 ISOUS = 1, NSOUS
  112. C
  113. C  on récupère l'information générale
  114. C
  115.         IMODEL = KMODEL(ISOUS)
  116.         SEGACT,IMODEL
  117.  
  118. C- Initialisations
  119.         IENVEL = 0
  120.         IPOGEO = 0
  121.  
  122.         IPT1 = IMAMOD
  123.         CONM = CONMOD
  124.         MELM = NEFMOD
  125. C
  126. C  création du tableau info
  127. C
  128.         iret = 1
  129.         CALL IDENT(IPT1,CONM,IPCHE1,0,INFOS,iret)
  130.         IF (iret.EQ.0) GOTO 1099
  131. C
  132. C  Determination de l'enveloppe du maillage massif du sous-modele
  133. C
  134.         CALL ECROBJ('MAILLAGE',IPT1)
  135.         IF (IDIM.EQ.3) THEN
  136.           CALL ENVELO
  137.         ELSE IF (IDIM.EQ.2) THEN
  138.           CALL PRCONT
  139.         ELSE
  140.           CALL ERREUR(5)
  141.         ENDIF
  142.         IF (IERR.NE.0) GOTO 1099
  143.         CALL LIROBJ('MAILLAGE',IENVEL,1,iret)
  144.         IF (IERR.NE.0) GOTO 1099
  145. C
  146. C  Elements de l'enveloppe IENVEL dans le maillage frontiere JPMAIL
  147. C
  148.         iret = 0
  149.         CALL INTERB(IENVEL,JPMAIL,iret,IPOGEO)
  150.         IF (iret.GT.0) GOTO 1099
  151.  
  152.         IPT3 = IPOGEO
  153.         SEGACT,IPT3
  154.         NBSOU3 = IPT3.LISOUS(/1)
  155.         IPT2 = IPT3
  156. C
  157. C  boucle sur les sous-zones de l'enveloppe
  158. C
  159.         DO 110 IB = 1, MAX(1,NBSOU3)
  160.  
  161. C-- Initialisations :
  162.           MOFORC = 0
  163.           MODEPL = 0
  164.           IPMINT = 0
  165.           MOMATR = 0
  166.           MOCARA = 0
  167.           MOTYPM = 0
  168.           MOTYPC = 0
  169.           ISUPM  = 0
  170.           ISUPC  = 0
  171.           IDESCR = 0
  172.  
  173. C-- Informations sur la (sous-zone de) l'enveloppe
  174.           IF (NBSOU3.NE.0) THEN
  175.             IPT2 = IPT3.LISOUS(IB)
  176.             SEGACT,IPT2
  177.           ENDIF
  178.           NBNOE2 = IPT2.NUM(/1)
  179.           NBELE2 = IPT2.NUM(/2)
  180.           LETYP  = IPT2.ITYPEL
  181. C-- Petit test sur le type
  182.           IF (LETYP.EQ.1) THEN
  183.             CALL ERREUR(16)
  184.             GOTO 1199
  185.           ENDIF
  186.           IPOGEO = IPT2
  187. C
  188. C-- On détermine la formulation associée à l'objet géométrique
  189. C-- elementaire de surface
  190.           CALL TYPFAC(MELM,NBNOE2,MELE)
  191. C
  192. C-- ERREUR : impossible d'utiliser FROABS pour les éléments
  193. C--          de formulation MELM
  194.           IF (MELE.EQ.0) THEN
  195.             MOTERR(1:8) = NOMTP(MELM)
  196.             CALL ERREUR(193)
  197.             GOTO 1199
  198.           ENDIF
  199.  
  200. C-- Information sur l'élément fini
  201.           CALL ELQUOI(MELE,0,3,IPINF,IMODEL)
  202.           IF (IERR.NE.0) GOTO 1199
  203. C
  204.           INFO  = IPINF
  205.           MFR   = INFELL(13)
  206.           LRE   = INFELL(9)
  207.           LW    = INFELL(7)
  208.           NDDL  = INFELL(15)
  209. c*          IELE  = INFELL(14)
  210.           IPPORE = 0
  211.           IPMINT = INFELL(11)
  212.           SEGSUP,INFO
  213.  
  214. C-- Recherche des inconnues primales et duales (DEPL-FORC)
  215.           CALL IDPRIM(IMODEL,MFR,MODEPL,NDEPL,ndum)
  216.           CALL IDDUAL(IMODEL,MFR,MOFORC,NFORC,ndum)
  217.  
  218.           IF (NDEPL.EQ.0 .OR. NFORC.EQ.0 .OR. NDEPL.NE.NFORC) THEN
  219.             CALL ERREUR(5)
  220.             GOTO 1199
  221.           ENDIF
  222.  
  223. C-- Remplissage du segment DESCRipteur
  224.           NLIGRP = LRE
  225.           NLIGRD = LRE
  226.           SEGINI,DESCR
  227.  
  228.           NCOMP = NDEPL
  229.           NBNNS = NBNOE2
  230.           IDDL=1
  231.           DO 1004 INOEUD=1,NBNNS
  232.             DO 1005 ICOMP=1,NCOMP
  233.               NOMID = MODEPL
  234.               LISINC(IDDL)=LESOBL(ICOMP)
  235.               NOMID = MOFORC
  236.               LISDUA(IDDL)=LESOBL(ICOMP)
  237.               NOELEP(IDDL)=INOEUD
  238.               NOELED(IDDL)=INOEUD
  239.               IDDL=IDDL+1
  240.  1005       CONTINUE
  241.  1004     CONTINUE
  242.  
  243.           IDESCR = DESCR
  244.  
  245. C-- Recuperation des noms de composantes MATERIAU
  246.           nbrobl = 0
  247.           nbrfac = 0
  248.           nomid  = 0
  249.           notype = 0
  250.  
  251. C  rho, E, nu pour les massifs
  252.           IF (MFR.EQ.1) THEN
  253.             nbrobl = 3
  254.             SEGINI,nomid
  255.             lesobl(1) = 'RHO '
  256.             lesobl(2) = 'YOUN'
  257.             lesobl(3) = 'NU  '
  258.  
  259.             nbtype = 1
  260.             SEGINI,notype
  261.             type(1) = 'REAL*8'
  262. C
  263. C  rho, cson, rhoref, cref, rlcar pour les liquides
  264.           ELSE IF (MFR.EQ.11.OR.MFR.EQ.41) THEN
  265.             nbrobl = 5
  266.             SEGINI,nomid
  267.             lesobl(1) = 'RHO '
  268.             lesobl(2) = 'CSON'
  269.             lesobl(3) = 'RORF'
  270.             lesobl(4) = 'CREF'
  271.             lesobl(5) = 'LCAR'
  272.  
  273.             nbtype = 1
  274.             SEGINI,notype
  275.             type(1) = 'REAL*8'
  276.           ENDIF
  277.  
  278.           MOMATR = nomid
  279.           MOTYPM = notype
  280.           NMATR = nbrobl
  281.           NMATF = nbrfac
  282.           NMATT = NMATR+NMATF
  283.  
  284. C--- Verification du support des composantes recherchées
  285.           IF (MOMATR.NE.0) THEN
  286.             CALL QUESUQ(IMODEL,IPCHE1,3,0,MOMATR,IPLAZ,ISUPM,iret)
  287.             IF (ISUPM.GT.1) GOTO 1199
  288.           ENDIF
  289. C
  290. C-- Recuperation des noms de composantes CARACTERISTIQUES
  291.           nbrobl = 0
  292.           nbrfac = 0
  293.           nomid  = 0
  294.           notype = 0
  295.  
  296. C Epaisseur du massif en contraintes planes
  297.           IF (MFR.EQ.1 .AND. IFOUR.EQ.-2) THEN
  298.             nbrfac = 1
  299.             SEGINI,nomid
  300.             lesfac(1) = 'DIM3'
  301.             nbtype = 1
  302.             SEGINI,notype
  303.             type(1) = 'REAL*8'
  304.           ENDIF
  305.  
  306.           MOCARA = nomid
  307.           MOTYPC = notype
  308.           NCARA = nbrobl
  309.           NCARF = nbrfac
  310.           NCARR = NCARA+NCARF
  311.  
  312. C--- Verification du support des composantes recherchées
  313.           IF (MOCARA.NE.0) THEN
  314.             CALL QUESUQ(IMODEL,IPCHE1,3,0,MOCARA,IPLAZ,ISUPC,iret)
  315.             IF (ISUPC.GT.1) GOTO 1199
  316.           ENDIF
  317.  
  318. C-- Segment d'integration MINTE
  319.           MINTE = IPMINT
  320.           SEGACT,MINTE
  321.           NBPGAU = POIGAU(/1)
  322.  
  323. C- Partionnement si necessaire de la matrice d'amortissement
  324. C- determinant ainsi le nombre d'objets elementaires de MRIGID
  325.           LTRK = oooval(1,4)
  326.           IF (LTRK.EQ.0) LTRK = oooval(1,1)
  327.       LTRK=MAX(LTRK,2**24)
  328. * Ajout a la taille en mots de la matrice des infos du segment
  329.           LSEG = LRE*LRE*NBELE2 + 16
  330.           NBLPRT = (LSEG-1)/LTRK + 1
  331.           NBLMAX = (NBELE2-1)/NBLPRT + 1
  332.           NBLPRT = (NBELE2-1)/NBLMAX + 1
  333. c*        write(ioimp,*) ' frvisq : nblprt nblmax = ',nblprt,nblmax,nbele2
  334. C*OF : Pour l'instant pas de partition pour FRVISQ
  335.           NBLPRT = 1
  336.  
  337. C-- Ajout de la matrice d'AMORTISSEMENT a la matrice globale
  338.           NRIGE0 = IRIGEL(/2)
  339.           NRIGEL = NRIGE0 + NBLPRT
  340.           SEGADJ,MRIGID
  341.  
  342.           descr  = IDESCR
  343.           meleme = IPOGEO
  344.           nbnn   = NBNOE2
  345.           nbelem = NBELE2
  346.           nbsous = 0
  347.           nbref  = 0
  348.  
  349.           DO 120 irige = 1, NBLPRT
  350.  
  351. C-- Mettre ici la partition du maillage IPOGEO
  352.             ipmail = meleme
  353.             ipdesc = descr
  354.  
  355. C- Initialisation de la matrice de rigidite elementaire (xmatri)
  356.             NELRIG = nbelem
  357.             SEGINI,xmatri
  358.             ipmatr = xmatri
  359.  
  360. C- Recuperation des valeurs des proprietes materiau et geometriques
  361. c* Note : les proprietes sont les valeurs au support des EF massifs
  362. c*        et non celles au niveau de l'enveloppe surfacique !
  363. c*        Cela ne marche que si les proprietes sont constantes. Dans
  364. c*        les autres cas, le resultat est... Pour eviter cela, on met
  365. c*        un test sur la constance du champ !
  366.             ivamat = 0
  367.             ivacar = 0
  368.             IF (MOMATR.NE.0) THEN
  369.               CALL KOMCHA(IPCHE1,IPT1,CONM,MOMATR,MOTYPM,1,
  370. c*              CALL KOMCHA(IPCHE1,ipmail,CONM,MOMATR,MOTYPM,1,
  371.      &                    INFOS,3, ivamat)
  372.               IF (IERR.NE.0) GOTO 1199
  373.               mptval = ivamat
  374.               do i = 1, NMATT
  375.                 if (ival(i).ne.0) then
  376.                   melval = IVAL(i)
  377.                   if (velche(/1).ne.1 .and. velche(/2).ne.1) then
  378.                     write(ioimp,*) 'Champ MATERIAU non constant'
  379.                     call erreur(21)
  380.                     goto 1199
  381.                   endif
  382.                 endif
  383.               enddo
  384.               IF (ISUPM.EQ.1) THEN
  385.                 CALL VALCHE(ivamat,NMATT,IPMINT,IPPORE,MOMATR,MELE)
  386.                 IF (IERR.NE.0) THEN
  387.                   ISUPM = 0
  388.                   GOTO 1199
  389.                 ENDIF
  390.               ENDIF
  391.             ENDIF
  392.             IF (MOCARA.NE.0) THEN
  393.               CALL KOMCHA(IPCHE1,IPT1,CONM,MOCARA,MOTYPC,1,
  394. c*              CALL KOMCHA(IPCHE1,ipmail,CONM,MOCARA,MOTYPC,1,
  395.      &                    INFOS,3, ivacar)
  396.               IF (IERR.NE.0) GOTO 1199
  397.               mptval = ivacar
  398.               do i = 1, NCARR
  399.                 if (ival(i).ne.0) then
  400.                   melval = IVAL(i)
  401.                   if (velche(/1).ne.1 .and. velche(/2).ne.1) then
  402.                     write(ioimp,*) 'Champ MATERIAU non constant'
  403.                     call erreur(21)
  404.                     goto 1199
  405.                   endif
  406.                 endif
  407.               enddo
  408.               IF (ISUPC.EQ.1) THEN
  409.                 CALL VALCHE(ivacar,NCARR,IPMINT,IPPORE,MOCARA,MELE)
  410.                 IF (IERR.NE.0)THEN
  411.                   ISUPC = 0
  412.                   GOTO 1199
  413.                 ENDIF
  414.               ENDIF
  415.             ENDIF
  416.  
  417. C  distinction des cas 2D et 3D
  418. C______________________________________________________________________C
  419. C                                                                      C
  420. C            CAS DES ELEMENTS MASSIFS BIDIMENSIONNELS                  C
  421. C                  FACES ASSOCIEES SEG2 OU SEG3                        C
  422. C______________________________________________________________________C
  423. C                                                                      C
  424.          IF (MELE.EQ.2.OR.MELE.EQ.3) THEN
  425. C
  426.             CALL FROA2D(ipmail,ipmatr,IPMINT,ivamat,ivacar,
  427.      1                  MELE,MFR,LRE,NDDL)
  428. C
  429. C______________________________________________________________________C
  430. C                                                                      C
  431. C               CAS DES ELEMENTS LIQUIDES 2D OU 3D                     C
  432. C               FACES ASSOCIEES LSE2, LTR3 OU LQU4                     C
  433. C______________________________________________________________________C
  434. C                                                                      C
  435. C
  436.          ELSE IF(MELE.EQ.97.OR.MELE.EQ.35.OR.MELE.EQ.36) THEN
  437. C
  438.             CALL LFROA(ipmail,ipmatr,IPMINT,ivamat,ivacar,
  439.      1                  MELE,MFR,LRE,NDDL)
  440. C
  441. C______________________________________________________________________C
  442. C                                                                      C
  443. C            CAS DES ELEMENTS MASSIFS TRIDIMENSIONNELS                 C
  444. C             FACES ASSOCIEES FAC3,FAC4,FAC6 OU FAC8                   C
  445. C______________________________________________________________________C
  446. C
  447.          ELSE IF(MELE.EQ.31.OR.MELE.EQ.32.OR.MELE.EQ.33.
  448.      1                                   OR.MELE.EQ.34)THEN
  449. C
  450.             CALL FROA3D(ipmail,ipmatr,IPMINT,ivamat,ivacar,
  451.      1                  MELE,MFR,LRE,NDDL)
  452. C
  453. C  erreur, l'élément n'est pas encore implémenté
  454. C
  455.          ELSE
  456. C
  457.             MOTERR(1:4)=NOMTP(MELE)
  458.             MOTERR(5:12)='FRVISQ'
  459.             CALL ERREUR (86)
  460.          ENDIF
  461. C
  462.          IF (ISUPM.EQ.1 .OR. nblprt.GT.1) THEN
  463.            CALL DTMVAL(ivamat,3)
  464.          ELSE
  465.            CALL DTMVAL(ivamat,1)
  466.          ENDIF
  467.          IF (ISUPC.EQ.1 .OR. nblprt.GT.1) THEN
  468.            CALL DTMVAL(ivacar,3)
  469.          ELSE
  470.            CALL DTMVAL(ivacar,1)
  471.          ENDIF
  472.  
  473.           xmatri = ipmatr
  474.           IF (NBLPRT.GT.1) THEN
  475.             meleme = ipmail
  476.           ENDIF
  477.  
  478. C- Sortie prematuree en cas d'erreur
  479.           IF (IERR.NE.0) GOTO 1199
  480.  
  481. C- Stockage de la matrice
  482.           jrige = NRIGE0 + irige
  483.           COERIG(jrige)   = 1.
  484.           IRIGEL(1,jrige) = ipmail
  485.           IRIGEL(2,jrige) = 0
  486.           IRIGEL(3,jrige) = ipdesc
  487.           IRIGEL(4,jrige) = ipmatr
  488.           IRIGEL(5,jrige) = NIFOUR
  489.           IRIGEL(6,jrige) = 0
  490.           IRIGEL(7,jrige) = 0
  491. * matrice non symetrique (forces sur pi seulement
  492. * qui dependent de p)
  493.           IF (MFR.EQ.11.OR.MFR.EQ.41) THEN
  494.             IRIGEL(7,jrige) = 2
  495.             xmatri.symre=2
  496.            ENDIF
  497.           SEGDES,xmatri
  498.           IRIGEL(8,jrige) = 0
  499.  
  500.  120      CONTINUE
  501. C- Fin de la boucle de partition maillage/rigidite
  502.  
  503.  1199     CONTINUE
  504.           IF (MOMATR.NE.0) THEN
  505.             nomid = MOMATR
  506.             SEGSUP,nomid
  507.             notype = MOTYPM
  508.             SEGSUP,notype
  509.           ENDIF
  510.           IF (MOCARA.NE.0) THEN
  511.             nomid = MOCARA
  512.             SEGSUP,nomid
  513.             notype = MOTYPC
  514.             SEGSUP,notype
  515.           ENDIF
  516. C
  517.  
  518. C- Sortie prematuree en cas d'erreur
  519.           IF (IERR.NE.0) GOTO 1098
  520.  
  521.  110    CONTINUE
  522. C- Fin de la boucle sur (les sous-zones de) l'enveloppe
  523. C
  524.  1098   CONTINUE
  525.  1099   CONTINUE
  526. C- Sortie prematuree en cas d'erreur
  527.         IF (IERR.NE.0) GOTO 999
  528. C
  529.  100  CONTINUE
  530. C- Fin de la boucle sur les modeles elementaires
  531. C
  532.       NRIGE0 = IRIGEL(/2)
  533.       IF (NRIGE0.EQ.0) THEN
  534.         CALL ERREUR(902)
  535.       ENDIF
  536.  
  537.  999  CONTINUE
  538.       IF (IERR.EQ.0) THEN
  539.         IPRIG = MRIGID
  540.         SEGDES,MRIGID
  541.       ELSE
  542.         IPRIG = 0
  543.         SEGSUP,MRIGID
  544.       ENDIF
  545.  
  546. c      RETURN
  547.       END
  548.  
  549.  
  550.  

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