Télécharger itrc.eso

Retour à la liste

Numérotation des lignes :

  1. C ITRC SOURCE BP208322 16/11/18 21:17:57 9177
  2. SUBROUTINE ITRC
  3. CA Ph. MIRANDA
  4. CD 30/11/2001
  5. CC--------------------------------------------------------------------------
  6. CC
  7. CC Le sous-programme ITRC calcule par une methode analytique exacte
  8. CC l'intersection entre des segments (connus par les coordonnees de
  9. CC leurs extremites) et les facettes triangulaires d'un maillage
  10. CC (connues par les coordonnees de leurs 3 extremites et leurs
  11. CC normales).
  12. CC
  13. CC ITRC reprend les fonctionnalites de la procedure GIBIANE @INTERC
  14. CC concue par R.Mitteau (CEA/DSM/DRFC) et B.Riou (CSSI) en oct 1998.
  15. CC
  16. CC version initiale : Ph.Miranda (CSSI) en nov 2001
  17. CC modification : T.Charras (CEA/DM2S) en oct 2003
  18. CC modification : A.Moal (CSSI) en jan 2004
  19. CC
  20. CC Syntaxe generale :
  21. CC
  22. CC CHDIST MINTER = ITRC CHOLD CHNEW TOL TAB1 ;
  23. CC
  24. CC--------------------------------------------------------------------------
  25. CE En entree :
  26. CE
  27. CE CHOLD : coordonnees des points extremites initiales des segments
  28. CE (type CHPOINT a 3 composantes)
  29. CE CHNEW : coordonnees des points extremites finales des segments
  30. CE (type CHPOINT a 3 composantes appuye sur le meme support
  31. CE que CH_OLD)
  32. CE TOL : tolerance (type FLOTTANT)
  33. CE TAB1 : Table qui doit contenir les parametres suivants en entree
  34. CE
  35. CE <chamx1 (type MCHAML) coordonnee x du premier noeud de l'element
  36. CE <chamy1 (type MCHAML) coordonnee y du premier noeud de l'element
  37. CE <chamz1 (type MCHAML) coordonnee z du premier noeud de l'element
  38. CE
  39. CE <chamx2 (type MCHAML) coordonnee x du deuxieme noeud de l'element
  40. CE <chamy2 (type MCHAML) coordonnee y du deuxieme noeud de l'element
  41. CE <chamz2 (type MCHAML) coordonnee z du deuxieme noeud de l'element
  42. CE
  43. CE <chamx3 (type MCHAML) coordonnee x du troisieme noeud de l'element
  44. CE <chamy3 (type MCHAML) coordonnee y du troisieme noeud de l'element
  45. CE <chamz3 (type MCHAML) coordonnee z du troisieme noeud de l'element
  46. CE
  47. CE <cosx (type MCHAML) cosinus directeur selon x
  48. CE <cosy (type MCHAML) cosinus directeur selon y
  49. CE <cosz (type MCHAML) cosinus directeur selon z
  50. CC
  51. CS En sortie :
  52. CS
  53. CS CHDIST : CHPOINT appuye sur le support de CHOLD et contenant
  54. CS le pas si pas d'intersection et la distance entre
  55. CS point initial du segment et le point d'intersection
  56. CS a la facette sinon (type CHPOINT)
  57. CS MINTER : contient les noeuds intersectes du support de CHOLD.
  58. CS (type MAILLAGE)
  59. CS
  60. CC Remarque : il est necessaire de travailler avec les champs par
  61. CC element du maillage interceptant et non pas avec les
  62. CC coordonnées des points de l'objet maillage car ces
  63. CC champs peuvent avoir subi des changements de repere
  64. CC
  65. CC--------------------------------------------------------------------------
  66. IMPLICIT INTEGER(I-N)
  67. IMPLICIT REAL*8 (A-H,O-Z)
  68. C
  69. -INC CCOPTIO
  70. -INC SMCHPOI
  71. -INC SMCHAML
  72. -INC SMELEME
  73. -INC SMCOORD
  74. -INC CCGEOME
  75. -INC CCREEL
  76. -INC SMTABLE
  77. C
  78. SEGMENT MTRAV1
  79. REAL*8 XPN(NZ1,NBP),YPN(NZ1,NBP),ZPN(NZ1,NBP)
  80. REAL*8 XPN1(NZ1,NBP),YPN1(NZ1,NBP),ZPN1(NZ1,NBP)
  81. REAL*8 CHPAS(NZ1,NBP),CHPAS0(NZ1,NBP)
  82. LOGICAL INTER(NZ1,NBP)
  83. ENDSEGMENT
  84. C
  85. SEGMENT MTRAV2
  86. REAL*8 COORFAC(NZ2,NEL,3,3)
  87. REAL*8 TXN(NZ2,NEL),TYN(NZ2,NEL),TZN(NZ2,NEL)
  88. ENDSEGMENT
  89. C
  90. POINTEUR CHOLD.MCHPOI,CHNEW.MCHPOI
  91. POINTEUR CHX1.MCHPOI,CHX2.MCHPOI,CHX3.MCHPOI
  92. POINTEUR CHY1.MCHPOI,CHY2.MCHPOI,CHY3.MCHPOI
  93. POINTEUR CHZ1.MCHPOI,CHZ2.MCHPOI,CHZ3.MCHPOI
  94. POINTEUR TAB1.MTABLE
  95. POINTEUR ICHOLD.MELEME
  96. POINTEUR CHCOSX.MCHAML,CHCOSY.MCHAML,CHCOSZ.MCHAML
  97. C
  98. CHARACTER*8 MTYPR
  99. CHARACTER*4 REP
  100. CHARACTER*9 MOTCLE(12)
  101. logical login,logre
  102. character*4 charre
  103. C
  104. C LECTURE DES OBJETS
  105. C
  106. CALL LIROBJ('CHPOINT',CHOLD,1,IRETOU)
  107. IF (IERR .NE. 0) RETURN
  108. CALL LIROBJ('CHPOINT',CHNEW,1,IRETOU)
  109. IF (IERR .NE. 0) RETURN
  110. CALL LIRREE(TOL1,1,IRETOU)
  111. IF (IERR .NE. 0) RETURN
  112. CALL LIROBJ('TABLE ',TAB1,1,IRETOU)
  113. IF (IERR .NE. 0) RETURN
  114. C
  115. C Dimensionnement des tableaux du segment MTRAV1
  116. C
  117. MCHPOI=CHOLD
  118. SEGACT MCHPOI
  119. NZ1=IPCHP(/1)
  120. DO 50 I=1,IPCHP(/1)
  121. MSOUPO=IPCHP(I)
  122. SEGACT MSOUPO
  123. MELEME=IGEOC
  124. SEGACT MELEME
  125. IF (I .EQ. 1) NBP=NUM(/2)
  126. IF (NUM(/2) .GT. NBP) NBP=NUM(/2)
  127. SEGDES MELEME
  128. SEGDES MSOUPO
  129. 50 CONTINUE
  130. SEGINI MTRAV1
  131. C
  132. C
  133. C EXTRACTION DES COORDONNEES DES POINTS Pn et Pn+1 a partir des CHPOINT
  134. C CHOLD et CHNEW
  135. C
  136. DO 100 I=1,IPCHP(/1)
  137. MSOUPO=IPCHP(I)
  138. SEGACT MSOUPO
  139. ICHOLD=IGEOC
  140. MPOVAL=IPOVAL
  141. SEGACT MPOVAL
  142. DO 200 J=1,NOCOMP(/2)
  143. IF (NOCOMP(J).EQ. 'X ')then
  144. do 331 k=1,vpocha(/1)
  145. XPN(I,K)=VPOCHA(K,J)
  146. 331 continue
  147. elseif(NOCOMP(J).EQ. 'Y ') then
  148. do 332 k=1,vpocha(/1)
  149. YPN(I,K)=VPOCHA(K,J)
  150. 332 continue
  151. elseif(NOCOMP(J).EQ. 'Z ')then
  152. do 333 k=1,vpocha(/1)
  153. ZPN(I,K)=VPOCHA(K,J)
  154. 333 continue
  155. endif
  156. 200 CONTINUE
  157. SEGDES MPOVAL
  158. SEGDES MELEME
  159. SEGDES MSOUPO
  160. 100 CONTINUE
  161. SEGDES MCHPOI
  162. C
  163. MCHPOI=CHNEW
  164. SEGACT MCHPOI
  165. DO 110 I=1,IPCHP(/1)
  166. MSOUPO=IPCHP(I)
  167. SEGACT MSOUPO
  168. MPOVAL=IPOVAL
  169. SEGACT MPOVAL
  170. DO 210 J=1,NOCOMP(/2)
  171. IF (NOCOMP(J).EQ. 'X ')then
  172. DO 310 K=1,VPOCHA(/1)
  173. XPN1(I,K)=VPOCHA(K,J)
  174. 310 CONTINUE
  175. ELSEIF (NOCOMP(J).EQ. 'Y ')then
  176. DO 311 K=1,VPOCHA(/1)
  177. YPN1(I,K)=VPOCHA(K,J)
  178. 311 CONTINUE
  179. ELSEIF (NOCOMP(J).EQ. 'Z ')then
  180. DO 312 K=1,VPOCHA(/1)
  181. ZPN1(I,K)=VPOCHA(K,J)
  182. 312 CONTINUE
  183. ENDIF
  184. 210 CONTINUE
  185. SEGDES MPOVAL
  186. SEGDES MSOUPO
  187. 110 CONTINUE
  188. SEGDES MCHPOI
  189. C
  190. MOTCLE(1) = '<COSX '
  191. MOTCLE(2) = '<COSY '
  192. MOTCLE(3) = '<COSZ '
  193. MOTCLE(4) = '<CHAMX1 '
  194. MOTCLE(5) = '<CHAMY1 '
  195. MOTCLE(6) = '<CHAMZ1 '
  196. MOTCLE(7) = '<CHAMX2 '
  197. MOTCLE(8) = '<CHAMY2 '
  198. MOTCLE(9) = '<CHAMZ2 '
  199. MOTCLE(10) = '<CHAMX3 '
  200. MOTCLE(11) = '<CHAMY3 '
  201. MOTCLE(12) = '<CHAMZ3 '
  202. C
  203. SEGACT TAB1
  204. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(1),LOGIN,IOBIN,
  205. . 'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHCOSX)
  206. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(2),LOGIN,IOBIN,
  207. . 'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHCOSY)
  208. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(3),LOGIN,IOBIN,
  209. . 'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHCOSZ)
  210. C
  211. MCHELM=CHCOSX
  212. SEGACT MCHELM
  213. C
  214. C Dimensionnement des tableaux du segment MTRAV2
  215. C
  216. NZ2 = ICHAML(/1)
  217. DO 60 I=1,ICHAML(/1)
  218. MCHAML=ICHAML(I)
  219. SEGACT MCHAML
  220. DO 70 IC=1,IELVAL(/1)
  221. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  222. MELVAL=IELVAL(IC)
  223. SEGACT MELVAL
  224. IF (I .EQ. 1) NEL=VELCHE(/2)
  225. IF (VELCHE(/2) .GT. NEL) NEL=VELCHE(/2)
  226. SEGDES MELVAL
  227. ENDIF
  228. 70 CONTINUE
  229. SEGDES MCHAML
  230. 60 CONTINUE
  231. SEGINI MTRAV2
  232. C
  233. C Extraction des trois cosinus directeurs (CHAMELEM) a partir de TAB1
  234. C Les elements du maillage associe sont ici des POI1
  235. C
  236. DO 120 I=1,ICHAML(/1)
  237. MCHAML=ICHAML(I)
  238. SEGACT MCHAML
  239. DO 220 IC=1,IELVAL(/1)
  240. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  241. MELVAL=IELVAL(IC)
  242. SEGACT MELVAL
  243. DO 320 IEL=1,VELCHE(/2)
  244. TXN(I,IEL) = VELCHE(1,IEL)
  245. 320 CONTINUE
  246. SEGDES MELVAL
  247. ENDIF
  248. 220 CONTINUE
  249. SEGDES MCHAML
  250. 120 CONTINUE
  251. SEGDES MCHELM
  252. C
  253. MCHELM=CHCOSY
  254. SEGACT MCHELM
  255. DO 130 I=1,ICHAML(/1)
  256. MCHAML=ICHAML(I)
  257. SEGACT MCHAML
  258. DO 230 IC=1,IELVAL(/1)
  259. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  260. MELVAL=IELVAL(IC)
  261. SEGACT MELVAL
  262. DO 330 IEL=1,VELCHE(/2)
  263. TYN(I,IEL) = VELCHE(1,IEL)
  264. 330 CONTINUE
  265. SEGDES MELVAL
  266. ENDIF
  267. 230 CONTINUE
  268. SEGDES MCHAML
  269. 130 CONTINUE
  270. SEGDES MCHELM
  271. C
  272. MCHELM=CHCOSZ
  273. SEGACT MCHELM
  274. DO 140 I=1,ICHAML(/1)
  275. MCHAML=ICHAML(I)
  276. SEGACT MCHAML
  277. DO 240 IC=1,IELVAL(/1)
  278. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  279. MELVAL=IELVAL(IC)
  280. SEGACT MELVAL
  281. DO 340 IEL=1,VELCHE(/2)
  282. TZN(I,IEL) = VELCHE(1,IEL)
  283. 340 CONTINUE
  284. SEGDES MELVAL
  285. ENDIF
  286. 240 CONTINUE
  287. SEGDES MCHAML
  288. 140 CONTINUE
  289. SEGDES MCHELM
  290. C
  291. C Extraction des coordonnees des 3 points de chaque facette
  292. C
  293. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(4),LOGIN,
  294. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHX1)
  295. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(5),LOGIN,
  296. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHY1)
  297. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(6),LOGIN,
  298. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHZ1)
  299. C
  300. MCHELM=CHX1
  301. SEGACT MCHELM
  302. DO 141 I=1,ICHAML(/1)
  303. MCHAML=ICHAML(I)
  304. SEGACT MCHAML
  305. DO 241 IC=1,IELVAL(/1)
  306. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  307. MELVAL=IELVAL(IC)
  308. SEGACT MELVAL
  309. DO 341 IEL=1,VELCHE(/2)
  310. COORFAC(I,IEL,1,1) = VELCHE(1,IEL)
  311. 341 CONTINUE
  312. SEGDES MELVAL
  313. ENDIF
  314. 241 CONTINUE
  315. SEGDES MCHAML
  316. 141 CONTINUE
  317. SEGDES MCHELM
  318. C
  319. MCHELM=CHY1
  320. SEGACT MCHELM
  321. DO 142 I=1,ICHAML(/1)
  322. MCHAML=ICHAML(I)
  323. SEGACT MCHAML
  324. DO 242 IC=1,IELVAL(/1)
  325. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  326. MELVAL=IELVAL(IC)
  327. SEGACT MELVAL
  328. DO 342 IEL=1,VELCHE(/2)
  329. COORFAC(I,IEL,1,2) = VELCHE(1,IEL)
  330. 342 CONTINUE
  331. SEGDES MELVAL
  332. ENDIF
  333. 242 CONTINUE
  334. SEGDES MCHAML
  335. 142 CONTINUE
  336. SEGDES MCHELM
  337. C
  338. MCHELM=CHZ1
  339. SEGACT MCHELM
  340. DO 143 I=1,ICHAML(/1)
  341. MCHAML=ICHAML(I)
  342. SEGACT MCHAML
  343. DO 243 IC=1,IELVAL(/1)
  344. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  345. MELVAL=IELVAL(IC)
  346. SEGACT MELVAL
  347. DO 343 IEL=1,VELCHE(/2)
  348. COORFAC(I,IEL,1,3) = VELCHE(1,IEL)
  349. 343 CONTINUE
  350. SEGDES MELVAL
  351. ENDIF
  352. 243 CONTINUE
  353. SEGDES MCHAML
  354. 143 CONTINUE
  355. SEGDES MCHELM
  356. C
  357. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(7),LOGIN,
  358. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHX2)
  359. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(8),LOGIN,
  360. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHY2)
  361. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(9),LOGIN,
  362. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHZ2)
  363. C
  364. MCHELM=CHX2
  365. SEGACT MCHELM
  366. DO 144 I=1,ICHAML(/1)
  367. MCHAML=ICHAML(I)
  368. SEGACT MCHAML
  369. DO 244 IC=1,IELVAL(/1)
  370. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  371. MELVAL=IELVAL(IC)
  372. SEGACT MELVAL
  373. DO 344 IEL=1,VELCHE(/2)
  374. COORFAC(I,IEL,2,1) = VELCHE(1,IEL)
  375. 344 CONTINUE
  376. SEGDES MELVAL
  377. ENDIF
  378. 244 CONTINUE
  379. SEGDES MCHAML
  380. 144 CONTINUE
  381. SEGDES MCHELM
  382. C
  383. MCHELM=CHY2
  384. SEGACT MCHELM
  385. DO 145 I=1,ICHAML(/1)
  386. MCHAML=ICHAML(I)
  387. SEGACT MCHAML
  388. DO 245 IC=1,IELVAL(/1)
  389. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  390. MELVAL=IELVAL(IC)
  391. SEGACT MELVAL
  392. DO 345 IEL=1,VELCHE(/2)
  393. COORFAC(I,IEL,2,2) = VELCHE(1,IEL)
  394. 345 CONTINUE
  395. SEGDES MELVAL
  396. ENDIF
  397. 245 CONTINUE
  398. SEGDES MCHAML
  399. 145 CONTINUE
  400. SEGDES MCHELM
  401. C
  402. MCHELM=CHZ2
  403. SEGACT MCHELM
  404. DO 146 I=1,ICHAML(/1)
  405. MCHAML=ICHAML(I)
  406. SEGACT MCHAML
  407. DO 246 IC=1,IELVAL(/1)
  408. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  409. MELVAL=IELVAL(IC)
  410. SEGACT MELVAL
  411. DO 346 IEL=1,VELCHE(/2)
  412. COORFAC(I,IEL,2,3) = VELCHE(1,IEL)
  413. 346 CONTINUE
  414. SEGDES MELVAL
  415. ENDIF
  416. 246 CONTINUE
  417. SEGDES MCHAML
  418. 146 CONTINUE
  419. SEGDES MCHELM
  420. C
  421. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(10),LOGIN,
  422. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHX3)
  423. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(11),LOGIN,
  424. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHY3)
  425. CALL ACCTAB(TAB1,'MOT ',IVALIN,XVALIN,MOTCLE(12),LOGIN,
  426. . IOBIN,'MCHAML ',IVALRE,XVALRE,CHARRE,LOGRE,CHZ3)
  427. C
  428. MCHELM=CHX3
  429. SEGACT MCHELM
  430. DO 147 I=1,ICHAML(/1)
  431. MCHAML=ICHAML(I)
  432. SEGACT MCHAML
  433. DO 247 IC=1,IELVAL(/1)
  434. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  435. MELVAL=IELVAL(IC)
  436. SEGACT MELVAL
  437. DO 347 IEL=1,VELCHE(/2)
  438. COORFAC(I,IEL,3,1) = VELCHE(1,IEL)
  439. 347 CONTINUE
  440. SEGDES MELVAL
  441. ENDIF
  442. 247 CONTINUE
  443. SEGDES MCHAML
  444. 147 CONTINUE
  445. SEGDES MCHELM
  446. C
  447. MCHELM=CHY3
  448. SEGACT MCHELM
  449. DO 148 I=1,ICHAML(/1)
  450. MCHAML=ICHAML(I)
  451. SEGACT MCHAML
  452. DO 248 IC=1,IELVAL(/1)
  453. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  454. MELVAL=IELVAL(IC)
  455. SEGACT MELVAL
  456. DO 348 IEL=1,VELCHE(/2)
  457. COORFAC(I,IEL,3,2) = VELCHE(1,IEL)
  458. 348 CONTINUE
  459. SEGDES MELVAL
  460. ENDIF
  461. 248 CONTINUE
  462. SEGDES MCHAML
  463. 148 CONTINUE
  464. SEGDES MCHELM
  465. C
  466. MCHELM=CHZ3
  467. SEGACT MCHELM
  468. DO 149 I=1,ICHAML(/1)
  469. MCHAML=ICHAML(I)
  470. SEGACT MCHAML
  471. DO 249 IC=1,IELVAL(/1)
  472. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  473. MELVAL=IELVAL(IC)
  474. SEGACT MELVAL
  475. DO 349 IEL=1,VELCHE(/2)
  476. COORFAC(I,IEL,3,3) = VELCHE(1,IEL)
  477. 349 CONTINUE
  478. SEGDES MELVAL
  479. ENDIF
  480. 249 CONTINUE
  481. SEGDES MCHAML
  482. 149 CONTINUE
  483. SEGDES MCHELM
  484. C
  485. SEGDES TAB1
  486. C
  487. C Initialisation du tableau des indicateurs d'interception
  488. C
  489. DO 302 I=1,NZ1
  490. DO 301 K=1,NBP
  491. INTER(I,K)=.FALSE.
  492. 301 CONTINUE
  493. 302 CONTINUE
  494. C
  495. C Boucle sur les facettes du maillage OMBRANT
  496. C
  497. MCHELM=CHX1
  498. SEGACT MCHELM
  499. DO 400 I=1,ICHAML(/1)
  500. MCHAML=ICHAML(I)
  501. SEGACT MCHAML
  502. DO 450 IC=1,IELVAL(/1)
  503. IF (TYPCHE(IC)(1:6) .EQ. 'REAL*8') THEN
  504. MELVAL=IELVAL(IC)
  505. SEGACT MELVAL
  506. DO 500 J=1,VELCHE(/2)
  507. XA1 = COORFAC(I,J,1,1)
  508. YA1 = COORFAC(I,J,1,2)
  509. ZA1 = COORFAC(I,J,1,3)
  510. C
  511. XB1 = COORFAC(I,J,2,1)
  512. YB1 = COORFAC(I,J,2,2)
  513. ZB1 = COORFAC(I,J,2,3)
  514. C
  515. XC1 = COORFAC(I,J,3,1)
  516. YC1 = COORFAC(I,J,3,2)
  517. ZC1 = COORFAC(I,J,3,3)
  518. C
  519. C Boucle sur les sous-zones du maillage OMBRE
  520. C
  521. DO 550 IZ=1,NZ1
  522. C
  523. C Boucle sur les lignes de champ (ou sur les elements POI1 de OMBRE)
  524. C
  525. DO 600 K=1,NBP
  526. IF (.NOT. INTER(IZ,K)) THEN
  527. C
  528. C Produits scalaires
  529. C
  530. PS1 = ((XPN(IZ,K)-XA1)*TXN(I,J)) +
  531. . ((YPN(IZ,K)-YA1)*TYN(I,J)) +
  532. . ((ZPN(IZ,K)-ZA1)*TZN(I,J))
  533. PS2 = ((XPN1(IZ,K)-XA1)*TXN(I,J)) +
  534. . ((YPN1(IZ,K)-YA1)*TYN(I,J)) +
  535. . ((ZPN1(IZ,K)-ZA1)*TZN(I,J))
  536. C -> -->
  537. C Calcul de PnPn+1 :
  538. C
  539. DEX = XPN1(IZ,K) - XPN(IZ,K)
  540. DEY = YPN1(IZ,K) - YPN(IZ,K)
  541. DEZ = ZPN1(IZ,K) - ZPN(IZ,K)
  542. C
  543. C calcul de la distance de connection par defaut
  544. C
  545. CHPAS0(IZ,K) = (DEX*DEX + DEY*DEY + DEZ*DEZ)**.5
  546. C
  547. C Test d'intersection avec le plan de la facette
  548. C
  549. IF ( (PS1*PS2 .LT. (-1*TOL1*TOL1))
  550. . .OR. (PS2 .EQ. 0.) ) THEN
  551. C
  552. C -> --> ->
  553. C Calcul de AD1 = PnPn+1 . N :
  554. C
  555. AD1 = DEX*TXN(I,J) + DEY*TYN(I,J) + DEZ*TZN(I,J)
  556. C
  557. C Test de detection des lignes de champ paralleles a la facette
  558. C
  559. IF (ABS(AD1) .GT. TOL1) THEN
  560. C
  561. C
  562. C Calcul des points d'intersection M entre les lignes de champ et
  563. C le plan de la facette
  564. C
  565. XM = ( (XA1*DEX*TXN(I,J)) -
  566. . ((((YPN(IZ,K)-YA1)*DEX)-
  567. . (XPN(IZ,K)*DEY))*TYN(I,J)) -
  568. . ((((ZPN(IZ,K)-ZA1)*DEX)-
  569. . (XPN(IZ,K)*DEZ))*TZN(I,J)) ) / AD1
  570. C
  571. YM = ( (YA1*DEY*TYN(I,J)) -
  572. . ((((ZPN(IZ,K)-ZA1)*DEY)-
  573. . (YPN(IZ,K)*DEZ))*TZN(I,J)) -
  574. . ((((XPN(IZ,K)-XA1)*DEY)-
  575. . (YPN(IZ,K)*DEX))*TXN(I,J)) ) / AD1
  576. C
  577. ZM = ( (ZA1*DEZ*TZN(I,J)) -
  578. . ((((XPN(IZ,K)-XA1)*DEZ)-
  579. . (ZPN(IZ,K)*DEX))*TXN(I,J)) -
  580. . ((((YPN(IZ,K)-YA1)*DEZ)-
  581. . (ZPN(IZ,K)*DEY))*TYN(I,J)) ) / AD1
  582. C
  583. C Calcul des vecteurs AM, BM, CM
  584. C
  585. XAM = XM - XA1
  586. XBM = XM - XB1
  587. XCM = XM - XC1
  588. C
  589. C
  590. YAM = YM - YA1
  591. YBM = YM - YB1
  592. YCM = YM - YC1
  593. C
  594. ZAM = ZM - ZA1
  595. ZBM = ZM - ZB1
  596. ZCM = ZM - ZC1
  597. C
  598. C Denominateurs
  599. C
  600. DX = (YCM*ZBM)-(ZCM*YBM)+(ZAM*YBM)-
    . (ZAM*YCM)-(YAM*ZBM)+(YAM*ZCM)
  601. C
  602. DY = (ZCM*XBM)-(XCM*ZBM)+(XAM*ZBM)-
  603. . (XAM*ZCM)-(ZAM*XBM)+(ZAM*XCM)
  604. C
  605. DZ = (YCM*XBM)-(XCM*YBM)+(XAM*YBM)-
  606. . (XAM*YCM)-(YAM*XBM)+(YAM*XCM)
  607. C
  608. IF (ABS(DZ) .GT. TOL1) THEN
  609. D1 = DZ
  610. A1 = YCM*XBM-XCM*YBM
  611. B1 = YAM*XCM-XAM*YCM
  612. C1 = YBM*XAM-XBM*YAM
  613. ELSEIF (ABS(DY) .GT. TOL1) THEN
  614. D1 = DY
  615. A1 = ZCM*XBM-XCM*ZBM
  616. B1 = ZAM*XCM-XAM*ZCM
  617. C1 = ZBM*XAM-XBM*ZAM
  618. ELSE
  619. D1 = DX
  620. A1 = YCM*ZBM-ZCM*YBM
  621. B1 = YAM*ZCM-ZAM*YCM
  622. C1 = YBM*ZAM-ZBM*YAM
  623. ENDIF
  624. C
  625. C Coordonnees barycentriques
  626. C
  627. ALPHA = A1 / D1
  628. BETA = B1 / D1
  629. GAMMA = C1 / D1
  630. C
  631. IF ( (ALPHA .GE. 0.) .AND.
  632. . (BETA .GE. 0.) .AND.
  633. . (GAMMA .GE. 0.) ) THEN
  634. C
  635. INTER(IZ,K)=.TRUE.
  636. C
  637. C Calcul de d(Pn,M) :
  638. C
  639. CHPAS(IZ,K) = ((XPN(IZ,K)-XM)**2. +
  640. . (YPN(IZ,K)-YM)**2. + (ZPN(IZ,K)-ZM)**2.)**.5
  641. ENDIF
  642. ENDIF
  643. ENDIF
  644. ENDIF
  645. 600 CONTINUE
  646. 550 CONTINUE
  647. 500 CONTINUE
  648. SEGDES MELVAL
  649. ENDIF
  650. 450 CONTINUE
  651. SEGDES MCHAML
  652. 400 CONTINUE
  653. SEGDES MCHELM
  654. C
  655. C--------------------------------------------------------
  656. C
  657. C Creation du maillage MINTER des points d'intersection
  658. C
  659. NBNN = 1
  660. NBELEM = NBP
  661. NBSOUS = 0
  662. NBREF = 0
  663. SEGINI MELEME
  664. MINTER = MELEME
  665. ITYPEL = 1
  666. C
  667. C Construction de MINTER
  668. C
  669. IPT1 = ICHOLD
  670. SEGACT IPT1
  671. NINTER = 0
  672. DO 701 I=1,NZ1
  673. DO 700 K=1,NBP
  674. IF (INTER(I,K)) THEN
  675. NINTER = NINTER+1
  676. NUM(1,NINTER) = IPT1.NUM(1,K)
  677. ENDIF
  678. 700 CONTINUE
  679. 701 CONTINUE
  680. NBELEM = NINTER
  681. SEGADJ MELEME
  682. SEGDES MELEME,IPT1
  683. C--------------------------------------------------------
  684. C
  685. C Creation du CHPOINT CHDIST
  686. C
  687. C segment MCHPOI
  688. NAT = 1
  689. NSOUPO = NZ1
  690. SEGINI MCHPOI
  691. ICHDIS = MCHPOI
  692. IFOPOI = 2
  693. JATTRI(1) = 2
  694. C segment MSOUPO
  695. NC = 1
  696. DO 901 I=1,NZ1
  697. SEGINI MSOUPO
  698. IPCHP(I)=MSOUPO
  699. NOCOMP(1)= 'SCAL'
  700. C segment MELEME : maillage associe a CHDIST
  701. NBNN = 1
  702. NBELEM = NBP
  703. NBSOUS = 0
  704. NBREF = 0
  705. C
  706. C LE CHAMPOINT CHDIST S'APPUIE SUR LE MEME MAILLAGE QUE CHOLD
  707. C
  708. IGEOC=ICHOLD
  709. C
  710. C Construction de MPOVAL
  711. C
  712. N = NBP
  713. SEGINI MPOVAL
  714. IPOVAL=MPOVAL
  715. C
  716. C On met la valeur de la distance de connection dans le chpoint CHDIST
  717. C Ce point n'existera plus dans le prochain maillage MINTER
  718. C
  719. DO 900 K=1,NBP
  720. IF (INTER(I,K)) THEN
  721. VPOCHA(K,1)=CHPAS(I,K)
  722. ELSE
  723. VPOCHA(K,1)=CHPAS0(I,K)
  724. ENDIF
  725. 900 CONTINUE
  726. C
  727. SEGDES MPOVAL
  728. SEGDES MSOUPO
  729. 901 CONTINUE
  730. SEGDES MCHPOI
  731. C
  732. C--------------------------------------------------------
  733. C
  734. C Ecriture des nouveaux objets en sortie
  735. C
  736. CALL ECROBJ('MAILLAGE',MINTER)
  737. C
  738. CALL ECROBJ('CHPOINT ',ICHDIS)
  739. C
  740. SEGSUP MTRAV1
  741. SEGSUP MTRAV2
  742. C
  743. RETURN
  744. END
  745.  
  746.  
  747.  
  748.  
  749.  
  750.  
  751.  
  752.  
  753.  
  754.  

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