Télécharger fpmass.eso

Retour à la liste

fpmass
C FPMASS    SOURCE    OF166741  25/02/06    21:15:06     12146          
 
C_____________________________________________________________________
C
C   CALCULE LES FORCES DE PRESSIONS APPLIQUEES SUR DES MASSIFS
C
C   ENTREES :
C   ---------
C
C    IPCHE1   CHPOINT CONTENANT LES VALEURS DES PRESSIONS AUX NOEUDS
C              DE LA FACE D UN MASSIF
C    IPCHM1   CHAMELEM CONTENANT LES VALEURS DES PRESSIONS AUX NOEUDS
C              DE LA FACE D UN MASSIF
C    IPMODL   OBJET MODELE  SUR LEQUEL S APPLIQUE LA PRESSION
C
C    JPMAIL   POINTEUR SUR LE MAILLAGE SI ON A LU UN FLOTTANT ET
C             UN MAILLAGE, SINON 0
C
C    XP       LA VALEUR DE LA PRESSION SI ON L'A LUE
C
C  SORTIES  :
C  ----------
C
C    IPTFP  = CHPOINT DES FORCES NODALES EQUIVALENTES
C    IRET   = 1 OU 0 SUIVANT SUCCES OU NON
C
C     REVISION JACQUELINE BROCHARD  SEPTEMBRE 86
C     MISE A JOUR P VERPEAUX MAI 88
C
C     PASSAGE AUX NOUVEAU CHAMELEM PAR JM CAMPENON LE 17 09 90
C_______________________________________________________________________
 
      SUBROUTINE FPMASS(IPCHE1,IPCHM1,IPMODL,IPTFP,JPMAIL,XP,IRET)
 
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
 
-INC PPARAM
-INC CCOPTIO
-INC CCREEL
-INC CCHAMP
 
-INC SMCOORD
-INC SMELEME
-INC SMMODEL
-INC SMCHAML
-INC SMCHPOI
-INC SMINTE
 
      SEGMENT INFO
        INTEGER INFELL(JG)
      ENDSEGMENT
 
      SEGMENT NOTYPE
        CHARACTER*16 TYPE(NBTYPE)
      ENDSEGMENT
 
      SEGMENT MPTVAL
        INTEGER IPOS(NS)    ,NSOF(NS)
        INTEGER      IVAL(NCOSOU)
        CHARACTER*16 TYVAL(NCOSOU)
      ENDSEGMENT
 
      segment netn(nonetn)
      segment ietn(letn)
 
      CHARACTER*4 MOSTRI,MOAPPU,MOGEOM
      CHARACTER*(NCONCH) CONM
      PARAMETER (NINF=3)
      INTEGER INFOS(NINF)
      LOGICAL LSUPFO,ltelq
 
      DATA MOAPPU/'APPU'/,MOSTRI/'STRI'/
      DATA MOGEOM/'GEOM'/
 
      IRET = 0
      IGEOM = 0
 
      NHRM=NIFOUR
 
C-----------------------------------------------------------------------
C LECTURE DU CHAMP DE CARACTERISTIQUES
C-----------------------------------------------------------------------
      IPCHE2 = 0
      ISUPCA = 0
C Prevoir la lecture en amont !
      CALL LIROBJ('MCHAML  ',IPCHE2,0,irt2)
      IF (IERR.NE.0) RETURN
      IF (IPCHE2.NE.0) THEN
        CALL ACTOBJ('MCHAML  ',IPCHE2,1)
        IF (IERR.NE.0) RETURN
 
        CALL REDUAF(IPCHE2,MODORI,ipche20,0,iretca,kerr)
        if (iretca.ne.1) call erreur(kerr)
        IF (IERR.NE.0) RETURN
        IPCHE2 = ipche20
C
C   Verification du lieu support du MCHAML de caracteristiques
C
        CALL QUESUP(IPMODL,IPCHE2,3,1,ISUPCA,iretca)
        IF (ISUPCA.GT.1) RETURN
      ENDIF
