Télécharger envvo2.eso

Retour à la liste

Numérotation des lignes :

  1. C ENVVO2 SOURCE GOUNAND 17/12/05 21:16:05 9645
  2. SUBROUTINE ENVVO2(icle)
  3. *
  4. IMPLICIT INTEGER(I-N)
  5. IMPLICIT REAL*8 (A-H,O-Z)
  6. *
  7. * 20160713 : SG nouvelle programmation de envvol, simplifiée, en vue
  8. * de traiter les faces TRI7/QUA9 nouvellement rajoutées
  9. * dans le bdata. La simplification sera egalement utile pour
  10. * les subroutines inspirees de envvol : faced2, faced, envel1,
  11. * dicho3
  12. * icle=0 : renvoie l'enveloppe (operateur ENVE)
  13. * icle=1 : renvoie les face (operateur CHAN FACE)
  14.  
  15. -INC CCOPTIO
  16. -INC SMELEME
  17. -INC CCGEOME
  18. -INC SMCOORD
  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 par type de face
  30. SEGMENT IPTFAC(NTYFAC)
  31. * Un segment pointant sur les IFACI
  32. SEGMENT IPOFAC(NTYFAC)
  33. * Segment IFACI contenant les noeuds et la couleurs des faces d'un
  34. * type donné
  35. SEGMENT IFACI(NNODE+1,NFACI)
  36. *
  37. SEGMENT IPPOL(NTPOL)
  38. SEGMENT NTFAC(NFAC)
  39. SEGMENT KFAK(NFAC)
  40. SEGMENT NAUX(max(2,NFAC))
  41. *SG
  42. * Logique loquaf : pour les faces TRI7 et QUA9, normalement, le
  43. * dernier noeud de la face est unique à la face : il peut donc
  44. * servir de clé de hachage et on peut éviter de vérifier l'égalité
  45. * de tous les autres noeuds lorsque l'on teste l'égalité des faces.
  46. * C'est ce qu'on fait si loquaf=vrai.
  47. *
  48. LOGICAL LOQUAF,LOPT
  49. PARAMETER (LOQUAF=.TRUE.)
  50. character*4 lnoid(1)
  51. * Pour chaque face dans KDFAC, le numéro d'élément associé
  52. * Ne se trouve pas dans CCGEOME, etonnant
  53. INTEGER ITYEL(NTYFAC)
  54. * T3, Q4, T6, Q8, ? , POLY, T7, Q9
  55. DATA ITYEL/4,8,6,10,0,0,7,11,12*0/
  56. data lnoid/'NOID'/
  57.  
  58. *dbg write(ioimp,*) 'coucou envvo2 IIMPI=',IIMPI
  59. c==== LECTURE ET OUVERTURE DU MELEME ==================================
  60. c et eventuelle boucle 10 sur les ojbets meleme elementaires
  61. CALL LIROBJ('MAILLAGE',MELEME,1,IRETOU)
  62. IF (IERR.NE.0) RETURN
  63. call lirmot(lnoid,1,inoid,0)
  64. c on compte le nombre d elements dont les faces sont de type 1 2 3 4
  65. c 7 8 dans NBFAC, on gère séparément les polygones
  66. SEGINI NBFAC
  67. NTPOL = 0
  68. IPT1=MELEME
  69. SEGACT MELEME
  70. DO 10 IOB=1,MAX(1,LISOUS(/1))
  71. IF (LISOUS(/1).NE.0) THEN
  72. IPT1=LISOUS(IOB)
  73. SEGACT IPT1
  74. ENDIF
  75. NBELEM=IPT1.NUM(/2)
  76. c LTEL,LDEL LFAC de CCGEOME remplis par bdata
  77. ILTEL=LTEL(1,IPT1.ITYPEL)
  78. IF (ILTEL.EQ.0) GOTO 10
  79. ILTAD=LTEL(2,IPT1.ITYPEL)
  80. c --- boucle sur les faces de chaque elements ---
  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. c --- fin de boucle sur les faces de chaque elements ---
  90. 10 CONTINUE
  91. c write(6,*) 'dimension de NFAC3,4,6,8=',NFAC3,NFAC4,NFAC6,NFAC8
  92.  
  93. c==== CREATION DES FACES ==============================================
  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(ITYFAC)=IFACI
  101. ENDIF
  102. ENDDO
  103. SEGINI IPPOL
  104. c NBFAC sert maintenant de compteur
  105. DO ITYFAC=1,NTYFAC
  106. NBFAC(ITYFAC)=0
  107. ENDDO
  108. NTPOL=0
  109. c eventuelle boucle 50 sur les objets meleme elementaires
  110. DO 50 IOB=1,MAX(1,LISOUS(/1))
  111. IF (LISOUS(/1).NE.0) IPT1=LISOUS(IOB)
  112. segact ipt1
  113. NBELEM=IPT1.NUM(/2)
  114. ILTEL=LTEL(1,IPT1.ITYPEL)
  115. IF (ILTEL.EQ.0) GOTO 52
  116. ILTAD=LTEL(2,IPT1.ITYPEL)
  117. c --- boucle 60 sur faces d'1 type d'element ------------------
  118. DO 60 IF=1,ILTEL
  119. ITYFAC=LDEL(1,ILTAD+IF-1)
  120. IAD=LDEL(2,ILTAD+IF-1)
  121. NNODE=KDFAC(1,ITYFAC)
  122. IF (NNODE.GT.0) THEN
  123. IFACI=IPOFAC(ITYFAC)
  124. c --- boucle 66 sur elements ---------------------------------
  125. DO 66 IEL=1,NBELEM
  126. NBFAC(ITYFAC)=NBFAC(ITYFAC)+1
  127. j=NBFAC(ITYFAC)
  128. DO i=1,NNODE
  129. IFACI(i,j)=IPT1.NUM(LFAC(IAD+i-1),IEL)
  130. ENDDO
  131. IFACI(NNODE+1,j)=IPT1.ICOLOR(IEL)
  132. 66 CONTINUE
  133. c --- fin de boucle 66 sur elements ---------------------------------
  134. ENDIF
  135. *SG Avant ce if était après le 52 CONTINUE mais alors ITYFAC etait
  136. * potentiellement non initialise
  137. IF (ITYFAC.EQ.6) THEN
  138. C Polygone
  139. NTPOL = NTPOL+1
  140. IPPOL(NTPOL)= IPT1
  141. ENDIF
  142. 60 CONTINUE
  143. c --- fin de boucle 60 sur faces d'1 type d'element ------------------
  144. 52 CONTINUE
  145. IF (LISOUS(/1).NE.0) SEGDES IPT1
  146. 50 CONTINUE
  147. SEGDES MELEME
  148.  
  149. C ======================================================================
  150. C IL FAUT MAINTENANT RETASSER ET CLASSER LES TABLEAUX DES FACES
  151. C PROBLEME ELLES NE SONT PAS FORCEMENT DECRITES DE LA MEME FACON
  152. C SG 20160712 NTFAC sert de cle de hachage, elle est égale à la
  153. C somme des numeros de noeuds de la face
  154. NFAC=0
  155. SEGINI IDXFAC
  156. IDXFAC(1)=NFAC+1
  157. DO ITYFAC=1,NTYFAC
  158. NFAC=NFAC+NBFAC(ITYFAC)
  159. IDXFAC(ITYFAC+1)=NFAC+1
  160. * write(ioimp,*) 'ityfac=',ityfac,' nbfac=',NBFAC(ITYFAC)
  161. ENDDO
  162. SEGINI NTFAC,KFAK
  163. IFAC=0
  164. DO ITYFAC=1,NTYFAC
  165. NNODE=KDFAC(1,ITYFAC)
  166. IF (NNODE.GT.0) THEN
  167. LOPT=LOQUAF.AND.(ITYFAC.EQ.7.OR.ITYFAC.EQ.8)
  168. IFACI=IPOFAC(ITYFAC)
  169. DO I=1,NBFAC(ITYFAC)
  170. IFAC=IFAC+1
  171. IF (LOPT) THEN
  172. NTFAC(IFAC)=IFACI(NNODE,I)
  173. ELSE
  174. DO J=1,NNODE
  175. NTFAC(IFAC)=NTFAC(IFAC)+IFACI(J,I)
  176. ENDDO
  177. ENDIF
  178. KFAK(IFAC)=I
  179. ENDDO
  180. ENDIF
  181. ENDDO
  182. C IL N'Y A PLUS QU'A TRIER ET RETASSER KFAK SUIVANT NTFAC
  183. SEGINI NAUX
  184. DO 300 ITYFAC=1,NTYFAC
  185. IDEB=IDXFAC(ITYFAC)
  186. IFIN=IDXFAC(ITYFAC+1)-1
  187. IF (IFIN.LE.IDEB) GOTO 300
  188. NAUX(1)=IDEB
  189. NAUX(2)=IFIN
  190. *
  191. IZ=2
  192. 208 IZ=IZ-1
  193. IF (IZ.LE.0) GOTO 209
  194. IPB=NAUX(IZ*2-1)
  195. IPH=NAUX(IZ*2)
  196. IF(IPB.GE.IPH) GOTO 208
  197. JPB=IPB-1
  198. JPH=IPH+1
  199. C CALCUL DU PIVOT
  200. NPV=0
  201. * DO 207 J=IPB,IPH
  202. * NPV=NPV+NTFAC(J)
  203. * 207 CONTINUE
  204. * NPV=NPV/(IPH-IPB+1)
  205. NPV=(NTFAC(IPB)+NTFAC(IPH))/2
  206. 242 JPB=JPB+1
  207. IF (JPH.EQ.JPB) GOTO 245
  208. IF (NTFAC(JPB).LE.NPV) GOTO 243
  209. GOTO 242
  210. 243 JPH=JPH-1
  211. IF (JPH.EQ.JPB) GOTO 245
  212. IF (NTFAC(JPH).GE.NPV) GOTO 244
  213. GOTO 243
  214. 244 IAUX=KFAK(JPB)
  215. KFAK(JPB)=KFAK(JPH)
  216. KFAK(JPH)=IAUX
  217. NTAUX=NTFAC(JPB)
  218. NTFAC(JPB)=NTFAC(JPH)
  219. NTFAC(JPH)=NTAUX
  220. GOTO 242
  221. 245 IF (JPB.GE.IPB) GOTO 247
  222. JPB=JPB+1
  223. JPH=JPH+2
  224. GOTO 248
  225. 247 IF (JPH.LE.IPH) GOTO 249
  226. JPB=JPB-2
  227. JPH=JPH-1
  228. GOTO 248
  229. 249 IF (NTFAC(JPB).LE.NPV) GOTO 250
  230. IF (JPH.EQ.IPH) GOTO 251
  231. 252 JPH=JPH+1
  232. GOTO 248
  233. 250 IF (JPB.EQ.IPB) GOTO 252
  234. 251 JPB=JPB-1
  235. 248 IF (JPB.EQ.IPB) GOTO 253
  236. NAUX(2*IZ)=JPB
  237. IZ=IZ+1
  238. 253 IF (JPH.EQ.IPH) GOTO 208
  239. NAUX(2*IZ)=IPH
  240. NAUX(2*IZ-1)=JPH
  241. IZ=IZ+1
  242. GOTO 208
  243. 209 CONTINUE
  244. 300 CONTINUE
  245.  
  246. C ======================================================================
  247. C LES FACES SONT CLASSEES DANS KFAK IL FAUT ELIMINER LES FACES EN DOUBL
  248. C ELLES SONT PAR TYPE LES UNES DERRIERES LES AUTRES LES PLUS HAUTES
  249. C DEVANT
  250. IF (IIMPI.NE.0) WRITE (IOIMP,9111) (KFAK(I),NTFAC(I),I=1,NFAC)
  251. 9111 FORMAT(5(2X,2I6))
  252. DO 400 ITYFAC=1,NTYFAC
  253. IDEB=IDXFAC(ITYFAC)
  254. IFIN=IDXFAC(ITYFAC+1)-1
  255. IF (IFIN.LE.IDEB) GOTO 400
  256. NNODE=KDFAC(1,ITYFAC)
  257. * A cette etape on doit avoir nnode.gt.0
  258. IF (NNODE.LE.0) THEN
  259. CALL ERREUR(5)
  260. RETURN
  261. ENDIF
  262. LOPT=LOQUAF.AND.(ITYFAC.EQ.7.OR.ITYFAC.EQ.8)
  263. IFACI=IPOFAC(ITYFAC)
  264. *
  265. IFINM=IFIN-1
  266. DO 450 I1=IDEB,IFINM
  267. NTI1=NTFAC(I1)
  268. IF (NTI1.EQ.0) GOTO 450
  269. IDEB1=I1+1
  270. DO 460 I2=IDEB1,IFIN
  271. NTI2=NTFAC(I2)
  272. IF (NTI2.EQ.0) GOTO 460
  273. IF (NTI2.NE.NTI1) GOTO 450
  274. IF (.NOT.LOPT) THEN
  275. IR1=KFAK(I1)
  276. IR2=KFAK(I2)
  277. DO 471 J1=1,NNODE
  278. INU=IFACI(J1,IR1)
  279. DO 472 J2=1,NNODE
  280. IF (INU.EQ.IFACI(J2,IR2)) GOTO 471
  281. 472 CONTINUE
  282. GOTO 460
  283. 471 CONTINUE
  284. ENDIF
  285. C DEUX FACES EGALES ON LES SUPPRIMENT (ENVE)
  286. if (icle.eq.0) then
  287. NTFAC(I1)=0
  288. NTFAC(I2)=0
  289. C DEUX FACES EGALES ON EN SUPPRIME UNE (CHAN FACE)
  290. elseif (icle.eq.1) then
  291. * il semble plus astucieux de supprimer la premiere pour le cas qui
  292. * ne doit pas arriver ici ou il y ait plus que deux faces egales
  293. NTFAC(I1)=0
  294. else
  295. call erreur(5)
  296. return
  297. endif
  298. GOTO 450
  299. 460 CONTINUE
  300. 450 CONTINUE
  301. 400 CONTINUE
  302. *
  303. IF (IIMPI.NE.0) WRITE (IOIMP,9111) (KFAK(I),NTFAC(I),I=1,NFAC)
  304.  
  305. SEGINI IPTFAC
  306. NBSOUS=0
  307. NBREF=0
  308. NBSOU2=0
  309. DO 600 ITYFAC=1,NTYFAC
  310. IDEB=IDXFAC(ITYFAC)
  311. IFIN=IDXFAC(ITYFAC+1)-1
  312. * write(ioimp,*) 'ityfac=',ityfac,' ideb=',ideb,' ifin=',ifin
  313. IF (IFIN.LT.IDEB) GOTO 600
  314. NNODE=KDFAC(1,ITYFAC)
  315. * A cette etape on doit avoir nnode.gt.0
  316. IF (NNODE.LE.0) THEN
  317. CALL ERREUR(5)
  318. RETURN
  319. ENDIF
  320. IFACI=IPOFAC(ITYFAC)
  321. NBELEM=0
  322. DO 611 I=IDEB,IFIN
  323. IF (NTFAC(I).NE.0) NBELEM=NBELEM+1
  324. 611 CONTINUE
  325. * write(ioimp,*) 'nbelem=',nbelem,' nnode=',nnode
  326. IF (NBELEM.EQ.0) GOTO 600
  327. NBSOU2=NBSOU2+1
  328. NBNN=NNODE
  329. SEGINI IPT1
  330. IPT1.ITYPEL=ITYEL(ITYFAC)
  331. JAUX=0
  332. DO 612 J=IDEB,IFIN
  333. IF (NTFAC(J).EQ.0) GOTO 612
  334. JAUX=JAUX+1
  335. IPT1.ICOLOR(JAUX)=IFACI(NNODE+1,KFAK(J))
  336. DO 613 I=1,NBNN
  337. IPT1.NUM(I,JAUX)=IFACI(I,KFAK(J))
  338. 613 CONTINUE
  339. 612 CONTINUE
  340. SEGDES IPT1
  341. IPTFAC(ITYFAC)=IPT1
  342. * write(ioimp,*) 'ipt1=',ipt1
  343. 600 CONTINUE
  344. * on rajoute les points et les segments qui pouvaient etre dans le
  345. * maillage initial
  346. ipt5=0
  347. segact meleme
  348. ipt6=meleme
  349. do 710 io=1,max(1,lisous(/1))
  350. if (lisous(/1).ne.0) ipt6=lisous(io)
  351. segact ipt6
  352. if (ipt6.itypel.le.3) then
  353. nbsou2=nbsou2+1
  354. ipt5=ipt6
  355. endif
  356. segdes ipt6
  357. 710 continue
  358. segdes meleme
  359. * write(ioimp,*) 'nbsou2=',nbsou2
  360. IF (NBSOU2.EQ.0.AND.NTPOL.EQ.0.AND.inoid.eq.0) CALL ERREUR(26)
  361. IF (NBSOU2.NE.1.OR.NTPOL.GT.0) THEN
  362. NBREF=0
  363. NBSOUS=NBSOU2+NTPOL
  364. NBNN=0
  365. NBELEM=0
  366. SEGINI IPT5
  367. I=0
  368. DO ITYFAC=1,NTYFAC
  369. IPT1=IPTFAC(ITYFAC)
  370. IF (IPT1.NE.0) THEN
  371. I=I+1
  372. IPT5.LISOUS(I)=IPT1
  373. ENDIF
  374. ENDDO
  375. segact meleme
  376. ipt1=meleme
  377. do 711 io=1,max(1,lisous(/1))
  378. if (lisous(/1).ne.0) ipt1=lisous(io)
  379. segact ipt1
  380. if (ipt1.itypel.le.3) then
  381. I=I+1
  382. IPT5.LISOUS(I)=IPT1
  383. endif
  384. segdes ipt1
  385. 711 continue
  386. DO 720, IO = 1, NTPOL
  387. I = I+1
  388. IPT5.LISOUS(I)=IPPOL(IO)
  389. 720 CONTINUE
  390. if (ipt5.lisous(/1).eq.1) then
  391. ipt6=ipt5
  392. ipt5=ipt6.lisous(1)
  393. segsup ipt6
  394. endif
  395. segdes meleme
  396. ELSE
  397. * Cas ou il n'y a qu'un meleme
  398. IPTF=0
  399. DO ITYFAC=1,NTYFAC
  400. IPTF=IPTF+IPTFAC(ITYFAC)
  401. * write(ioimp,*) 'iptf=',iptf
  402. ENDDO
  403. IPTF=IPTF+IPT5
  404. IPT5=IPTF
  405. ENDIF
  406. SEGDES IPT5
  407. CALL ECROBJ('MAILLAGE',IPT5)
  408.  
  409. SEGSUP IPTFAC
  410. SEGSUP NAUX
  411. SEGSUP NTFAC,KFAK
  412. SEGSUP IDXFAC
  413. SEGSUP IPPOL
  414. DO ITYFAC=1,NTYFAC
  415. IFACI=IPOFAC(ITYFAC)
  416. IF (IFACI.NE.0) THEN
  417. SEGSUP IFACI
  418. ENDIF
  419. ENDDO
  420. SEGSUP IPOFAC
  421. SEGSUP NBFAC
  422. RETURN
  423. END
  424.  
  425.  
  426.  
  427.  
  428.  
  429.  
  430.  
  431.  
  432.  
  433.  
  434.  
  435.  
  436.  
  437.  
  438.  
  439.  
  440.  

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