Sources Esope | Cast3M

solpri
C SOLPRI    SOURCE    PV        17/09/12    21:15:03     9542           
C---------------------------------------------------------------------|
C                                                                     |
       LOGICAL FUNCTION SOLPRI(IF1,IF2,IF3,IF4,IF5)
C                                                                     |
C      CETTE FONCTION LOGIQUE TESTE SI LE PRISME DECRIT PAR LES       |
C      FACETTES IF1..IF5 EST VIDE DE POINTS.                          |
C    LES FACETTES SUPPRIMEES SONT ORIENTES DANS L'AUTRE SENS
C      SOLPRI EST VRAI SI LE PRISME EST VIDE (DONC VALIDE)            |
C      SOLPRI EST FAUX SI LE PRISME CONTIENT UN POINT (ET EST         |
C      DONC INVALIDE)                                                 |
C                                                                     |
C---------------------------------------------------------------------|
C
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
 
-INC PPARAM
-INC CCOPTIO
-INC TDEMAIT
      LOGICAL FACET2,vervol
C
C   CHERCHER SI LES FACETTES SONT NOUVELLES OU ANCIENNES
      dimension ift(5)
      ift(1)=if1
      ift(2)=if2
      ift(3)=if3
      ift(4)=if4
      ift(5)=if5
        N1=-1
        N2=-1
        N3=-1
        N4=-1
        N5=-1
             IF (IFAT(IF1).NE.0) N1=1
             IF (IFAT(IF2).NE.0) N2=1
             IF (IFAT(IF3).NE.0) N3=1
             IF (IFAT(IF4).NE.0) N4=1
             IF (IFAT(IF5).NE.0) N5=1
*  encadrer l'hexaedre
       xmin=1e30
       xmax=-1e30
       ymin=1e30
       ymax=-1e30
       zmin=1e30
       zmax=-1e30
       do 92 j=1,5
       do 90 i=1,3
       xmin=min(xmin,xyz(1,nfc(i,ift(j))))
       xmax=max(xmax,xyz(1,nfc(i,ift(j))))
       ymin=min(ymin,xyz(2,nfc(i,ift(j))))
       ymax=max(ymax,xyz(2,nfc(i,ift(j))))
       zmin=min(zmin,xyz(3,nfc(i,ift(j))))
       zmax=max(zmax,xyz(3,nfc(i,ift(j))))
  90   continue
  92   continue
       iteste=0
*      DO 100 I=1,NPTMAX
       DO 100 I=1,0
          IF (NPF(1,I).EQ.0) GOTO 100
          DO 110 J=1,4
             IF (I.EQ.NFC(J,IF1)) GOTO 100
             IF (I.EQ.NFC(J,IF2)) GOTO 100
             IF (I.EQ.NFC(J,IF3)) GOTO 100
             IF (I.EQ.NFC(J,IF4)) GOTO 100
             IF (I.EQ.NFC(J,IF5)) GOTO 100
110       CONTINUE
       if (iteste.ne.0) then
       if (xyz(1,i).lt.xmin.or.xyz(1,i).gt.xmax) goto 100
       if (xyz(2,i).lt.ymin.or.xyz(2,i).gt.ymax) goto 100
       if (xyz(3,i).lt.zmin.or.xyz(3,i).gt.zmax) goto 100
       endif
       iteste=1
C
       V1=N1*VOL(I,NFC(1,IF1),NFC(2,IF1),NFC(3,IF1))
       V2=N1*VOL(I,NFC(1,IF1),NFC(3,IF1),NFC(4,IF1))
       V3=N2*VOL(I,NFC(1,IF2),NFC(2,IF2),NFC(3,IF2))
       V4=N2*VOL(I,NFC(1,IF2),NFC(3,IF2),NFC(4,IF2))
       V5=N3*VOL(I,NFC(1,IF3),NFC(2,IF3),NFC(3,IF3))
       V6=N4*VOL(I,NFC(1,IF4),NFC(2,IF4),NFC(3,IF4))
       V7=N5*VOL(I,NFC(1,IF5),NFC(2,IF5),NFC(3,IF5))
       V8=N5*VOL(I,NFC(1,IF5),NFC(3,IF5),NFC(4,IF5))
       VV=V1+V2+V3+V4+V5+V6+V7+V8
       IF (VV.LE.0.) then
        if (IVERB.EQ.1) write (6,*) ' prisme volume negatif '
        GOTO 120
       endif
       IF (V1.LE.-0.01*VV) GOTO 100
       IF (V2.LE.-0.01*VV) GOTO 100
       IF (V3.LE.-0.01*VV) GOTO 100
       IF (V4.LE.-0.01*VV) GOTO 100
       IF (V5.LE.-0.01*VV) GOTO 100
       IF (V6.LE.-0.01*VV) GOTO 100
       IF (V7.LE.-0.01*VV) GOTO 100
       IF (V8.LE.-0.01*VV) GOTO 100
       GOTO 120
