Télécharger inclu4.eso

Retour à la liste

Numérotation des lignes :

  1. C INCLU4 SOURCE PASCAL 09/01/22 21:15:35 6263
  2. SUBROUTINE INCLU4(IPT1,IPT2,IPEX,XCRIT)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8(A-H,O-Z)
  5. -INC CCREEL
  6. -INC CCOPTIO
  7. -INC SMELEME
  8. -INC SMCOORD
  9. *
  10. * IL EST SUPPOSER QUE LE IPT2 ne contiennent quer des TET4 et que ipt1
  11. * soit des poi1,
  12. *
  13. * IPEX(I)=1 veut dire que le noeud I est interne à un element de ipt2
  14. *
  15. * POUR DECIDER SI UN POINT EST A L'INTERIEUR D'UN ELEMENT ON CALCULE
  16. * LES COORDONNEES BARYCENTRIQUES DU POINT ET IL FAUT QU'ELLES SOIENT
  17. * TOUTES POSITIVES OU QUE CELLES QUI SOIENT NEGATIVES SOIENT D'UN
  18. * ORDRE DE GRANDEUR TRES INFERIEUR AUX AUTRES ( 0.0001 FOIS ).
  19. *
  20. SEGMENT ISEG1
  21. REAL*8 XLIM(2,NBEL),YLIM(2,NBEL),ZLIM(2,NBEL)
  22. ENDSEGMENT
  23. SEGMENT ISEG3
  24. INTEGER NIZO(NZO+1)
  25. ENDSEGMENT
  26. SEGMENT ISEG4
  27. INTEGER NUMZO(NZO)
  28. ENDSEGMENT
  29. SEGMENT ISEG5
  30. INTEGER NNMEL(ILON),IDEJ(NZO)
  31. ENDSEGMENT
  32. SEGMENT ISEG6
  33. REAL*8 AM(4,4),P(4),AL(4),A(4),B(4),C(4),D(4)
  34. REAL*8 XA(3,4),XPU(3)
  35. ENDSEGMENT
  36. SEGMENT IPEX(XCOOR(/1)/(IDIM+1))
  37. *
  38. * ON CALCULE LA TAILLE MAXI D'UN ELEMENT DANS TOUTES LES DIRECTIONS
  39. * AFIN DE CREER UN ZONAGE DE L'ESPACE. EN MEME TEMPS ON CALCULE
  40. * LA DIMENSION HORS TOUT DU MAILLAGE ET ON COMMENCE A CREER
  41. * UNE NUMEROTATION LOCALE DES POINTS DU MAILLAGE IPT
  42. *
  43. IDIM1=IDIM+1
  44. XPREC=-XCRIT
  45. MELEME= IPT2
  46. SEGACT MELEME
  47. IF(ITYPEL.NE.23) THEN
  48. CALL ERREUR(5)
  49. RETURN
  50. ENDIF
  51. NBEL = NUM(/2)
  52. NBNN=NUM(/1)
  53. * WRITE(6,FMT='('' NBEL NBNN '',2I6)') NBEL,NBNN
  54. SEGINI ISEG1
  55. ILOC=0
  56. XZO=0.D0
  57. YZO=0.D0
  58. ZZO=0.D0
  59. XZA=XGRAND
  60. YZA=XGRAND
  61. ZZA=XGRAND
  62. XTOMI=XGRAND
  63. XTOMA=-XGRAND
  64. YTOMI=XGRAND
  65. YTOMA=-XGRAND
  66. ZTOMI=XGRAND
  67. ZTOMA=-XGRAND
  68. DO 1 I1=1,NBEL
  69. XMI=XGRAND
  70. YMI=XGRAND
  71. ZMI=XGRAND
  72. YMA=-XGRAND
  73. XMA=-XGRAND
  74. ZMA=-XGRAND
  75. DO 2 I2 = 1,NBNN
  76. IB=NUM(I2,I1)
  77. IA=(IB-1)*IDIM1
  78. IF(XCOOR(IA+1).LT.XMI) XMI=XCOOR(IA+1)
  79. IF(XCOOR(IA+1).GT.XMA) XMA=XCOOR(IA+1)
  80. IF(XCOOR(IA+2).LT.YMI) YMI=XCOOR(IA+2)
  81. IF(XCOOR(IA+2).GT.YMA) YMA=XCOOR(IA+2)
  82. * IF( IDIM.EQ.3 ) THEN
  83. IF(XCOOR(IA+3).LT.ZMI) ZMI=XCOOR(IA+3)
  84. IF(XCOOR(IA+3).GT.ZMA) ZMA=XCOOR(IA+3)
  85. * ENDIF
  86. 2 CONTINUE
  87. XLIM(1,I1)=XMI
  88. XLIM(2,I1)=XMA
  89. YLIM(1,I1)=YMI
  90. YLIM(2,I1)=YMA
  91. XZO=MAX (XZO,XMA-XMI)
  92. YZO=MAX (YZO,YMA-YMI)
  93. XZA=MIN(XZA,XMA-XMI)
  94. YZA=MIN(YZA,YMA-YMI)
  95. IF(XMI.LT.XTOMI) XTOMI=XMI
  96. IF(XMA.GT.XTOMA) XTOMA=XMA
  97. IF(YMI.LT.YTOMI) YTOMI=YMI
  98. IF(YMA.GT.YTOMA) YTOMA=YMA
  99. * IF(IDIM.EQ.3) THEN
  100. ZLIM(1,I1)=ZMI
  101. ZLIM(2,I1)=ZMA
  102. ZZO=MAX(ZZO,ZMA-ZMI)
  103. ZZA=MIN(ZZA,ZMA-ZMI)
  104. IF(ZMI.LT.ZTOMI) ZTOMI=ZMI
  105. IF(ZMA.GT.ZTOMA) ZTOMA=ZMA
  106. * ENDIF
  107. 1 CONTINUE
  108. * WRITE(6,FMT='(''XZO YZO '',4E12.5)') XZO,YZO
  109. * WRITE(6,FMT='(''XZA YZA '',4E12.5)') XZA,YZA
  110. * WRITE(6,FMT='(''XTOMI XTOMA '',4E12.5)') XTOMI,XTOMA
  111. * WRITE(6,FMT='(''YTOMI YTOMA '',4E12.5)') YTOMI,YTOMA
  112. XPR=MIN(XZO*1.D-2,(XTOMA-XTOMI)/2.D+4)
  113. YPR=MIN(YZO*1.D-2,(YTOMA-YTOMI)/2.D+4)
  114. C XZO=XZO*1.1
  115. C YZO=YZO*1.1
  116. XZA=XZA*0.97
  117. YZA=YZA*0.97
  118. XTOMI= XTOMI - (XTOMA-XTOMI)/1.D+4
  119. XTOMA= XTOMA + (XTOMA-XTOMI)/1.D+4
  120. YTOMI= YTOMI - (YTOMA-YTOMI)/1.D+4
  121. YTOMA= YTOMA + (YTOMA-YTOMI)/1.D+4
  122. C XZO=MIN ( XZO, XTOMA-XTOMI)
  123. C YZO=MIN ( YZO, YTOMA-YTOMI)
  124. XZA=MIN ( XZA, XTOMA-XTOMI)
  125. YZA=MIN ( YZA, YTOMA-YTOMI)
  126. NXZO=INT((XTOMA-XTOMI)/XZA) + 1
  127. NYZO=INT((YTOMA-YTOMI)/YZA) + 1
  128. XZO=XZA
  129. YZO=YZA
  130. NZZO=1
  131. * WRITE(6,FMT='('' NXZO NYZO'',2I7)') NXZO,NYZO
  132. * IF(IDIM.EQ.3) THEN
  133. ZPR=MIN(ZZO*1.D-2,(ZTOMA-ZTOMI)/2.D+4)
  134. C ZZO=ZZO*1.1
  135. ZZA=ZZA*0.97
  136. ZTOMI= ZTOMI - (ZTOMA-ZTOMI)/1.D+4
  137. ZTOMA= ZTOMA + (ZTOMA-ZTOMI)/1.D+4
  138. C ZZO=MIN ( ZZO, ZTOMA-ZTOMI)
  139. ZZA=MIN ( ZZA, ZTOMA-ZTOMI)
  140. NZZO=INT((ZTOMA-ZTOMI)/ZZA)+ 1
  141. ZZO=ZZA
  142. * WRITE(6,FMT='('' zz0,zzA,ztomi,ztoma'',4e12.5)')
  143. * $ xzo,xza,ztomi,ztoma
  144. * ENDIF
  145. * WRITE(6,FMT='('' XTOMI XTOMA YTOMI YTOMA '',4E12.5 )')
  146. * $ XTOMI, XTOMA, YTOMI ,YTOMA
  147. NXDEP=MIN(NXZO,10)
  148. NYDEP=MIN(NYZO,10)
  149. * IF(IDIM.EQ.2) THEN
  150. * IF(FLOAT(NXZO)*FLOAT(NYZO).GT.10000.) THEN
  151. * XY=SQRT(FLOAT(NXZO)*FLOAT(NYZO))/90
  152. * NXZO=MAX(INT(NXZO/XY),NXDEP)
  153. * NYZO=MAX(INT(NYZO/XY),NYDEP)
  154. * IF(FLOAT(NXZO)*FLOAT(NYZO).GT.10000.) THEN
  155. * XY=SQRT(FLOAT(NXZO)*FLOAT(NYZO))/60
  156. * NXZO=MAX(INT(NXZO/XY),NXDEP)
  157. * NYZO=MAX(INT(NYZO/XY),NYDEP)
  158. * ENDIF
  159. * XZO=(XTOMA-XTOMI)/NXZO
  160. * YZO=(YTOMA-YTOMI)/NYZO
  161. * NXZO=(XTOMA-XTOMI)/XZO +1
  162. * NYZO=(YTOMA-YTOMI)/YZO +1
  163. * ENDIF
  164. C WRITE(6,FMT='('' XZO NXZO YZO NYZO '' , E12.5,I5,E12.5,I5)')
  165. C $ XZO ,NXZO, YZO, NYZO
  166. * ELSE
  167. NZDEP=MIN(NZZO,10)
  168. * WRITE(6,FMT='('' XZO NXZO YZO NYZO ZZO NZZO'' , E12.5,I7,/,
  169. * $ E12.5,I7,E12.5,I7)')
  170. C $ XZO ,NXZO, YZO, NYZO,ZZO,NZZO
  171. IF(IIMPI.NE.0)WRITE(IOIMP,FMT='('' NXZO NYZO NZZO ''
  172. $,4I7) ') NXZO,NYZO,NZZO
  173. IF(FLOAT(NXZO)*FLOAT(NYZO)*FLOAT(NZZO).GT.25000.) THEN
  174. XYZ =(FLOAT(NXZO)*FLOAT(NYZO)*FLOAT(NZZO))**0.3333/25.
  175. NXZO=MAX(INT(FLOAT(NXZO)/XYZ),NXDEP)
  176. NYZO=MAX(INT(FLOAT(NYZO)/XYZ),NYDEP)
  177. NZZO=MAX(INT(FLOAT(NZZO)/XYZ),NZDEP)
  178. IF(IIMPI.NE.0)WRITE(IOIMP,FMT='('' NXZO NYZO NZZO ''
  179. $,4I7) ') NXZO,NYZO,NZZO
  180. IF(FLOAT(NXZO)*FLOAT(NYZO)*FLOAT(NZZO).GT.20000.) THEN
  181. XYZ =(FLOAT(NXZO)*FLOAT(NYZO)*FLOAT(NZZO))**0.3333/25.
  182. NXZO=MAX(INT(FLOAT(NXZO)/XYZ),NXDEP)
  183. NYZO=MAX(INT(FLOAT(NYZO)/XYZ),NYDEP)
  184. NZZO=MAX(INT(FLOAT(NZZO)/XYZ),NZDEP)
  185. IF(IIMPI.NE.0)WRITE(IOIMP,FMT='('' NXZO NYZO NZZO ''
  186. $,4I7) ') NXZO,NYZO,NZZO
  187. IF(FLOAT(NXZO)*FLOAT(NYZO)*FLOAT(NZZO).GT.20000.) THEN
  188. XYZ =(FLOAT(NXZO)*FLOAT(NYZO)*FLOAT(NZZO))**0.3333/25.
  189. NXZO=MAX(INT(FLOAT(NXZO)/XYZ),NXDEP)
  190. NYZO=MAX(INT(FLOAT(NYZO)/XYZ),NYDEP)
  191. NZZO=MAX(INT(FLOAT(NZZO)/XYZ),NZDEP)
  192. ENDIF
  193. ENDIF
  194. XZO=(XTOMA-XTOMI)/FLOAT(NXZO)
  195. YZO=(YTOMA-YTOMI)/FLOAT(NYZO)
  196. ZZO=(ZTOMA-ZTOMI)/FLOAT(NZZO)
  197. NXZO=INT((XTOMA-XTOMI)/XZO)+1
  198. NYZO=INT((YTOMA-YTOMI)/YZO)+1
  199. NZZO=INT((ZTOMA-ZTOMI)/ZZO)+1
  200. ENDIF
  201. * ENDIF
  202. *
  203. * ON VEUT CONSTRUIRE LA LISTE DES ELEMENTS TOUCHANT UNE ZONE
  204. * POUR CELA ON COMMENCE PAR COMPTER COMBIEN D'ELEMENT TOUCHENT
  205. * CHAQUE ZONE ET EN MEME TEMPS ON STOCKE LES ZONES TOUCHEES
  206. * PAR CHAQUE ELEMENT ET LEUR NOMBRE
  207. *
  208.  
  209. NZO=NXZO*NYZO*NZZO
  210. IF(IIMPI.NE.0)WRITE(IOIMP,FMT='('' NZO NXZO NYZO NZZO ''
  211. $,4I7) ') NZO,NXZO,NYZO,NZZO
  212. NXYZO=NXZO*NYZO
  213. * IDI=4
  214. * IF(IDIM.EQ.3) THEN
  215. * IDI=8
  216. * ENDIF
  217. SEGINI ISEG3
  218. SEGINI ISEG4
  219. DO 3 I1=1,NBEL
  220. NIZ1X=INT((XLIM(1,I1)-XTOMI-XPR)/XZO) +1
  221. NIZ1Y=INT((YLIM(1,I1)-YTOMI-YPR)/YZO) +1
  222. NIZ2X=INT((XLIM(2,I1)-XTOMI+XPR)/XZO) +1
  223. NIZ2Y=INT((YLIM(2,I1)-YTOMI+YPR)/YZO) +1
  224. * IF(IDIM.EQ.3) THEN
  225. NIZ1Z=INT((ZLIM(1,I1)-ZTOMI-ZPR)/ZZO) +1
  226. NIZ2Z=INT((ZLIM(2,I1)-ZTOMI+ZPR)/ZZO) +1
  227. DO 200 L3=NIZ1Z,NIZ2Z
  228. DO 200 L1=NIZ1Y,NIZ2Y
  229. DO 200 L2=NIZ1X,NIZ2X
  230. NIZA = L2 + ( L1-1) * NXZO + ( L3-1)*NXYZO
  231. NUMZO(NIZA) = NUMZO(NIZA) +1
  232. 200 CONTINUE
  233. * ELSE
  234. * DO 201 L1=NIZ1Y,NIZ2Y
  235. * DO 201 L2=NIZ1X,NIZ2X
  236. * NIZA = L2 + ( L1-1) * NXZO
  237. * NUMZO(NIZA) = NUMZO(NIZA) +1
  238. * 201 CONTINUE
  239. * ENDIF
  240. 3 CONTINUE
  241. *
  242. * CONSTRUCTION DU TABLEAU D'ADRESSAGE DU TABLEAU DONNANT LES
  243. * ELEMENTS CONCERNEES PAR UNE ZONE
  244. *
  245. ILON=0
  246. NIZO(1)=1
  247. DO 202 L1=1,NZO
  248. NIZO(L1+1)=NIZO(L1)+NUMZO(L1)
  249. ILON=ILON+ NUMZO(L1)
  250. 202 CONTINUE
  251. * WRITE(6,FMT='('' ILON '',I5)') ILON
  252. * WRITE(6,109) (KKK,NUMZO(KKK),(NELZO(KI,KKK),KI=1,4),KKK=1,NBEL)
  253. * 109 FORMAT(I6,I5,4I5)
  254. * WRITE(6,110)( NIZO(KI),KI=1,NZO+1)
  255. 110 FORMAT(16I5)
  256. SEGINI ISEG5
  257. DO 5 I1=1,NBEL
  258. NIZ1X=INT((XLIM(1,I1)-XTOMI-XPR)/XZO) +1
  259. NIZ1Y=INT((YLIM(1,I1)-YTOMI-YPR)/YZO) +1
  260. NIZ2X=INT((XLIM(2,I1)-XTOMI+XPR)/XZO) +1
  261. NIZ2Y=INT((YLIM(2,I1)-YTOMI+YPR)/YZO) +1
  262. * IF(IDIM.EQ.3) THEN
  263. NIZ1Z=INT((ZLIM(1,I1)-ZTOMI-ZPR)/ZZO) +1
  264. NIZ2Z=INT((ZLIM(2,I1)-ZTOMI+ZPR)/ZZO) +1
  265. DO 205 L3=NIZ1Z,NIZ2Z
  266. DO 205 L1=NIZ1Y,NIZ2Y
  267. DO 205 L2=NIZ1X,NIZ2X
  268. NIZA = L2 + ( L1-1) * NXZO + ( L3-1)*NXYZO
  269. IAD=NIZO(NIZA)+IDEJ(NIZA)
  270. NNMEL(IAD)=I1
  271. IDEJ(NIZA)=IDEJ(NIZA)+1
  272. 205 CONTINUE
  273. * ELSE
  274. * DO 203 L1=NIZ1Y,NIZ2Y
  275. * DO 203 L2=NIZ1X,NIZ2X
  276. * NIZA = L2 + ( L1-1) * NXZO
  277. * IAD=NIZO(NIZA)+IDEJ(NIZA)
  278. * NNMEL(IAD)=I1
  279. * IDEJ(NIZA)=IDEJ(NIZA)+1
  280. * 203 CONTINUE
  281. * ENDIF
  282. 5 CONTINUE
  283. *
  284. * IL NE RESTE PLUS QU'A FAIRE LE TRAVAIL PROPREMENT DIT POUR CHAQUE
  285. * POINT DE L'OBJETR MAILLAGE IPT1, ON COMMENCE PAR LE METTRE SOUS
  286. * FORME D'ELEMENTS DE TYPE POI1
  287. *
  288. SEGSUP ISEG1,ISEG4
  289. SEGACT IPT1
  290. IF(IPT1.ITYPEL.NE.1) THEN
  291. CALL CHANGE(IPT1,1)
  292. ENDIF
  293. SEGINI IPEX
  294. SEGINI ISEG6
  295. C WRITE(6,FMT='('' AVANT BOUCLE 10'')')
  296. DO 10 I=1,IPT1.NUM(/2)
  297. IP=IPT1.NUM(1,I)
  298. XPU(1)=XCOOR((IP-1)*IDIM1+1)
  299. XPU(2)=XCOOR((IP-1)*IDIM1+2)
  300. XPU(3)=XCOOR((IP-1)*IDIM1+3)
  301. * write(6,fmt='('' point x y '',i5,2e12.5)') iP,XPU(1),xpu(2)
  302. IF(XPU(1).LT.XTOMI.OR.XPU(1).GT.XTOMA) GO TO 10
  303. IF(XPU(2).LT.YTOMI.OR.XPU(2).GT.YTOMA) GO TO 10
  304. * IF(IDIM.EQ.3) THEN
  305. IF(XPU(3).LT.ZTOMI.OR.XPU(3).GT.ZTOMA) GO TO 10
  306. * ENDIF
  307. INDZO=INT((XPU(1)-XTOMI)/XZO)+ 1 +INT((XPU(2)-YTOMI)/YZO)*NXZO
  308. # +INT((XPU(3)-ZTOMI)/ZZO)*NXZO*NYZO
  309. * IF(IDIM.EQ.3) INDZO=INDZO+INT((XPU(3)-ZTOMI)/ZZO)*NXZO*NYZO
  310. IDEB=NIZO(INDZO)
  311. IFIN=NIZO(INDZO+1)-1
  312. * write(6,fmt='('' ideb ifin'',2i5)') ideb,ifin
  313. IF(IDEB.GT.IFIN) GO TO 10
  314. IEL=0
  315. DO 11 KK=IDEB,IFIN
  316. *
  317. * ON CALCULE LES COORDONNEES BARYCENTRIQUES ( AU PLUS 4 )
  318. *
  319. K=NNMEL(KK)
  320. * write(6,fmt='('' k '',i5)') k
  321. J1=NUM(1,K)
  322. J2=NUM(2,K)
  323. J3=NUM(3,K)
  324. J1IDIM=(J1-1)*IDIM1
  325. J2IDIM=(J2-1)*IDIM1
  326. J3IDIM=(J3-1)*IDIM1
  327. * IF(IDIM.EQ.3) THEN
  328. J4=NUM(4,K)
  329. J4IDIM=(J4-1)*IDIM1
  330. * ENDIF
  331. P(1)=1.D0
  332. DO 12 K1=1,NBNN
  333. AM(1,K1)=1.D0
  334. 12 CONTINUE
  335. DO 13 K1=1,IDIM
  336. P(K1+1)=XPU(K1)
  337. AM(K1+1,1)=XCOOR(J1IDIM+K1)
  338. AM(K1+1,2)=XCOOR(J2IDIM+K1)
  339. AM(K1+1,3)=XCOOR(J3IDIM+K1)
  340. * IF(IDIM.EQ.3) THEN
  341. AM(K1+1,4)=XCOOR(J4IDIM+K1)
  342. * ENDIF
  343. 13 CONTINUE
  344. * IF(IDIM.EQ.2) THEN
  345. * X1=AM(2,1)
  346. * X2=AM(2,2)
  347. * X3=AM(2,3)
  348. * Y1=AM(3,1)
  349. * Y2=AM(3,2)
  350. * Y3=AM(3,3)
  351. * X=P(2)
  352. * Y=P(3)
  353. * DETAM=X1*Y2+X2*Y3+X3*Y1-Y1*X2-Y2*X3-Y3*X1
  354. * A(1)=X2*Y3-X3*Y2
  355. * A(2)=X3*Y1-X1*Y3
  356. * A(3)=X1*Y2-X2*Y1
  357. * B(1)=Y2-Y3
  358. * B(2)=Y3-Y1
  359. * B(3)=Y1-Y2
  360. * C(1)=X3-X2
  361. * C(2)=X1-X3
  362. * C(3)=X2-X1
  363. * DO 14 IK=1,NBNN
  364. * AL(IK)=(A(IK)+B(IK)*X+C(IK)*Y)/DETAM
  365. * 14 CONTINUE
  366. * AL(4)=1.D0
  367. * ELSE
  368. X1=AM(2,1)
  369. X2=AM(2,2)
  370. X3=AM(2,3)
  371. X4=AM(2,4)
  372. Y1=AM(3,1)
  373. Y2=AM(3,2)
  374. Y3=AM(3,3)
  375. Y4=AM(3,4)
  376. Z1=AM(4,1)
  377. Z2=AM(4,2)
  378. Z3=AM(4,3)
  379. Z4=AM(4,4)
  380. X=P(2)
  381. Y=P(3)
  382. Z=P(4)
  383. DETAM=X2*Y3*Z4+X3*Y4*Z2+X4*Y2*Z3-X4*Y3*Z2-X2*Y4*Z3-X3*Y2*Z4-X1*Y3*
  384. 1Z4-X3*Y4*Z1-X4*Y1*Z3+X4*Y3*Z1+X3*Y1*Z4+X1*Y4*Z3+X1*Y2*Z4+X4*Y1*Z2+
  385. 2X2*Y4*Z1-X4*Y2*Z1-X2*Y1*Z4-X1*Y4*Z2-X1*Y2*Z3-X3*Y1*Z2-X2*Y3*Z1+X3*
  386. 3Y2*Z1+X2*Y1*Z3+X1*Y3*Z2
  387. A(1)=X2*Y3*Z4+X3*Y4*Z2+X4*Y2*Z3-X4*Y3*Z2-X2*Y4*Z3-X3*Y2*Z4
  388. A(2)=X4*Y3*Z1+X3*Y1*Z4+X1*Y4*Z3-X1*Y3*Z4-X3*Y4*Z1-X4*Y1*Z3
  389. A(3)=X1*Y2*Z4+X4*Y1*Z2+X2*Y4*Z1-X4*Y2*Z1-X2*Y1*Z4-X1*Y4*Z2
  390. A(4)=X3*Y2*Z1+X2*Y1*Z3+X1*Y3*Z2-X1*Y2*Z3-X3*Y1*Z2-X2*Y3*Z1
  391. B(1)=Y4*Z3-Y3*Z4+Y2*Z4-Y4*Z2+Y3*Z2-Y2*Z3
  392. B(2)=Y3*Z4-Y4*Z3+Y4*Z1-Y1*Z4+Y1*Z3-Y3*Z1
  393. B(3)=Y4*Z2-Y2*Z4+Y1*Z4-Y4*Z1+Y2*Z1-Y1*Z2
  394. B(4)=Y2*Z3-Y3*Z2+Y3*Z1-Y1*Z3+Y1*Z2-Y2*Z1
  395. C(1)=X3*Z4-X4*Z3+X4*Z2-X2*Z4+X2*Z3-X3*Z2
  396. C(2)=X4*Z3-X3*Z4+X1*Z4-X4*Z1+X3*Z1-X1*Z3
  397. C(3)=X2*Z4-X4*Z2+X4*Z1-X1*Z4+X1*Z2-X2*Z1
  398. C(4)=X3*Z2-X2*Z3+X1*Z3-X3*Z1+X2*Z1-X1*Z2
  399. D(1)=X4*Y3-X3*Y4+X2*Y4-X4*Y2+X3*Y2-X2*Y3
  400. D(2)=X3*Y4-X4*Y3+X4*Y1-X1*Y4+X1*Y3-X3*Y1
  401. D(3)=X4*Y2-X2*Y4+X1*Y4-X4*Y1+X2*Y1-X1*Y2
  402. D(4)=X2*Y3-X3*Y2+X3*Y1-X1*Y3+X1*Y2-X2*Y1
  403. DO 15 II=1,NBNN
  404. AL(II)=(A(II)+B(II)*X+C(II)*Y+D(II)*Z)/DETAM
  405. 15 CONTINUE
  406. * ENDIF
  407. * write(6,fmt='(''K al '',I5,4e12.3)')K,
  408. * $ al(1),al(2),al(3),al(4)
  409. IF( AL(1).GT.XPREC.AND.AL(2).GT.XPREC.AND.AL(3).GT.XPREC.
  410. $ AND.AL(4).GT.XPREC) THEN
  411. *
  412. * LE POINT EST INTERNE A L'ELEMENT
  413. *
  414. * IEL=IEL+1
  415. IPEX(IP)=1
  416. * WRITE(6,*) IP
  417. GO TO 10
  418. ENDIF
  419. 11 CONTINUE
  420. 10 CONTINUE
  421. SEGSUP ISEG5,ISEG6,ISEG3
  422. SEGDES IPT1,MELEME
  423. RETURN
  424. END
  425.  
  426.  
  427.  
  428.  
  429.  
  430.  
  431.  

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