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

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