Télécharger soltet.eso

Retour à la liste

Numérotation des lignes :

  1. C SOLTET SOURCE PV 17/09/12 21:15:04 9542
  2. C---------------------------------------------------------------------|
  3. C |
  4. LOGICAL FUNCTION SOLTET(IF1,IF2,IF3,IF4,IPIN)
  5. C |
  6. C CETTE FONCTION LOGIQUE TESTE SI LE TETRAEDRE DECRIT PAR LES |
  7. C FACETTES IF1..IF4 EST VIDE DE POINT. |
  8. C LES ENTIERS N1..N4 INDIQUENT L'ORIENTATION DES FACETTES |
  9. C SOLTET EST VRAI SI LE TETRAEDRE EST VIDE (DONC VALIDE) |
  10. C SOLTET EST FAUX SI LE TETRAEDRE CONTIENT UN POINT (ET EST |
  11. C DONC INVALIDE) |
  12. C |
  13. C---------------------------------------------------------------------|
  14. C
  15. IMPLICIT INTEGER(I-N)
  16. IMPLICIT REAL*8(A-H,O-Z)
  17. -INC CCOPTIO
  18. -INC TDEMAIT
  19. LOGICAL VAL,FACET2,FACET3,vervol
  20. C
  21. dimension ift(4)
  22. real*8 xo(3),xa(3),xb(3),xc(3),eps
  23. IPIN=0
  24. ift(1)=if1
  25. ift(2)=if2
  26. ift(3)=if3
  27. ift(4)=if4
  28. ip1=nfc(1,if1)
  29. ip2=nfc(2,if1)
  30. ip3=nfc(3,if1)
  31. ip4=nfc(1,if2)
  32. if (ip4.eq.ip1.or.ip4.eq.ip2.or.ip4.eq.ip3) ip4=nfc(2,if2)
  33. if (ip4.eq.ip1.or.ip4.eq.ip2.or.ip4.eq.ip3) ip4=nfc(3,if2)
  34. C
  35. C CHERCHER SI LES FACETTES SONT NOUVELLES OU ANCIENNES
  36. N1=-1
  37. N2=-1
  38. N3=-1
  39. N4=-1
  40. IF (IFAT(IF1).NE.0) N1=1
  41. IF (IFAT(IF2).NE.0) N2=1
  42. IF (IFAT(IF3).NE.0) N3=1
  43. IF (IFAT(IF4).NE.0) N4=1
  44. * encadrer le tetraedre
  45. xmin=1e30
  46. xmax=-1e30
  47. ymin=1e30
  48. ymax=-1e30
  49. zmin=1e30
  50. zmax=-1e30
  51. xmin=min(xyz(1,ip1),xyz(1,ip2),xyz(1,ip3),xyz(1,ip4))
  52. xmax=max(xyz(1,ip1),xyz(1,ip2),xyz(1,ip3),xyz(1,ip4))
  53. xmin=min(xyz(2,ip1),xyz(2,ip2),xyz(2,ip3),xyz(2,ip4))
  54. xmax=max(xyz(2,ip1),xyz(2,ip2),xyz(2,ip3),xyz(2,ip4))
  55. xmin=min(xyz(3,ip1),xyz(3,ip2),xyz(3,ip3),xyz(3,ip4))
  56. xmax=max(xyz(3,ip1),xyz(3,ip2),xyz(3,ip3),xyz(3,ip4))
  57. xd=(xmax-xmax)/2
  58. yd=(ymax-ymax)/2
  59. zd=(zmax-zmax)/2
  60. td=max(xd,yd,zd)
  61. xmin=xmin-td
  62. xmax=xmax+td
  63. ymin=ymin-td
  64. ymax=ymax+td
  65. zmin=zmin-td
  66. zmax=zmax+td
  67. iteste=0
  68. V1=N1*VOL(1,NFC(1,IF1),NFC(2,IF1),NFC(3,IF1))
  69. V2=N2*VOL(1,NFC(1,IF2),NFC(2,IF2),NFC(3,IF2))
  70. V3=N3*VOL(1,NFC(1,IF3),NFC(2,IF3),NFC(3,IF3))
  71. V4=N4*VOL(1,NFC(1,IF4),NFC(2,IF4),NFC(3,IF4))
  72. VV=V1+V2+V3+V4
  73. tm=(xyz(4,ip1)+xyz(4,ip2)+xyz(4,ip3)+xyz(4,ip4))/4
  74. tv=tm**3*1e-6
  75. IF (VV.LE.tv) then
  76. if(IVERB.EQ.1) write (6,*) ' tetraedre volume negatif '
  77. GOTO 120
  78. endif
  79. * DO 100 I=1,NPTMAX
  80. DO 100 I=1,-1
  81. IF (NPF(1,I).EQ.0) GOTO 100
  82. IF (I.EQ.IP1) GOTO 100
  83. IF (I.EQ.IP2) GOTO 100
  84. IF (I.EQ.IP3) GOTO 100
  85. IF (I.EQ.IP4) GOTO 100
  86. if (xyz(1,i).lt.xmin.or.xyz(1,i).gt.xmax) goto 100
  87. if (xyz(2,i).lt.ymin.or.xyz(2,i).gt.ymax) goto 100
  88. if (xyz(3,i).lt.zmin.or.xyz(3,i).gt.zmax) goto 100
  89. C
  90. V1=N1*VOL(I,NFC(1,IF1),NFC(2,IF1),NFC(3,IF1))
  91. V2=N2*VOL(I,NFC(1,IF2),NFC(2,IF2),NFC(3,IF2))
  92. V3=N3*VOL(I,NFC(1,IF3),NFC(2,IF3),NFC(3,IF3))
  93. V4=N4*VOL(I,NFC(1,IF4),NFC(2,IF4),NFC(3,IF4))
  94. IF (V1.LE.-volcri*VV) GOTO 100
  95. IF (V2.LE.-volcri*VV) GOTO 100
  96. IF (V3.LE.-volcri*VV) GOTO 100
  97. IF (V4.LE.-volcri*VV) GOTO 100
  98. IF (IVERB.EQ.1) WRITE(6,1010)I
  99. 1010 FORMAT(' LE POINT ',I5,' EST INTERNE AU tetraedre CREE |')
  100. IF (IVERB.EQ.1) write (6,*) ' vv,v1,v2,v3,v4 ',vv,v1,v2,v3,v4
  101. ipin=i
  102. GOTO 120
  103. C
  104. C
  105. C
  106. 100 CONTINUE
  107. C
  108. C IL N'EXISTE PAS DE POINTS INTERNES AU VOLUME
  109. SOLTET=.TRUE.
  110. C RAJOUT PV TEST ON NE CREE PAS D'ARETE A PLUS DE DEUX FACETTES
  111. * IF (N1.EQ.1) SOLTET=SOLTET.AND.FACET2(IF1)
  112. * IF (.NOT.SOLTET) RETURN
  113. * IF (N2.EQ.1) SOLTET=SOLTET.AND.FACET2(IF2)
  114. * IF (.NOT.SOLTET) RETURN
  115. * IF (N3.EQ.1) SOLTET=SOLTET.AND.FACET2(IF3)
  116. * IF (.NOT.SOLTET) RETURN
  117. * IF (N4.EQ.1) SOLTET=SOLTET.AND.FACET2(IF4)
  118. * IF (.NOT.SOLTET) RETURN
  119. C
  120. C on teste maintenant qu'on ne recouvre pas un autre volume
  121. * write (6,*) ' soltet sommets ',ip1,ip2,ip3,ip4
  122. * write (6,*) ' soltet coordonnees '
  123. * write (6,*) xyz(1,ip1),xyz(2,ip1),xyz(3,ip1)
  124. * write (6,*) xyz(1,ip2),xyz(2,ip2),xyz(3,ip2)
  125. * write (6,*) xyz(1,ip3),xyz(2,ip3),xyz(3,ip3)
  126. * write (6,*) xyz(1,ip4),xyz(2,ip4),xyz(3,ip4)
  127. xbar=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3)+xyz(1,ip4))/4.
  128. ybar=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3)+xyz(2,ip4))/4.
  129. zbar=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3)+xyz(3,ip4))/4.
  130. * write (6,*) ' bary ',xbar,ybar,zbar
  131. xfac=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip3))/3.
  132. yfac=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip3))/3
  133. zfac=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip3))/3.
  134. xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
  135. xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
  136. xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
  137. soltet=soltet.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,-1,-1,-1,-1)
  138. if (.not.soltet) return
  139. xfac=(xyz(1,ip1)+xyz(1,ip2)+xyz(1,ip4))/3.
  140. yfac=(xyz(2,ip1)+xyz(2,ip2)+xyz(2,ip4))/3
  141. zfac=(xyz(3,ip1)+xyz(3,ip2)+xyz(3,ip4))/3.
  142. xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
  143. xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
  144. xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
  145. soltet=soltet.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,-1,-1,-1,-1)
  146. if (.not.soltet) return
  147. xfac=(xyz(1,ip2)+xyz(1,ip3)+xyz(1,ip4))/3.
  148. yfac=(xyz(2,ip2)+xyz(2,ip3)+xyz(2,ip4))/3
  149. zfac=(xyz(3,ip2)+xyz(3,ip3)+xyz(3,ip4))/3.
  150. xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
  151. xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
  152. xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
  153. soltet=soltet.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,-1,-1,-1,-1)
  154. if (.not.soltet) return
  155. xfac=(xyz(1,ip3)+xyz(1,ip1)+xyz(1,ip4))/3.
  156. yfac=(xyz(2,ip3)+xyz(2,ip1)+xyz(2,ip4))/3
  157. zfac=(xyz(3,ip3)+xyz(3,ip1)+xyz(3,ip4))/3.
  158. xyz(1,nptmax+1)=xfac+0.1*(xbar-xfac)
  159. xyz(2,nptmax+1)=yfac+0.1*(ybar-yfac)
  160. xyz(3,nptmax+1)=zfac+0.1*(zbar-zfac)
  161. soltet=soltet.and.vervol(nptmax+1,ip1,ip2,ip3,ip4,-1,-1,-1,-1)
  162. if (.not.soltet) return
  163. * test sur les facettes touchees
  164. * do 200 i=1,4
  165. * if (i.eq.1) ip=ip1
  166. * if (i.eq.2) ip=ip2
  167. * if (i.eq.3) ip=ip3
  168. * if (i.eq.4) ip=ip4
  169. * do 210 j=1,40
  170. * if=npf(j,ip)
  171. * if (if.eq.0) goto 200
  172. * iprob=0
  173. * do 220 k=1,3
  174. * jp=nfc(k,if)
  175. * if (npf(4,jp).ne.0.and.npf(5,jp).eq.0) iprob=iprob+1
  176. *220 continue
  177. * if (iprob.eq.3) soltet=.false.
  178. *210 continue
  179. *200 continue
  180. * if (.not.soltet) write (6,*) ' soltet nouveau test echoue'
  181. * test sur qualite
  182. * if (qualt(ip1,ip2,ip3,ip4).lt.5e-2) then
  183. * write (6,*) ' soltet qualite insuffisante '
  184. * soltet=.false.
  185. * endif
  186. RETURN
  187. C
  188. 120 CONTINUE
  189. C LE POINT I EST INTERNE AU VOLUME
  190. IF (IVERB.EQ.1) write (6,*) ' facettes ',if1,if2,if3,if4
  191. SOLTET=.FALSE.
  192. RETURN
  193. END
  194.  
  195.  
  196.  
  197.  
  198.  

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