Télécharger envel1.eso

Retour à la liste

Numérotation des lignes :

envel1
  1. C ENVEL1 SOURCE CB215821 19/08/20 21:17:02 10287
  2. * copier sur envvol avec gestion du chamelem de valeurs associes
  3. * utilise dans trac cham
  4. C
  5. C SG 2016/07/20 Programmation comme faced2, envvo2 pour gerer les faces TRI7/QUA9
  6. *
  7. SUBROUTINE ENVEL1(MELEME,MELRES,mcoup)
  8. IMPLICIT INTEGER(I-N)
  9. IMPLICIT REAL*8 (A-H,O-Z)
  10.  
  11.  
  12.  
  13. -INC PPARAM
  14. -INC CCOPTIO
  15. -INC SMELEME
  16. -INC CCGEOME
  17. -INC SMCHAML
  18. *
  19. * Type de faces prises en compte: T3, Q4, T6, Q8, T7, Q9
  20. * Numero dans KDFAC 1 2 3 4 7 8
  21. * Pour ne pas se prendre la tête, on numerote pareil que KDFAC
  22. * Pour les 5 (non utilisé), 6 (polygone) et >8, ca restera à 0
  23. * NTYFAC=20, codé en dur dans CCGEOME pour KDFAC
  24. PARAMETER (NTYFAC=20)
  25. * Nb de faces de chaque type, sert également de compteur
  26. SEGMENT NBFAC(NTYFAC)
  27. * Tableau d'index de début fin dans les tableaux ???(NFAC)
  28. SEGMENT IDXFAC(NTYFAC+1)
  29. * Pointeurs sur des segments MELEME et MELVAL par type de face
  30. SEGMENT IPTFAC(NTYFAC)
  31. SEGMENT MLVFAC(NTYFAC)
  32. * Un segment pointant sur les IFACI et les XFACI
  33. SEGMENT IPOFAC(2,NTYFAC)
  34. * Segment IFACI contenant les noeuds, la couleur et si la face d'un
  35. * type donné est vue ou non
  36. SEGMENT IFACI(NNODE+2,NFACI)
  37. * Segment XFACI contenant les coordonnees noeuds, la couleur et si la face d'un
  38. * type donné est vue ou non
  39. SEGMENT XFACI(NNODE,NFACI)
  40. *
  41. SEGMENT IPPOL(NTPOL)
  42. SEGMENT IPREF(NTPOL)
  43. SEGMENT NTFAC(NFAC)
  44. SEGMENT KFAK(NFAC)
  45. SEGMENT NAUX(max(2,NFAC))
  46. *SG
  47. * Logique loquaf : pour les faces TRI7 et QUA9, normalement, le
  48. * dernier noeud de la face est unique à la face : il peut donc
  49. * servir de clé de hachage et on peut éviter de vérifier l'égaliteé
  50. * de tous les autres noeuds lorsque l'on teste l'égalité des faces.
  51. * C'est ce qu'on fait si loquaf=vrai.
  52. *
  53. LOGICAL LOQUAF,LOPT
  54. PARAMETER (LOQUAF=.TRUE.)
  55. * Pour chaque face dans KDFAC, le numéro d'élément associé
  56. * Ne se trouve pas dans CCGEOME, etonnant
  57. INTEGER ITYEL(NTYFAC)
  58. * T3, Q4, T6, Q8, ? , POLY, T7, Q9
  59. DATA ITYEL/4,8,6,10,0,0,7,11,12*0/
  60.  
  61. *dbg write(ioimp,*) 'coucou envel1'
  62. n2ptel=0
  63. n2el=0
  64. SEGACT MELEME
  65.  
  66. c on compte le nombre d elements dont les faces sont de type 1 2 3 4
  67. c 7 8 dans NBFAC
  68. SEGINI NBFAC
  69. NTPOL=0
  70. IPT1=MELEME
  71. SEGACT MELEME
  72. nbsour=lisous(/1)
  73. if (mcoup.ne.0) nbsour=nbsour-1
  74. DO 10 IOB=1,nbsour
  75. IPT1=LISOUS(IOB)
  76. SEGACT IPT1
  77. NBELEM=IPT1.NUM(/2)
  78. ILTEL=LTEL(1,IPT1.ITYPEL)
  79. IF (ILTEL.EQ.0) GOTO 12
  80. ILTAD=LTEL(2,IPT1.ITYPEL)
  81. DO 13 IF=1,ILTEL
  82. IFT=LDEL(1,ILTAD+IF-1)
  83. IF (IFT.EQ.6) THEN
  84. NTPOL=NTPOL+1
  85. ELSE
  86. NBFAC(IFT)=NBFAC(IFT)+NBELEM
  87. ENDIF
  88. 13 CONTINUE
  89. 12 CONTINUE
  90. 10 CONTINUE
  91.  
  92. c==== CREATION DES FACES ==============================================
  93. * Initialisation des IFACI,XFACI
  94. SEGINI IPOFAC
  95. DO ITYFAC=1,NTYFAC
  96. NNODE=KDFAC(1,ITYFAC)
  97. IF (NNODE.GT.0) THEN
  98. NFACI=NBFAC(ITYFAC)
  99. SEGINI IFACI
  100. IPOFAC(1,ITYFAC)=IFACI
  101. SEGINI XFACI
  102. IPOFAC(2,ITYFAC)=XFACI
  103. ENDIF
  104. ENDDO
  105. SEGINI IPPOL,IPREF
  106. c NBFAC sert maintenant de compteur
  107. DO ITYFAC=1,NTYFAC
  108. NBFAC(ITYFAC)=0
  109. ENDDO
  110. NTPOL=0
  111. DO 50 IOB=1,nbsour
  112. IPT1=LISOUS(IOB)
  113. * si objet en double on saute
  114. do 51 io2=1,iob-1
  115. if (ipt1.eq.lisous(io2)) goto 50
  116. 51 continue
  117. SEGACT IPT1
  118. IELIM=1
  119. IF (KSURF(IPT1.ITYPEL).EQ.IPT1.ITYPEL) THEN
  120. * face non eliminable (pas un volume)
  121. IELIM=0
  122. ENDIF
  123. melval=lisref(iob)
  124. if (melval.eq.0) goto 50
  125. segact melval
  126. lval=velche(/1)
  127. ival=velche(/2)
  128. NBELEM=IPT1.NUM(/2)
  129. ILTEL=LTEL(1,IPT1.ITYPEL)
  130. IF (ILTEL.EQ.0) GOTO 52
  131. ILTAD=LTEL(2,IPT1.ITYPEL)
  132. DO 60 IF=1,ILTEL
  133. ITYFAC=LDEL(1,ILTAD+IF-1)
  134. IAD=LDEL(2,ILTAD+IF-1)
  135. NNODE=KDFAC(1,ITYFAC)
  136. IF (NNODE.GT.0) THEN
  137. IFACI=IPOFAC(1,ITYFAC)
  138. XFACI=IPOFAC(2,ITYFAC)
  139. DO 66 IEL=1,NBELEM
  140. ielr=min(ival,iel)
  141. NBFAC(ITYFAC)=NBFAC(ITYFAC)+1
  142. j=NBFAC(ITYFAC)
  143. DO i=1,NNODE
  144. IFACI(i,j)=IPT1.NUM(LFAC(IAD+i-1),IEL)
  145. XFACI(i,j)=velche(min(lval,LFAC(IAD+i-1)),ielr)
  146. ENDDO
  147. IFACI(NNODE+1,j)=IPT1.ICOLOR(IEL)
  148. IFACI(NNODE+2,j)=IELIM
  149. 66 CONTINUE
  150. ENDIF
  151. * Avant ce if était après le 52 CONTINUE mais alors ITYFAC etait
  152. * potentiellement non initialise
  153. IF (ITYFAC.EQ.6) THEN
  154. C Polygone
  155. NTPOL = NTPOL+1
  156. IPPOL(NTPOL)= IPT1
  157. SEGINI,MELVA1 = MELVAL
  158. IPREF(NTPOL) = MELVA1
  159. ENDIF
  160. 60 CONTINUE
  161. 52 CONTINUE
  162. 50 CONTINUE
  163. C IF FAUT MAINTENANT RETASSER ET CLASSER LES TABLEAUX DES FACES
  164. C PROBLEME ELLES NE SONT PAS FORCEMENT DECRITES DE LA MEME FACON
  165. C SG 20160712 NTFAC sert de cle de hachage, elle est égale à la
  166. C somme des numeros de noeuds de la face
  167. NFAC=0
  168. SEGINI IDXFAC
  169. IDXFAC(1)=NFAC+1
  170. DO ITYFAC=1,NTYFAC
  171. NFAC=NFAC+NBFAC(ITYFAC)
  172. IDXFAC(ITYFAC+1)=NFAC+1
  173. * write(ioimp,*) 'ityfac=',ityfac,' nbfac=',NBFAC(ITYFAC)
  174. ENDDO
  175. SEGINI NTFAC,KFAK
  176. IFAC=0
  177. DO ITYFAC=1,NTYFAC
  178. NNODE=KDFAC(1,ITYFAC)
  179. IF (NNODE.GT.0) THEN
  180. LOPT=LOQUAF.AND.(ITYFAC.EQ.7.OR.ITYFAC.EQ.8)
  181. IFACI=IPOFAC(1,ITYFAC)
  182. DO I=1,NBFAC(ITYFAC)
  183. IFAC=IFAC+1
  184. IF (LOPT) THEN
  185. NTFAC(IFAC)=IFACI(NNODE,I)
  186. ELSE
  187. DO J=1,NNODE
  188. NTFAC(IFAC)=NTFAC(IFAC)+IFACI(J,I)
  189. ENDDO
  190. ENDIF
  191. KFAK(IFAC)=I
  192. ENDDO
  193. ENDIF
  194. ENDDO
  195. C IL N'Y A PLUS QU'A TRIER ET RETASSER KFAK SUIVANT NTFAC
  196. SEGINI NAUX
  197. DO 300 ITYFAC=1,NTYFAC
  198. IDEB=IDXFAC(ITYFAC)
  199. IFIN=IDXFAC(ITYFAC+1)-1
  200. IF (IFIN.LE.IDEB) GOTO 300
  201. NAUX(1)=IDEB
  202. NAUX(2)=IFIN
  203. IZ=2
  204. 208 IZ=IZ-1
  205. IF (IZ.LE.0) GOTO 209
  206. IPB=NAUX(IZ*2-1)
  207. IPH=NAUX(IZ*2)
  208. IF(IPB.GE.IPH) GOTO 208
  209. JPB=IPB-1
  210. JPH=IPH+1
  211. C CALCUL DU PIVOT
  212. NPV=0
  213. * DO 207 J=IPB,IPH
  214. * NPV=NPV+NTFAC(J)
  215. *207 CONTINUE
  216. * NPV=NPV/(IPH-IPB+1)
  217. NPV=(NTFAC(IPB)+NTFAC(IPH))/2
  218. 242 JPB=JPB+1
  219. IF (JPH.EQ.JPB) GOTO 245
  220. IF (NTFAC(JPB).LE.NPV) GOTO 243
  221. GOTO 242
  222. 243 JPH=JPH-1
  223. IF (JPH.EQ.JPB) GOTO 245
  224. IF (NTFAC(JPH).GE.NPV) GOTO 244
  225. GOTO 243
  226. 244 IAUX=KFAK(JPB)
  227. KFAK(JPB)=KFAK(JPH)
  228. KFAK(JPH)=IAUX
  229. NTAUX=NTFAC(JPB)
  230. NTFAC(JPB)=NTFAC(JPH)
  231. NTFAC(JPH)=NTAUX
  232. GOTO 242
  233. 245 IF (JPB.GE.IPB) GOTO 247
  234. JPB=JPB+1
  235. JPH=JPH+2
  236. GOTO 248
  237. 247 IF (JPH.LE.IPH) GOTO 249
  238. JPB=JPB-2
  239. JPH=JPH-1
  240. GOTO 248
  241. 249 IF (NTFAC(JPB).LE.NPV) GOTO 250
  242. IF (JPH.EQ.IPH) GOTO 251
  243. 252 JPH=JPH+1
  244. GOTO 248
  245. 250 IF (JPB.EQ.IPB) GOTO 252
  246. 251 JPB=JPB-1
  247. 248 IF (JPB.EQ.IPB) GOTO 253
  248. NAUX(2*IZ)=JPB
  249. IZ=IZ+1
  250. 253 IF (JPH.EQ.IPH) GOTO 208
  251. NAUX(2*IZ)=IPH
  252. NAUX(2*IZ-1)=JPH
  253. IZ=IZ+1
  254. GOTO 208
  255. 209 CONTINUE
  256. 300 CONTINUE
  257. C LES FACES SONT CLASSEES DANS KFAK IL FAUT ELIMINER LES FACES EN DOUBL
  258. C ELLES SONT PAR TYPE LES UNES DERRIERES LES AUTRES LES PLUS HAUTES
  259. C DEVANT
  260. IF (IIMPI.NE.0) WRITE (IOIMP,9111) (KFAK(I),NTFAC(I),I=1,NFAC)
  261. 9111 FORMAT(5(2X,2I6))
  262. DO 400 ITYFAC=1,NTYFAC
  263. IDEB=IDXFAC(ITYFAC)
  264. IFIN=IDXFAC(ITYFAC+1)-1
  265. IF (IFIN.LE.IDEB) GOTO 400
  266. NNODE=KDFAC(1,ITYFAC)
  267. * A cette etape on doit avoir nnode.gt.0
  268. IF (NNODE.LE.0) THEN
  269. CALL ERREUR(5)
  270. RETURN
  271. ENDIF
  272. LOPT=LOQUAF.AND.(ITYFAC.EQ.7.OR.ITYFAC.EQ.8)
  273. IFACI=IPOFAC(1,ITYFAC)
  274. *
  275. IFINM=IFIN-1
  276. DO 450 I1=IDEB,IFINM
  277. NTI1=NTFAC(I1)
  278. IF (NTI1.EQ.0) GOTO 450
  279. IDEB1=I1+1
  280. DO 460 I2=IDEB1,IFIN
  281. NTI2=NTFAC(I2)
  282. IF (NTI2.EQ.0) GOTO 460
  283. IF (NTI2.NE.NTI1) GOTO 450
  284. IR1=KFAK(I1)
  285. IR2=KFAK(I2)
  286. IF (IFACI(NNODE+2,IR1).EQ.0) GOTO 460
  287. IF (IFACI(NNODE+2,IR2).EQ.0) GOTO 460
  288. IF (.NOT.LOPT) THEN
  289. DO 471 J1=1,NNODE
  290. INU=IFACI(J1,IR1)
  291. DO 472 J2=1,NNODE
  292. IF (INU.EQ.IFACI(J2,IR2)) GOTO 471
  293. 472 CONTINUE
  294. GOTO 460
  295. 471 CONTINUE
  296. ENDIF
  297. C DEUX FACES EGALES ON LES SUPPRIMENT
  298. NTFAC(I1)=0
  299. NTFAC(I2)=0
  300. GOTO 450
  301. 460 CONTINUE
  302. 450 CONTINUE
  303. 400 CONTINUE
  304. *
  305. IF (IIMPI.NE.0) WRITE (IOIMP,9111) (KFAK(I),NTFAC(I),I=1,NFAC)
  306.  
  307. SEGINI IPTFAC,MLVFAC
  308. NBSOUS=0
  309. NBREF=0
  310. NBSOU2=0
  311. DO 600 ITYFAC=1,NTYFAC
  312. IDEB=IDXFAC(ITYFAC)
  313. IFIN=IDXFAC(ITYFAC+1)-1
  314. * write(ioimp,*) 'ityfac=',ityfac,' ideb=',ideb,' ifin=',ifin
  315. IF (IFIN.LT.IDEB) GOTO 600
  316. NNODE=KDFAC(1,ITYFAC)
  317. * A cette etape on doit avoir nnode.gt.0
  318. IF (NNODE.LE.0) THEN
  319. CALL ERREUR(5)
  320. RETURN
  321. ENDIF
  322. IFACI=IPOFAC(1,ITYFAC)
  323. XFACI=IPOFAC(2,ITYFAC)
  324. NBELEM=0
  325. DO 611 I=IDEB,IFIN
  326. IF (NTFAC(I).NE.0) NBELEM=NBELEM+1
  327. 611 CONTINUE
  328. * write(ioimp,*) 'nbelem=',nbelem,' nnode=',nnode
  329. IF (NBELEM.EQ.0) GOTO 600
  330. NBSOU2=NBSOU2+1
  331. NBNN=NNODE
  332. SEGINI IPT1
  333. IPT1.ITYPEL=ITYEL(ITYFAC)
  334. n1ptel=nnode
  335. n1el=nbelem
  336. segini melva1
  337. JAUX=0
  338. DO 612 J=IDEB,IFIN
  339. IF (NTFAC(J).EQ.0) GOTO 612
  340. JAUX=JAUX+1
  341. IPT1.ICOLOR(JAUX)=IFACI(NNODE+1,KFAK(J))
  342. DO 613 I=1,NBNN
  343. IPT1.NUM(I,JAUX)=IFACI(I,KFAK(J))
  344. melva1.velche(I,JAUX)=XFACI(I,KFAK(J))
  345. 613 CONTINUE
  346. 612 CONTINUE
  347. IPTFAC(ITYFAC)=IPT1
  348. * write(ioimp,*) 'ipt1=',ipt1
  349. MLVFAC(ITYFAC)=melva1
  350. 600 CONTINUE
  351. * on rajoute les points et les segments qui pouvaient etre dans le
  352. * maillage initial
  353. ipt5=0
  354. segact meleme
  355. ipt6=meleme
  356. do 710 io=1,max(1,nbsour)
  357. if (nbsour.ne.0) ipt6=lisous(io)
  358. segact ipt6
  359. if (ipt6.itypel.le.3) then
  360. nbsou2=nbsou2+1
  361. ipt5=ipt6
  362. endif
  363. 710 continue
  364. * write(ioimp,*) 'nbsou2=',nbsou2
  365. NBSOUS=NBSOU2+NTPOL
  366. if (mcoup.ne.0) nbsous=nbsous+1
  367. IF (NBSOUS.EQ.0) CALL ERREUR(26)
  368. IF (IERR.NE.0) RETURN
  369. NBREF=nbsous
  370. NBNN=0
  371. NBELEM=0
  372. SEGINI IPT5
  373. I=0
  374. DO ITYFAC=1,NTYFAC
  375. IPT1=IPTFAC(ITYFAC)
  376. melva1=MLVFAC(ITYFAC)
  377. IF (IPT1.NE.0) THEN
  378. if (melva1.eq.0) then
  379. call erreur(5)
  380. return
  381. endif
  382. I=I+1
  383. IPT5.LISOUS(I)=IPT1
  384. IPT5.LISref(I)=melva1
  385. ENDIF
  386. ENDDO
  387. segact meleme
  388. ipt1=meleme
  389. do 711 io=1,max(1,nbsour)
  390. if (nbsour.ne.0) ipt1=lisous(io)
  391. segact ipt1
  392. if (ipt1.itypel.le.3) then
  393. I=I+1
  394. IPT5.LISOUS(I)=IPT1
  395. IPT5.LISref(I)=lisref(io)
  396. endif
  397. 711 continue
  398. DO 720, IO = 1, NTPOL
  399. I= I+1
  400. IPT5.LISOUS(I) = IPPOL(IO)
  401. IPT5.LISREF(I) = IPREF(IO)
  402. 720 CONTINUE
  403. if (mcoup.ne.0) then
  404. I= I+1
  405. IPT5.LISOUS(I) = lisous(nbsour+1)
  406. IPT5.LISREF(I) = lisref(nbsour+1)
  407. endif
  408. melres=ipt5
  409.  
  410. SEGSUP IPTFAC,MLVFAC,NAUX,NTFAC,KFAK,IDXFAC,IPPOL,IPREF
  411. DO ITYFAC=1,NTYFAC
  412. IFACI=IPOFAC(1,ITYFAC)
  413. IF (IFACI.NE.0) THEN
  414. SEGSUP IFACI
  415. ENDIF
  416. XFACI=IPOFAC(2,ITYFAC)
  417. IF (XFACI.NE.0) THEN
  418. SEGSUP XFACI
  419. ENDIF
  420. ENDDO
  421. SEGSUP IPOFAC,NBFAC
  422. END
  423.  
  424.  
  425.  

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