C
100    CONTINUE
C
C      IL N'EXISTE PAS DE POINTS INTERNES AU VOLUME
       SOLPRI=.TRUE.
C  RAJOUT PV TEST INTERSECTION DES FACETTES
       IF (N1.EQ.1) SOLPRI=SOLPRI.AND.FACET2(IF1)
       IF (N2.EQ.1) SOLPRI=SOLPRI.AND.FACET2(IF2)
       IF (N3.EQ.1) SOLPRI=SOLPRI.AND.FACET2(IF3)
       IF (N4.EQ.1) SOLPRI=SOLPRI.AND.FACET2(IF4)
       IF (N5.EQ.1) SOLPRI=SOLPRI.AND.FACET2(IF5)
       IF (.NOT.SOLPRI.AND.IVERB.EQ.1)
     &       WRITE (6,*) ' SOLPRI INTERSECTION TROUVEE'
C
C  on teste maintenant qu'on ne recouvre pas un autre volume
       ip1=nfc(1,if3)
       ip2=nfc(2,if3)
       ip3=nfc(3,if3)
       ip4=nfc(1,if4)
       ip5=nfc(2,if4)
       ip6=nfc(3,if4)
*      write (6,*) ' solpyr sommets ',ip1,ip2,ip3,ip4,ip5
*      write (6,*) ' solpyr coordonnees '
*      write (6,*) xyz(1,ip1),xyz(2,ip1),xyz(3,ip1)
*      write (6,*) xyz(1,ip2),xyz(2,ip2),xyz(3,ip2)
*      write (6,*) xyz(1,ip3),xyz(2,ip3),xyz(3,ip3)
*      write (6,*) xyz(1,ip4),xyz(2,ip4),xyz(3,ip4)
*      write (6,*) xyz(1,ip5),xyz(2,ip5),xyz(3,ip5)
       xbar=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3)+xyz(1,ip4)+
     *   xyz(1,ip5)+xyz(1,ip6))/6.
       ybar=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3)+xyz(2,ip4)+
     *   xyz(2,ip5)+xyz(2,ip6))/6.
       zbar=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3)+xyz(3,ip4)+
     *   xyz(3,ip5)+xyz(3,ip6))/6.
*      write (6,*) ' bary ',xbar,ybar,zbar
       xfac=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3))/3.
       yfac=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3))/3
       zfac=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3))/3.
       xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
       xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
       xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
       solpri=solpri.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,ip5,ip6,-1,-1)
       if (.not.solpri) return
       xfac=(xyz(1,ip4)+xyz(1,ip5)+xyz(1,ip6))/3.
       yfac=(xyz(2,ip4)+xyz(2,ip5)+xyz(2,ip6))/3
       zfac=(xyz(3,ip4)+xyz(3,ip5)+xyz(3,ip6))/3.
       xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
       xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
       xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
       solpri=solpri.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,ip5,ip6,-1,-1)
       if (.not.solpri) return
       ip1=nfc(1,if1)
       ip2=nfc(2,if1)
       ip3=nfc(3,if1)
       ip4=nfc(4,if1)
       xfac=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3)+xyz(1,ip3))/4.
       yfac=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3)+xyz(2,ip3))/4.
       zfac=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3)+xyz(3,ip3))/4.
       xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
       xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
       xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
       solpri=solpri.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,ip5,ip6,-1,-1)
       if (.not.solpri) return
       ip1=nfc(1,if2)
       ip2=nfc(2,if2)
       ip3=nfc(3,if2)
       ip4=nfc(4,if2)
       xfac=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3)+xyz(1,ip3))/4.
       yfac=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3)+xyz(2,ip3))/4.
       zfac=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3)+xyz(3,ip3))/4.
       xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
       xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
       xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
       solpri=solpri.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,ip5,ip6,-1,-1)
       if (.not.solpri) return
       ip1=nfc(1,if5)
       ip2=nfc(2,if5)
       ip3=nfc(3,if5)
       ip4=nfc(4,if5)
       xfac=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3)+xyz(1,ip3))/4.
       yfac=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3)+xyz(2,ip3))/4.
       zfac=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3)+xyz(3,ip3))/4.
       xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
       xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
       xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
       solpri=solpri.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,ip5,ip6,-1,-1)
       if (.not.solpri) return
       RETURN
C
120    CONTINUE
C      LE POINT I EST INTERNE AU VOLUME
       WRITE(6,1010)I
 1010   FORMAT(' LE POINT ',I5,' EST INTERNE AU prisme CREE |')
        IF (IVERB.EQ.1) write (6,*) xyz(1,i),xyz(2,i),xyz(3,i)
       SOLPRI=.FALSE.
       RETURN
       END