C-----------------------------------------------------------------------
C   CAS OU UN CHPOINT EST FOURNI
C   ON CREE L OBJET GEOMETRIQUE CONTENANT TOUS LES PTS SI BESOIN
      IF (JPMAIL.EQ.0.AND.IPCHM1.EQ.0) THEN
        MCHPOI=IPCHE1
        ltelq=.false.
        DO I=1,IPCHP(/1)
           MSOUPO=IPCHP(I)
           IF (I.GT.1) THEN
              CALL FUSE(IGEOM,IGEOC,IPPT,ltelq)
              IGEOM=IPPT
           ELSE
              IGEOM=IGEOC
           ENDIF
        ENDDO
        IF (IERR.NE.0) RETURN
      ENDIF
C-----------------------------------------------------------------------
C   CAS OU UN CHAMELEM EST FOURNI
C   ON CREE L OBJET GEOMETRIQUE CONTENANT TOUS LES PTS SI BESOIN
      IF (IPCHM1.NE.0) THEN
        MCHEL2 = IPCHM1
        ltelq=.false.
        DO I=1,MCHEL2.IMACHE(/1)
          IMTMP=MCHEL2.IMACHE(I)
          IF (I.GT.1) THEN
            CALL FUSE(IGEOM,IMTMP,IPPT,ltelq)
            IGEOM=IPPT
          ELSE
            IGEOM=IMTMP
          ENDIF
        ENDDO
        IF (IERR.NE.0) RETURN
      ENDIF
 
C= Cas des modes de calculs en DEFORMATIONS GENERALISEES
      IF (IFOUR.EQ.-3) THEN
        NDPGE=3
      ELSE IF (IFOUR.EQ.11) THEN
        NDPGE=2
      ELSE IF (IFOUR.EQ. 7.OR.IFOUR.EQ. 8.OR.IFOUR.EQ. 9.OR.
     &         IFOUR.EQ.10.OR.IFOUR.EQ.14) THEN
        NDPGE=1
      ELSE
        NDPGE=0
      ENDIF
 
C- Un petit segment toujours utile :
      nbtype = 1
      SEGINI,notype
      notype.TYPE(1) = 'REAL*8  '
      MOTYR8 = notype
 
c- Segments utiles pour accelerer la recherche des elements touchant un noeud.
      SEGACT,mcoord
      nonetn = nbpts+1
      netn = 0
      ietn = 0
C
C     TRAITEMENT DU MODELE
C
      MMODEL = IPMODL
      NSOUS  = mmodel.KMODEL(/1)
 
      IRRT=0
      DO 100 ISOUS = 1, NSOUS
C
        ISOK = 0
        MOCARA = 0
        IVACAR = 0
C
C        TRAITEMENT DU MODELE
C
        IMODEL = mmodel.KMODEL(ISOUS)
        IPMAIL = imodel.IMAMOD
        CONM   = imodel.CONMOD
        MELM   = imodel.NEFMOD
 
C*         write(*,*) ISOUS,'/',NSOUS,' : ',IMODEL,'NEFMOD=',MELM
         if ((melm .eq. 22).OR.(melm .eq. 259)) then
C      ... Ici sous modele de multiplicateur de lagrange on
C          incrémente le compteur et on passe à la zone suivante ...
            IRRT=IRRT+1
            GOTO 100
         endif
C
C        ON RECUPERE LES ELTS DE L ENVELOPPE DU MASSIF APPUYES
C        STRICTEMENT SUR LE CHPOINT DE PRESSIONS OU appartenant au
C        MAILLAGE DONNE
C
         CALL ECROBJ('MAILLAGE',IPMAIL)
         IF (IDIM.EQ.2) THEN
           CALL PRCONT
         ELSE IF (IDIM.EQ.3) THEN
           CALL ENVELO
         ELSE IF (IDIM.EQ.1) THEN
           CALL PREX1D
         ENDIF
         IF (IERR.NE.0) GOTO 9900
         CALL LIROBJ('MAILLAGE',ienvel,1,iretou)
         IF (IERR .NE. 0) GOTO 9900
         CALL ACTOBJ('MAILLAGE',ienvel,1)
 
