Télécharger kalpre.eso

Retour à la liste

Numérotation des lignes :

  1. C KALPRE SOURCE CB215821 16/04/21 21:17:23 8920
  2. SUBROUTINE KALPRE (MYMOD,INFOEL,SKFACE,XDEC,NELD)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8 (A-H,O-Z)
  5. C----------------------------------------------------------------------
  6. C Facteurs de forme OPTION CACHE 3D
  7. C Calcul du decoupage des faces et initialisations
  8. C SP appelle par facgen
  9. C entree:
  10. C MYMOD : pointeur de l'objet modèle
  11. C INFOEL : informations concernant le type des éléments des maillages
  12. C SKFACE : pointeur sur l'objet 'faces' pour le calcul des FF
  13. C SHC3D : pointeur sur la surface de projection
  14. C XDEC : paramŠtre pour le decoupage des faces
  15. C NELD : nombre d'éléments virtuels
  16. C (un élément COQ -> 2 facteurs de forme <=> 2 éléments virtuels)
  17. C sortie:
  18. C SKFACE : pointeur sur l'objet 'faces' pour le calcul des FF
  19. C SHC3D : pointeur sur la surface de projection
  20. C----------------------------------------------------------------------
  21. C
  22. -INC CCOPTIO
  23. -INC SMCOORD
  24. -INC SMMODEL
  25. -INC SMELEME
  26. -INC TFFOR3D
  27. C
  28. C -----------------------------------
  29. SEGMENT SDECOU
  30. INTEGER KDECOU(MFACE)
  31. ENDSEGMENT
  32. POINTEUR MYMOD.MMODEL
  33. C -----------------------------------
  34. C Stockage d'informations concernant le type des éléments des maillages
  35. SEGMENT ,INFOEL
  36. LOGICAL KCOQ(N1),KQUAD(N1)
  37. ENDSEGMENT
  38. C -----------------------------------
  39. DIMENSION A(3,4),UU(4),IN(3)
  40. DIMENSION A1(3,3),A2(3,3),U1(4),U2(4),G1(3),G2(3)
  41. DIMENSION GG(3,400),SS(400)
  42. LOGICAL ICOQ
  43. C
  44. C--------------------------------------------------------------------
  45. C
  46. C
  47. DMIN=1.E-5
  48. DMAX=1./DMIN
  49. C... nombre max de triangles de subdivision
  50. NPM=20
  51. C... quadratiques
  52. NSPA=1
  53.  
  54. IF(INFOEL.EQ.0) THEN
  55. ICOQ=.FALSE.
  56. ELSE
  57. ICOQ=.TRUE.
  58. SEGACT INFOEL
  59. ENDIF
  60. C
  61. C>>> CREATION DE L'OBJET FACE
  62. C
  63. NFACE = 0
  64. NELD = 0
  65. C
  66. SEGACT MYMOD
  67. MTYP = MYMOD.KMODEL(/1)
  68. SEGACT INFOEL
  69. C
  70. C On va compter le nombre d'éléments virtuels (égal au nombre de
  71. C facteurs de forme par ligne) : NELD , et le nombre de faces : NFACE .
  72. DO 10 ITYP=1,MTYP
  73. IMODE1 = MYMOD.KMODEL(ITYP)
  74. SEGACT IMODE1
  75. IPT1 = IMODE1.IMAMOD
  76. SEGDES IMODE1
  77. SEGACT IPT1
  78. NEL = IPT1.NUM(/2)
  79. NPGEO = IPT1.NUM(/1)
  80. SEGDES IPT1
  81. C
  82. C un élément COQ -> 2 facteurs de formes
  83. IF (KCOQ(ITYP)) NEL=2*NEL
  84. NELD = NELD + NEL
  85. C
  86. C un élément carré -> 2 éléments triangles
  87. IF(NPGEO.EQ.4.OR.NPGEO.EQ.8) NEL = 2*NEL
  88. NFACE=NFACE+NEL
  89. 10 CONTINUE
  90. C
  91. IF (KIMP.GE.1) THEN
  92. WRITE( 6,*) ' NOMBRE TOTAL D ELEMENTS ',NELD,' DE FACES ',NFACE
  93. ENDIF
  94. C
  95. C
  96. MFACE = NFACE
  97. MS = 3
  98. SEGINI SKFACE
  99. SEGINI SDECOU
  100. C
  101. C KEL : INDICE ELEMENT DECOUPAGE (-> numéro global de la facette)
  102. C IEL : INDICE ELEMENT D'ORIGINE (-> numero global de l'élément virtuel)
  103. C Ces indices sont différents car on découpe les quadrangles en triangles
  104. C
  105. KEL = 1
  106. IEL = 0
  107. C
  108. DO 100 ITYP = 1, MTYP
  109. IMODE1 = MYMOD.KMODEL(ITYP)
  110. SEGACT IMODE1
  111. IPT1 = IMODE1.IMAMOD
  112. SEGDES IMODE1
  113. SEGACT IPT1
  114.  
  115. NSGEO = IPT1.NUM(/1)
  116. IF(ICOQ.AND.KQUAD(ITYP)) THEN
  117. NS = NSGEO/2
  118. ELSE
  119. NS=NSGEO
  120. ENDIF
  121. IF(ICOQ.AND.KQUAD(ITYP)) NSPA=2
  122. NEL = IPT1.NUM(/2)
  123. C
  124. C... TRI3 ou TRI6
  125. IF (NS.EQ.3) THEN
  126. *////////
  127. * WRITE (6,*) 'Le maillage ',IPT1,
  128. * # 'est constitué d éléments triangles'
  129. *////////
  130. DO 110 I = 1,NEL
  131. IF (KCOQ(ITYP)) THEN
  132. IEL1 = IEL + (2*I-1)
  133. ELSE
  134. IEL1 = IEL + I
  135. ENDIF
  136. KCORR(KEL) = IEL1
  137.  
  138. DO 111 IS = 1,NSGEO,NSPA
  139. LS=IS
  140. IF(ICOQ.AND.KQUAD(ITYP)) LS=(IS+1)/2
  141. C** NUK(IS,KEL) = IPT1.NUM(IS,I)
  142. NUK(LS,KEL) = IPT1.NUM(IS,I)
  143. IREF = (IDIM+1)*(NUK(LS,KEL)-1)
  144. DO 112 K = 1,KES
  145. C** A(K,LS) = XCOOR(IREF+K)
  146. A(K,LS) = XCOOR(IREF+K)
  147. 112 CONTINUE
  148. *////////
  149. * WRITE (6,*) 'noeud numero',IPT1.NUM(IS,I)
  150. * WRITE (6,*) 'Coordonnées :',(A(K,IS),K=1,KES)
  151. *////////
  152. 111 CONTINUE
  153. IF (KIMP.GE.3) THEN
  154. WRITE (6,*) 'La facette ',KEL,'repose sur les points',
  155. # (IPT1.NUM(IS,I),IS=1,NS)
  156. ENDIF
  157. CALL KNORM(KES,A,NS,UU,S(KEL))
  158. DO 113 L = 1,4
  159. U(L,KEL) = UU(L)
  160. 113 CONTINUE
  161. C WRITE(6,*) ' KEL S U ',KEL,S(KEL),(U(K,NEL),K=1,4)
  162.  
  163. KEL = KEL + 1
  164. *////////
  165. IF (KCOQ(ITYP)) THEN
  166. C On remplit KCOOR , NUK , U et S pour l'élément inverse
  167. KCORR(KEL) = IEL1 + 1
  168. NUK(1,KEL) = NUK(3,KEL-1)
  169. NUK(2,KEL) = NUK(2,KEL-1)
  170. NUK(3,KEL) = NUK(1,KEL-1)
  171. DO 114 L = 1,4
  172. U(L,KEL) = - UU(L)
  173. 114 CONTINUE
  174. S(KEL) = S(KEL-1)
  175. KEL = KEL + 1
  176. ENDIF
  177. C
  178. 110 CONTINUE
  179. C
  180. IF (KCOQ(ITYP)) THEN
  181. IEL = IEL + NEL*2
  182. ELSE
  183. IEL = IEL + NEL
  184. ENDIF
  185. C
  186. C... QUA4 ou QUA8
  187. ELSE
  188. *////////
  189. * WRITE (6,*) 'Le maillage ',IPT1,
  190. * # 'est constitué d éléments carrés'
  191. *////////
  192. C
  193. DO 120 I = 1,NEL
  194. C
  195. IF (KCOQ(ITYP)) THEN
  196. IEL1 = IEL + (2*I-1)
  197. ELSE
  198. IEL1 = IEL + I
  199. ENDIF
  200. DO 121 IS = 1,NSGEO,NSPA
  201. LS=IS
  202. IF(ICOQ.AND.KQUAD(ITYP)) LS=(IS+1)/2
  203. IREF = (IDIM+1)*(IPT1.NUM(IS,I)-1)
  204. *////////
  205. * WRITE (6,*) 'noeud numero',IPT1.NUM(IS,I)
  206. * WRITE (6,*) 'IREF =',IREF
  207. *////////
  208. DO 122 K = 1,KES
  209. C** A(K,IS) = XCOOR(IREF+K)
  210. A(K,LS) = XCOOR(IREF+K)
  211. 122 CONTINUE
  212. *////////
  213. * WRITE (6,*) 'noeud numero',IPT1.NUM(IS,I)
  214. * WRITE (6,*) 'Coordonnées :',(A(K,IS),K=1,KES)
  215. *////////
  216. 121 CONTINUE
  217. C
  218. C On regarde comment on va découper le quadrangle
  219. DIA1 = 0.
  220. DO 123 K = 1,KES
  221. DIA1 = DIA1 + (A(K,1)-A(K,3))**2
  222. 123 CONTINUE
  223. DIA2 = 0.
  224. DO 124 K = 1,KES
  225. DIA2 = DIA2 + (A(K,2)-A(K,4))**2
  226. 124 CONTINUE
  227. *//////
  228. * WRITE(6,*) ' DIA1 DIA2 ',DIA1,DIA2
  229. *//////
  230.  
  231. C stockage des connectivités des traiangles
  232.  
  233. IF(ICOQ.AND.KQUAD(ITYP)) THEN
  234.  
  235. * cas des qua8
  236.  
  237. IF (DIA1.LE.DIA2) THEN
  238. NUK(1,KEL) = IPT1.NUM(1,I)
  239. NUK(2,KEL) = IPT1.NUM(3,I)
  240. NUK(3,KEL) = IPT1.NUM(5,I)
  241. KEL = KEL + 1
  242. IF (KCOQ(ITYP)) THEN
  243. NUK(1,KEL) = IPT1.NUM(5,I)
  244. NUK(2,KEL) = IPT1.NUM(3,I)
  245. NUK(3,KEL) = IPT1.NUM(1,I)
  246. KEL = KEL + 1
  247. ENDIF
  248. NUK(1,KEL) = IPT1.NUM(5,I)
  249. NUK(2,KEL) = IPT1.NUM(7,I)
  250. NUK(3,KEL) = IPT1.NUM(1,I)
  251. KEL = KEL + 1
  252. IF (KCOQ(ITYP)) THEN
  253. NUK(1,KEL) = IPT1.NUM(1,I)
  254. NUK(2,KEL) = IPT1.NUM(7,I)
  255. NUK(3,KEL) = IPT1.NUM(5,I)
  256. KEL = KEL + 1
  257. ENDIF
  258. ELSE
  259. NUK(1,KEL) = IPT1.NUM(3,I)
  260. NUK(2,KEL) = IPT1.NUM(5,I)
  261. NUK(3,KEL) = IPT1.NUM(7,I)
  262. KEL = KEL + 1
  263. IF (KCOQ(ITYP)) THEN
  264. NUK(1,KEL) = IPT1.NUM(7,I)
  265. NUK(2,KEL) = IPT1.NUM(5,I)
  266. NUK(3,KEL) = IPT1.NUM(3,I)
  267. KEL = KEL + 1
  268. ENDIF
  269. NUK(1,KEL) = IPT1.NUM(7,I)
  270. NUK(2,KEL) = IPT1.NUM(1,I)
  271. NUK(3,KEL) = IPT1.NUM(3,I)
  272. KEL = KEL + 1
  273. IF (KCOQ(ITYP)) THEN
  274. NUK(1,KEL) = IPT1.NUM(3,I)
  275. NUK(2,KEL) = IPT1.NUM(1,I)
  276. NUK(3,KEL) = IPT1.NUM(7,I)
  277. KEL = KEL + 1
  278. ENDIF
  279. ENDIF
  280.  
  281. ELSE
  282.  
  283. * cas des qua4
  284.  
  285. IF (DIA1.LE.DIA2) THEN
  286. NUK(1,KEL) = IPT1.NUM(1,I)
  287. NUK(2,KEL) = IPT1.NUM(2,I)
  288. NUK(3,KEL) = IPT1.NUM(3,I)
  289. KEL = KEL + 1
  290. IF (KCOQ(ITYP)) THEN
  291. NUK(1,KEL) = IPT1.NUM(3,I)
  292. NUK(2,KEL) = IPT1.NUM(2,I)
  293. NUK(3,KEL) = IPT1.NUM(1,I)
  294. KEL = KEL + 1
  295. ENDIF
  296. NUK(1,KEL) = IPT1.NUM(3,I)
  297. NUK(2,KEL) = IPT1.NUM(4,I)
  298. NUK(3,KEL) = IPT1.NUM(1,I)
  299. KEL = KEL + 1
  300. IF (KCOQ(ITYP)) THEN
  301. NUK(1,KEL) = IPT1.NUM(1,I)
  302. NUK(2,KEL) = IPT1.NUM(4,I)
  303. NUK(3,KEL) = IPT1.NUM(3,I)
  304. KEL = KEL + 1
  305. ENDIF
  306. ELSE
  307. NUK(1,KEL) = IPT1.NUM(2,I)
  308. NUK(2,KEL) = IPT1.NUM(3,I)
  309. NUK(3,KEL) = IPT1.NUM(4,I)
  310. KEL = KEL + 1
  311. IF (KCOQ(ITYP)) THEN
  312. NUK(1,KEL) = IPT1.NUM(4,I)
  313. NUK(2,KEL) = IPT1.NUM(3,I)
  314. NUK(3,KEL) = IPT1.NUM(2,I)
  315. KEL = KEL + 1
  316. ENDIF
  317. NUK(1,KEL) = IPT1.NUM(4,I)
  318. NUK(2,KEL) = IPT1.NUM(1,I)
  319. NUK(3,KEL) = IPT1.NUM(2,I)
  320. KEL = KEL + 1
  321. IF (KCOQ(ITYP)) THEN
  322. NUK(1,KEL) = IPT1.NUM(2,I)
  323. NUK(2,KEL) = IPT1.NUM(1,I)
  324. NUK(3,KEL) = IPT1.NUM(4,I)
  325. KEL = KEL + 1
  326. ENDIF
  327. ENDIF
  328.  
  329. ENDIF
  330.  
  331. C
  332. DO 125 IL = 1,2
  333. C
  334. IF (KCOQ(ITYP)) THEN
  335. INDICE = KEL - IL*2
  336. ELSE
  337. INDICE = KEL - IL
  338. ENDIF
  339. KCORR(INDICE) = IEL1
  340. C
  341. DO 126 IS = 1,3
  342. IREF = (IDIM+1)*(NUK(IS,INDICE)-1)
  343. DO 127 K = 1,KES
  344. A(K,IS) = XCOOR(IREF+K)
  345. C WRITE(6,*) 'A ',K,IS,A(K,IS)
  346. 127 CONTINUE
  347. 126 CONTINUE
  348. CALL KNORM(KES,A,3,UU,S(INDICE))
  349. DO 128 L = 1,4
  350. U(L,INDICE) = UU(L)
  351. 128 CONTINUE
  352. CC WRITE(6,*) ' KEL S U ',INDICE,S(INDICE),(U(K,INDICE),K=1,4)
  353.  
  354. 125 CONTINUE
  355. C
  356. IF (KCOQ(ITYP)) THEN
  357. C On remplit KCOOR , U et S pour les éléments inverses
  358. DO 135 IL = 1,2
  359. INDICE = KEL - IL*2 +1
  360. KCORR(INDICE) = IEL1 +1
  361. DO 138 L = 1,4
  362. U(L,INDICE) = - U(L,INDICE-1)
  363. 138 CONTINUE
  364. S(INDICE)=S(INDICE-1)
  365. 135 CONTINUE
  366. ENDIF
  367. C
  368. 120 CONTINUE
  369. C
  370. IF (KCOQ(ITYP)) THEN
  371. IEL = IEL + NEL*2
  372. ELSE
  373. IEL = IEL + NEL
  374. ENDIF
  375. C
  376. C
  377. ENDIF
  378. C
  379. SEGDES IPT1
  380. C
  381. 100 CONTINUE
  382. C
  383. C On a pris toute l'information nécessaire concernant le maillage
  384. SEGDES MYMOD
  385. SEGDES INFOEL
  386. *////////
  387. IF(KIMP.GE.2) THEN
  388. WRITE (6,*) 'Subdivision'
  389. ENDIF
  390. *////////
  391. C-----------------------------------------------------------
  392.  
  393. DO 200 K1=1,MFACE
  394.  
  395. DO 202 K = 1,KES
  396. G1(K)=0.
  397. DO 201 ISS = 1,MS
  398. IREF = (IDIM+1)*(NUK(ISS,K1)-1)
  399. A1(K,ISS) = XCOOR(IREF+K)
  400. G1(K) = G1(K) + A1(K,ISS)
  401. 201 CONTINUE
  402. G1(K)=G1(K)/3.
  403. 202 CONTINUE
  404. DO 203 K=1,KES+1
  405. U1(K) = U(K,K1)
  406. 203 CONTINUE
  407. C
  408. IF (XDEC.GE.0.01) THEN
  409. DK1 = DMAX
  410. DO 400 K2 = 1,MFACE
  411. C
  412. DO 412 K=1,KES+1
  413. U2(K) = U(K,K2)
  414. 412 CONTINUE
  415.  
  416. DO 212 K = 1,KES
  417. G2(K)=0.
  418. DO 211 ISS = 1,MS
  419. IREF = (IDIM+1)*(NUK(ISS,K2)-1)
  420. A2(K,ISS) = XCOOR(IREF+K)
  421. G2(K)=G2(K)+A2(K,ISS)
  422. 211 CONTINUE
  423. G2(K)=G2(K)/3.
  424. 212 CONTINUE
  425.  
  426. C>>> VISIBILITE A PRIORI
  427. C
  428.  
  429. CALL KPRIOR(KES,MS,MS,A1,A2,U1,U2,KVU)
  430. IF (KVU.NE.0) THEN
  431. D1=G1(1)-G2(1)
  432. D2=G1(2)-G2(2)
  433. D3=G1(3)-G2(3)
  434. DKK1 = SQRT(D1*D1+D2*D2+D3*D3)
  435. C WRITE(6,*) ' DKK1 CK1 ',DKK1,CK1
  436. IF(DKK1.GE.DMIN) THEN
  437. CK1 = ABS(U1(1)*D1+ U1(2)*D2 + U1(3)*D3)/DKK1
  438. IF(CK1.GE.0.01) THEN
  439. DK1 = DMIN1(DK1,DKK1)
  440. ENDIF
  441. ENDIF
  442. ENDIF
  443. 400 CONTINUE
  444. DR = DK1/XDEC
  445. ELSE
  446. DR = DMAX
  447. ENDIF
  448.  
  449. C
  450. 901 CONTINUE
  451. CALL KCREPA(DR,A1,S(K1),KES,MS,NPAT,GG,SS)
  452.  
  453. IF(NPAT.GT.NPM) THEN
  454. DR=DR*2.
  455.  
  456. IF (KIMP.GE.2) THEN
  457. WRITE(6,900) K1,NPAT
  458. 900 FORMAT(1X,' FACE ',I4,' : SUBDIVISION EXCESSIVE ',I4)
  459. ENDIF
  460. GOTO 901
  461. ENDIF
  462. NPATCH = NPAT
  463. KDECOU(K1)=NPAT
  464. MSP = KES
  465. SEGINI IPATCH
  466.  
  467. DO 130 IP = 1,NPAT
  468. SP(IP) = SS(IP)/S(K1)
  469. DO 129 K=1,KES
  470. GP(K,IP) = GG(K,IP)
  471. 129 CONTINUE
  472. C WRITE(6,*) ' SP GP ',SP(IP),GP(1,IP),GP(2,IP),GP(3,IP)
  473. 130 CONTINUE
  474.  
  475. KPATCH(K1) = IPATCH
  476. SEGDES IPATCH
  477.  
  478. 200 CONTINUE
  479. C-----------------------------------------------------------
  480.  
  481. IF(KIMP.GE.2) THEN
  482. WRITE(6,*)
  483. WRITE(6,*) 'NOMBRE DE SUBDIVISIONS PAR FACE '
  484. WRITE(6,1000) (KDECOU(I),I=1,MFACE)
  485. 1000 FORMAT(1X,10(I4))
  486. ENDIF
  487. C
  488. SEGDES SKFACE
  489. SEGSUP SDECOU
  490. C
  491. RETURN
  492. END
  493.  
  494.  
  495.  
  496.  
  497.  
  498.  
  499.  
  500.  
  501.  
  502.  
  503.  

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