C     ... si un CHPOINT a été donné, on va chercher la partie de
C         l'enveloppe s'appuyant strictement sur le support du CHPOINT ...
         IF (JPMAIL.EQ.0) THEN
              CALL ECROBJ('MAILLAGE',IGEOM)
              CALL ECRCHA(MOSTRI)
              CALL ECRCHA(MOAPPU)
              CALL ECROBJ('MAILLAGE',ienvel)
              CALL EXTREL(IRR,0,IBNOR)
         ELSE
C      ... sinon, on va chercher l'intersection de l'enveloppe avec
C          le maillage fourni ...
              CALL INTERB(ienvel,jpmail,irr,IPOGEO)
         ENDIF
C      ... Ici on teste si intersection est vide, si OUI on
C          incrémente le compteur et on passe à la zone suivante ...
         IF (irr.gt.0) THEN
              IRRT=IRRT+1
              GOTO 100
         ENDIF
         IF (JPMAIL.EQ.0) THEN
             CALL LIROBJ('MAILLAGE',IPOGEO,1,iret)
             IF (IERR.NE.0) GOTO 9900
             CALL ACTOBJ('MAILLAGE',IPOGEO,1)
         ENDIF
 
C  pour accelerer la recherche, utilisation d'un tableau des elements touchant un noeud.
        if (netn.EQ.0) THEN
          segini,netn
        else
          do i = 1, nonetn
            netn(i) = 0
          enddo
        endif
        IPT1 = IPMAIL
        nbnn1 = ipt1.num(/1)
        nbel1 = ipt1.num(/2)
        do j = 1, nbel1
          do i = 1, nbnn1
            ino = ipt1.num(i,j)
            netn(ino) = netn(ino)+1
          enddo
        enddo
        do i = 2, nonetn
          netn(i) = netn(i) + netn(i-1)
        enddo
        letn = netn(nonetn)
        if (ietn.eq.0) then
          segini,ietn
        else
          if (letn.gt.ietn(/1)) segadj,ietn
          do i = 1, letn
            ietn(i) = 0
          enddo
        endif
        do j = 1, nbel1
          do i = 1, nbnn1
            ino = ipt1.num(i,j)
            ietn(netn(ino)) = j
            netn(ino) = netn(ino)-1
          enddo
        enddo
        ietn1 = ietn
        netn1 = netn
 
C_______________________________________________________________________
C
C TRAITEMENT DES CHAMPS DE CARACTERISTIQUES
C_______________________________________________________________________
        NBROBL = 0
        NBRFAC = 0
        IF (IPCHE2.NE.0 .AND. IFOUR.EQ.-2) THEN
C
C       CREATION DU TABLEAU INFOS
C
          CALL IDENT(IPMAIL,CONM,IPCHE2,0,INFOS,IRTD)
          IF (IRTD.EQ.0) GOTO 9900
C
          NBRFAC=1
          SEGINI,NOMID
          LESFAC(1)='DIM3'
          MOCARA = NOMID
C
          IF (ISUPCA.NE.1) THEN
            CALL KOMCHA(IPCHE2,IPMAIL,CONM,MOCARA,MOTYR8,0,
     &                                     INFOS,3,IVACAR)
            IF (IERR.NE.0) GOTO 9900
          ENDIF
        ENDIF
        NCARA  = NBROBL
        NCARF  = NBRFAC
        NCARR  = NCARA+NCARF
C
C     ON DETERMINE LA FORMULATION ASSOCIEE A L OBJET
C     GEOMETRIQUE ELEMENTAIRE DE SURFACE
C
        IPT3 = IPOGEO
        NBSOU3 = IPT3.LISOUS(/1)
        IPT2=IPT3
C
C     BOUCLE SUR LES SOUS ZONES DE L ENVELOPPE
C
        DO 110 IB=1,MAX(1,NBSOU3)
 
          MOFORC = 0
          IVAFOR = 0
          IVACA1 = 0
          IPMOD1 = 0
          IPTVPR = 0
          lsupfo = .false.
          ISOK = 0
 
          IF (NBSOU3.NE.0) THEN
            IPT2=IPT3.LISOUS(IB)
          ENDIF
          IPOGEO=IPT2
 
          NBNN  = IPT2.NUM(/1)
          LETYP = IPT2.ITYPEL
C
C     PETIT TEST SUR LE TYPE
          IF (LETYP.EQ.1.AND.IDIM.NE.1) THEN
            CALL ERREUR(16)
            GOTO 9990
          ENDIF
          CALL TYPFAC(MELM,NBNN,MELE)
C          write(*,*) 'TYPFAC --> MELE=',MELE
C  ERREUR : IMPOSSIBLE D UTILISER L OPERATEUR PRESSI POUR
C           LES ELEMENTS DE FORMULATION MELM
          IF (MELE.EQ.0) THEN
            MOTERR(1:8)=NOMTP(MELM)
            CALL ERREUR(193)
            GOTO 9990
          ENDIF
C
C     CAS OU UN CHAMP PAR POINT A ETE FOURNI
C     ON CREE L OBJET MODEL ASSOCIE A LA SURFACE ELEMENTAIRE
C     ON TRANSFORME LE CHPOINT DE PRESSION EN CHELEM ELEMENTAIRE
          IF (JPMAIL.EQ.0.AND.IPCHM1.EQ.0) THEN
            N1 = 1
            SEGINI,MMODE1
            NFOR = imodel.FORMOD(/2)
            NMAT = imodel.MATMOD(/2)
c*            MN3  = imodel.INFMOD(/1)
            MN3    = 1
            NPARMO = 0
            NOBMOD = 0
            SEGINI,IMODE1
            imode1.IMAMOD = IPOGEO
            imode1.NEFMOD = MELE
            imode1.CONMOD = imodel.CONMOD
            DO i = 1, NFOR
              imode1.FORMOD(i) = imodel.FORMOD(i)
            ENDDO
            DO i = 1, NMAT
              imode1.MATMOD(i) = imodel.MATMOD(i)
            ENDDO
c*            DO i = 1, MN3
c*              imode1.INFMOD(i) = imodel.INFMOD(i)
c*            ENDDO
c*            lzero = 0
c*            call inomid(imode1,lzero,lzero,lzero,lzero)
c*            call prquoi(imode1)
            mmode1.KMODEL(1) = IMODE1
            IPMOD1 = MMODE1
c* Il faut redefinir a chaque fois IPMOD1 pour eviter rappel du 
c* preconditionnement dans CHAME1 qui ne cree pas 
            CALL CHAME1(0,IPMOD1,IPCHE1,' ',ICHELP,3)
            IF (IERR.NE.0) GOTO 9990
            MCHEL1=ICHELP
            MCHAM1=MCHEL1.ICHAML(1)
            IPTVPR=MCHAM1.IELVAL(1)
          ENDIF
C
C      INFORMATION SUR L'ELEMENT FINI
C
         CALL ELQUOI(MELE,0,3,IPINF,IMODEL)
Cbp : on aurait voulu faire CALL ELQUOI(MELE,0,3,IPINF,IMODE1),
C   : mais cela ne marche evidemment pas bien...
         IF (IERR.NE.0) GOTO 9990
         INFO=IPINF
         IPTINT=INFELL(11)
         MFR   =INFELL(13)
C*OF En DIMEnsion 1, on force FORMULATION MASSIVE pour POI1
         IF (IDIM.EQ.1.AND.MELE.EQ.45) MFR=1
         IPPORE=0
         IF (MFR.EQ.33) IPPORE=NBNN
C Destruction immediate du segment
         SEGSUP,INFO
C_______________________________________________________________________
C
C RECHERCHE DES NOMS DE COMPOSANTES
C_______________________________________________________________________
          MOFORC = imodel.lnomid(2)
          if (moforc.ne.0) then
            lsupfo = .false.
            nomid = moforc
            nfor = lesobl(/2)
            nfac = 0
C            write(*,*) 'nomid deja existant dans IMODEL',IMODEL
          else
            lsupfo = .true.
            CALL IDFORC(MFR,IFOUR,MOFORC,NFOR,NFAC)
            write(ioimp,*) 'FPMASS : appel a IDFORC pour creer nomid'
          endif
          NCOMP=NFOR-NDPGE
          NOMID=MOFORC
 
Cbp     on verifie qu on a suffisamment de composantes d'effort
          NFO=0
          IF (MELE.EQ.2  .OR. MELE.EQ.3) NFO=2
          IF (MELE.EQ.31 .OR. MELE.EQ.32 .OR. MELE.EQ.33 .OR.
     &        MELE.EQ.34) NFO=3
          IF (MELE.EQ.45) NFO=1
          IF (NFO.ne.0) THEN
            IF (NCOMP.lt.NFO) GOTO 444
            DO ICOMP=1,NFO
              IF(LESOBL(ICOMP)(1:1).NE.'F') GOTO 444
            ENDDO
            GOTO 440
          ENDIF
C         -erreur
 444      CONTINUE
          write(IOIMP,*) 'on attend un MODELE avec au moins',NFO,
     &                   'composantes de FORCES !'
          write(IOIMP,*) 'Ici, on a :',(LESOBL(i),i=1,NCOMP)
          MOTERR(1:16)='MECANIQUE, ...  '
          CALL ERREUR(719)
          GOTO 9990
C         -pas d'erreur
 440      CONTINUE
C
C        CAS OU UN CHAMP PAR ELEMENT A ETE FOURNI
C        -> Verification de son support
C
         IF (IPCHM1.NE.0) THEN
           CALL QUESUP(0,IPCHM1,0,0,ISUP1,ISUP2)
           MCHEL2=IPCHM1
           MCHAM2 = MCHEL2.ICHAML(1)
           IF (ISUP2.NE.3) THEN
             IF (ISUP2.EQ.4) THEN
               CALL ERREUR(609)
               GOTO 9990
             ELSE IF (ISUP2.EQ.5) THEN
               IPTVPR = MCHAM2.IELVAL(1)
             ELSE IF (ISUP2.EQ.1.OR.ISUP2.EQ.2) THEN
               IVPRES = MCHAM2.IELVAL(1)
               CALL VALMEL(IVPRES,IPTINT,IPTVPR)
             ENDIF
           ELSE
             IPTVPR = MCHAM2.IELVAL(1)
           ENDIF
         ENDIF
C
C        INITIALISATION DU CHELEM ELEMENTAIRE DES FORCES NODALES
C
         N1=1
         L1=6
         N3=6
         SEGINI MCHELM
         TITCHE='FORCES'
         IFOCHE=IFOUR
         IPCHEL=MCHELM
C
         IMACHE(1)=IPOGEO
         INFCHE(1,1)=0
         INFCHE(1,2)=0
         INFCHE(1,3)=NHRM
         INFCHE(1,4)=IPTINT
         INFCHE(1,5)=0
         INFCHE(1,6)=3
C
C        RECHERCHE DE LA TAILLE DES MELVALS
C
         MELEME=IPOGEO
         N1PTEL=NUM(/1)
         N1EL  =NUM(/2)
         N2PTEL=0
         N2EL  =0
C
C        CREATION DU MCHAML DE LA SOUS ZONE
C
         N2=NCOMP
         SEGINI MCHAML
         ICHAML(1)=MCHAML
         NS=1
         NCOSOU=NCOMP
         SEGINI MPTVAL
         IVAFOR=MPTVAL
         nomid = MOFORC
 
         DO ICOMP=1,NCOMP
            NOMCHE(ICOMP)=LESOBL(ICOMP)
            TYPCHE(ICOMP)='REAL*8'
            SEGINI MELVAL
            IELVAL(ICOMP)=MELVAL
            IVAL(ICOMP)=MELVAL
         ENDDO
C____________________________________________________________________
C
C TRAITEMENT DES CHAMPS DE CARACTERISTIQUES
C____________________________________________________________________
          IF (MOCARA.NE.0) THEN
            IVACA1 = IVACAR
            IF (ISUPCA.EQ.1) THEN
              CALL KOMCHA(IPCHE2,IPMAIL,CONM,MOCARA,MOTYR8,0,
     &                                       INFOS,3,IVACA1)
              IF (IERR.NE.0) GOTO 9900
              CALL VALCHE(IVACA1,NCARR,IPTINT,IPPORE,MOCARA,MELE)
              IF (IERR.NE.0) THEN
                ISUPCA = 0
                GOTO 9990
              ENDIF
            ENDIF
          ENDIF
C
C        CALCUL DES FORCES NODALES EQUIVALENTES
C        DEBRANCHEMENT SUIVANT LE TYPE DES ELEMENTS
C
C           CAS DES ELEMENTS MASSIFS BIDIMENSIONNELS
C           FACES ASSOCIEES SEG2 OU SEG3
C
         IF (MELE.EQ.2.OR.MELE.EQ.3) THEN
            CALL FPMA2D(IPTVPR,IPOGEO,ipt1,IPTINT,IVAFOR,IVACA1,XP
     +                   ,netn1,ietn1)
C
C           CAS DES ELEMENTS MASSIFS TRIDIMENSIONNELS
C           FACES ASSOCIEES FAC3,FAC4,FAC6 OU FAC8
C
         ELSE IF(MELE.EQ.31.OR.MELE.EQ.32.OR.MELE.EQ.33.
     +                                   OR.MELE.EQ.34)THEN
            CALL FPMA3D(IPTVPR,IPOGEO,ipt1,IPTINT,IVAFOR,XP
     +                  ,netn1,ietn1)
C
C= Cas des elements MASSIFs UNIDIMENSIONNELs (1D)
C= Face associee : POI1 (numero 45)
         ELSE IF (MELE.EQ.45) THEN
            CALL FPMA1D(IPTVPR,IPOGEO,ipt1,IPTINT,IVAFOR,XP
     +                  ,netn1,ietn1)
C
C           ERREUR L ELEMENT N EST PAS ENCORE IMPLEMENTE
C
         ELSE
            MOTERR(1:4)=NOMTP(MELE)
            MOTERR(5:12)='FPMASS'
            CALL ERREUR (86)
            GOTO 9990
         ENDIF
C
C        ON TRANSFORME LE CHAM/ELEM EN CHAM/POIN
C        ET ON ADDITIONNE LES CHAM/POIN ELEMENTAIRES
C
         CALL CHAMPO(IPCHEL,0,IPCHPO,IPPT)
         CALL DTCHAM(IPCHEL)
         IF (IPPT.EQ.0) THEN
            GOTO 9990
         ENDIF
         IF ((ISOUS-IRRT).GT.1.OR.IB.GT.1) THEN
            CALL ADCHPO(IPCHPO,IPTFP,IPPT,1D0,1D0)
C            CALL ECRCHA(MOGEOM)
            CALL DTCHPO(IPCHPO)
C            CALL ECRCHA(MOGEOM)
            CALL DTCHPO(IPTFP)
            IF (IPPT.EQ.0) GOTO 9990
            IPTFP=IPPT
         ELSE
            IPTFP=IPCHPO
         ENDIF
 
          ISOK = 1
 
 9990     CONTINUE
          mptval = IVAFOR
          IF (IVAFOR.NE.0) SEGSUP,mptval
          nomid = MOFORC
          if (MOFORC.NE.0 .and. lsupfo) SEGSUP,nomid
          IF (IVACA1.NE.0 .AND. ISUPCA.EQ.1) THEN
            CALL DTMVAL(IVACA1,3)
          ENDIF
          IF (IPMOD1.NE.0) CALL DTMODL(IPMOD1)
          IF (ISOK.EQ.0) GOTO 9900
 
 110    CONTINUE
C- Fin de la boucle sur les sous zones de l'enveloppe
 
 9900   CONTINUE
        nomid = MOCARA
        IF (MOCARA.NE.0) SEGSUP,nomid
        IF (IVACAR.NE.0) CALL DTMVAL(IVACAR,1)
        IF (ISOK.EQ.0) GOTO 9000
 
 100  CONTINUE
 
      IF (IRRT.EQ.NSOUS) THEN
         IRET = 0
         CALL ERREUR(395)
      ELSE
         IRET = 1
      ENDIF
 
 9000 CONTINUE
      notype = MOTYR8
      SEGSUP,notype
      if (netn.ne.0) SEGSUP,netn
      if (ietn.ne.0) SEGSUP,ietn
 
      RETURN
      END
 
 

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