Télécharger envvo2.eso

Retour à la liste

Numérotation des lignes :

  1. C ENVVO2 SOURCE CB215821 19/08/20 21:17:04 10287
  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. CALL ACTOBJ('MAILLAGE',MELEME,1)
  63. IF (IERR.NE.0) RETURN
  64. call lirmot(lnoid,1,inoid,0)
  65. c on compte le nombre d elements dont les faces sont de type 1 2 3 4
  66. c 7 8 dans NBFAC, on gère séparément les polygones
  67. SEGINI NBFAC
  68. NTPOL = 0
  69. IPT1=MELEME
  70. SEGACT MELEME
  71. DO 10 IOB=1,MAX(1,LISOUS(/1))
  72. IF (LISOUS(/1).NE.0) THEN
  73. IPT1=LISOUS(IOB)
  74. SEGACT IPT1
  75. ENDIF
  76. NBELEM=IPT1.NUM(/2)
  77. c LTEL,LDEL LFAC de CCGEOME remplis par bdata
  78. ILTEL=LTEL(1,IPT1.ITYPEL)
  79. IF (ILTEL.EQ.0) GOTO 10
  80. ILTAD=LTEL(2,IPT1.ITYPEL)
  81. c --- boucle sur les faces de chaque elements ---
  82. DO 13 IF=1,ILTEL
  83. IFT=LDEL(1,ILTAD+IF-1)
  84. IF (IFT.EQ.6) THEN
  85. NTPOL=NTPOL+1
  86. ELSE
  87. NBFAC(IFT)=NBFAC(IFT)+NBELEM
  88. ENDIF
  89. 13 CONTINUE
  90. c --- fin de boucle sur les faces de chaque elements ---
  91. 10 CONTINUE
  92. c write(6,*) 'dimension de NFAC3,4,6,8=',NFAC3,NFAC4,NFAC6,NFAC8
  93.  
  94. c==== CREATION DES FACES ==============================================
  95. SEGINI IPOFAC
  96. DO ITYFAC=1,NTYFAC
  97. NNODE=KDFAC(1,ITYFAC)
  98. IF (NNODE.GT.0) THEN
  99. NFACI=NBFAC(ITYFAC)
  100. SEGINI IFACI
  101. IPOFAC(ITYFAC)=IFACI
  102. ENDIF
  103. ENDDO
  104. SEGINI IPPOL
  105. c NBFAC sert maintenant de compteur
  106. DO ITYFAC=1,NTYFAC
  107. NBFAC(ITYFAC)=0
  108. ENDDO
  109. NTPOL=0
  110. c eventuelle boucle 50 sur les objets meleme elementaires
  111. DO 50 IOB=1,MAX(1,LISOUS(/1))
  112. IF (LISOUS(/1).NE.0) IPT1=LISOUS(IOB)
  113. segact ipt1
  114. NBELEM=IPT1.NUM(/2)
  115. ILTEL=LTEL(1,IPT1.ITYPEL)
  116. IF (ILTEL.EQ.0) GOTO 52
  117. ILTAD=LTEL(2,IPT1.ITYPEL)
  118. c --- boucle 60 sur faces d'1 type d'element ------------------
  119. DO 60 IF=1,ILTEL
  120. ITYFAC=LDEL(1,ILTAD+IF-1)
  121. IAD=LDEL(2,ILTAD+IF-1)
  122. NNODE=KDFAC(1,ITYFAC)
  123. IF (NNODE.GT.0) THEN
  124. IFACI=IPOFAC(ITYFAC)
  125. c --- boucle 66 sur elements ---------------------------------
  126. DO 66 IEL=1,NBELEM
  127. NBFAC(ITYFAC)=NBFAC(ITYFAC)+1
  128. j=NBFAC(ITYFAC)
  129. DO i=1,NNODE
  130. IFACI(i,j)=IPT1.NUM(LFAC(IAD+i-1),IEL)
  131. ENDDO
  132. IFACI(NNODE+1,j)=IPT1.ICOLOR(IEL)
  133. 66 CONTINUE
  134. c --- fin de boucle 66 sur elements ---------------------------------
  135. ENDIF
  136. *SG Avant ce if était après le 52 CONTINUE mais alors ITYFAC etait
  137. * potentiellement non initialise
  138. IF (ITYFAC.EQ.6) THEN
  139. C Polygone
  140. NTPOL = NTPOL+1
  141. IPPOL(NTPOL)= IPT1
  142. ENDIF
  143. 60 CONTINUE
  144. c --- fin de boucle 60 sur faces d'1 type d'element ------------------
  145. 52 CONTINUE
  146. 50 CONTINUE
  147.  
  148. C ======================================================================
  149. C IL FAUT MAINTENANT RETASSER ET CLASSER LES TABLEAUX DES FACES
  150. C PROBLEME ELLES NE SONT PAS FORCEMENT DECRITES DE LA MEME FACON
  151. C SG 20160712 NTFAC sert de cle de hachage, elle est égale à la
  152. C somme des numeros de noeuds de la face
  153. NFAC=0
  154. SEGINI IDXFAC
  155. IDXFAC(1)=NFAC+1
  156. DO ITYFAC=1,NTYFAC
  157. NFAC=NFAC+NBFAC(ITYFAC)
  158. IDXFAC(ITYFAC+1)=NFAC+1
  159. * write(ioimp,*) 'ityfac=',ityfac,' nbfac=',NBFAC(ITYFAC)
  160. ENDDO
  161. SEGINI NTFAC,KFAK
  162. IFAC=0
  163. DO ITYFAC=1,NTYFAC
  164. NNODE=KDFAC(1,ITYFAC)
  165. IF (NNODE.GT.0) THEN
  166. LOPT=LOQUAF.AND.(ITYFAC.EQ.7.OR.ITYFAC.EQ.8)
  167. IFACI=IPOFAC(ITYFAC)
  168. DO I=1,NBFAC(ITYFAC)
  169. IFAC=IFAC+1
  170. IF (LOPT) THEN
  171. NTFAC(IFAC)=IFACI(NNODE,I)
  172. ELSE
  173. DO J=1,NNODE
  174. NTFAC(IFAC)=NTFAC(IFAC)+IFACI(J,I)
  175. ENDDO
  176. ENDIF
  177. KFAK(IFAC)=I
  178. ENDDO
  179. ENDIF
  180. ENDDO
  181. C IL N'Y A PLUS QU'A TRIER ET RETASSER KFAK SUIVANT NTFAC
  182. SEGINI NAUX
  183. DO 300 ITYFAC=1,NTYFAC
  184. IDEB=IDXFAC(ITYFAC)
  185. IFIN=IDXFAC(ITYFAC+1)-1
  186. IF (IFIN.LE.IDEB) GOTO 300
  187. NAUX(1)=IDEB
  188. NAUX(2)=IFIN
  189. *
  190. IZ=2
  191. 208 IZ=IZ-1
  192. IF (IZ.LE.0) GOTO 209
  193. IPB=NAUX(IZ*2-1)
  194. IPH=NAUX(IZ*2)
  195. IF(IPB.GE.IPH) GOTO 208
  196. JPB=IPB-1
  197. JPH=IPH+1
  198. C CALCUL DU PIVOT
  199. NPV=0
  200. * DO 207 J=IPB,IPH
  201. * NPV=NPV+NTFAC(J)
  202. * 207 CONTINUE
  203. * NPV=NPV/(IPH-IPB+1)
  204. NPV=(NTFAC(IPB)+NTFAC(IPH))/2
  205. 242 JPB=JPB+1
  206. IF (JPH.EQ.JPB) GOTO 245
  207. IF (NTFAC(JPB).LE.NPV) GOTO 243
  208. GOTO 242
  209. 243 JPH=JPH-1
  210. IF (JPH.EQ.JPB) GOTO 245
  211. IF (NTFAC(JPH).GE.NPV) GOTO 244
  212. GOTO 243
  213. 244 IAUX=KFAK(JPB)
  214. KFAK(JPB)=KFAK(JPH)
  215. KFAK(JPH)=IAUX
  216. NTAUX=NTFAC(JPB)
  217. NTFAC(JPB)=NTFAC(JPH)
  218. NTFAC(JPH)=NTAUX
  219. GOTO 242
  220. 245 IF (JPB.GE.IPB) GOTO 247
  221. JPB=JPB+1
  222. JPH=JPH+2
  223. GOTO 248
  224. 247 IF (JPH.LE.IPH) GOTO 249
  225. JPB=JPB-2
  226. JPH=JPH-1
  227. GOTO 248
  228. 249 IF (NTFAC(JPB).LE.NPV) GOTO 250
  229. IF (JPH.EQ.IPH) GOTO 251
  230. 252 JPH=JPH+1
  231. GOTO 248
  232. 250 IF (JPB.EQ.IPB) GOTO 252
  233. 251 JPB=JPB-1
  234. 248 IF (JPB.EQ.IPB) GOTO 253
  235. NAUX(2*IZ)=JPB
  236. IZ=IZ+1
  237. 253 IF (JPH.EQ.IPH) GOTO 208
  238. NAUX(2*IZ)=IPH
  239. NAUX(2*IZ-1)=JPH
  240. IZ=IZ+1
  241. GOTO 208
  242. 209 CONTINUE
  243. 300 CONTINUE
  244.  
  245. C ======================================================================
  246. C LES FACES SONT CLASSEES DANS KFAK IL FAUT ELIMINER LES FACES EN DOUBL
  247. C ELLES SONT PAR TYPE LES UNES DERRIERES LES AUTRES LES PLUS HAUTES
  248. C DEVANT
  249. IF (IIMPI.NE.0) WRITE (IOIMP,9111) (KFAK(I),NTFAC(I),I=1,NFAC)
  250. 9111 FORMAT(5(2X,2I6))
  251. DO 400 ITYFAC=1,NTYFAC
  252. IDEB=IDXFAC(ITYFAC)
  253. IFIN=IDXFAC(ITYFAC+1)-1
  254. IF (IFIN.LE.IDEB) GOTO 400
  255. NNODE=KDFAC(1,ITYFAC)
  256. * A cette etape on doit avoir nnode.gt.0
  257. IF (NNODE.LE.0) THEN
  258. CALL ERREUR(5)
  259. RETURN
  260. ENDIF
  261. LOPT=LOQUAF.AND.(ITYFAC.EQ.7.OR.ITYFAC.EQ.8)
  262. IFACI=IPOFAC(ITYFAC)
  263. *
  264. IFINM=IFIN-1
  265. DO 450 I1=IDEB,IFINM
  266. NTI1=NTFAC(I1)
  267. IF (NTI1.EQ.0) GOTO 450
  268. IDEB1=I1+1
  269. DO 460 I2=IDEB1,IFIN
  270. NTI2=NTFAC(I2)
  271. IF (NTI2.EQ.0) GOTO 460
  272. IF (NTI2.NE.NTI1) GOTO 450
  273. IF (.NOT.LOPT) THEN
  274. IR1=KFAK(I1)
  275. IR2=KFAK(I2)
  276. DO 471 J1=1,NNODE
  277. INU=IFACI(J1,IR1)
  278. DO 472 J2=1,NNODE
  279. IF (INU.EQ.IFACI(J2,IR2)) GOTO 471
  280. 472 CONTINUE
  281. GOTO 460
  282. 471 CONTINUE
  283. ENDIF
  284. C DEUX FACES EGALES ON LES SUPPRIMENT (ENVE)
  285. if (icle.eq.0) then
  286. NTFAC(I1)=0
  287. NTFAC(I2)=0
  288. C DEUX FACES EGALES ON EN SUPPRIME UNE (CHAN FACE)
  289. elseif (icle.eq.1) then
  290. * il semble plus astucieux de supprimer la premiere pour le cas qui
  291. * ne doit pas arriver ici ou il y ait plus que deux faces egales
  292. NTFAC(I1)=0
  293. else
  294. call erreur(5)
  295. return
  296. endif
  297. GOTO 450
  298. 460 CONTINUE
  299. 450 CONTINUE
  300. 400 CONTINUE
  301. *
  302. IF (IIMPI.NE.0) WRITE (IOIMP,9111) (KFAK(I),NTFAC(I),I=1,NFAC)
  303.  
  304. SEGINI IPTFAC
  305. NBSOUS=0
  306. NBREF=0
  307. NBSOU2=0
  308. DO 600 ITYFAC=1,NTYFAC
  309. IDEB=IDXFAC(ITYFAC)
  310. IFIN=IDXFAC(ITYFAC+1)-1
  311. * write(ioimp,*) 'ityfac=',ityfac,' ideb=',ideb,' ifin=',ifin
  312. IF (IFIN.LT.IDEB) GOTO 600
  313. NNODE=KDFAC(1,ITYFAC)
  314. * A cette etape on doit avoir nnode.gt.0
  315. IF (NNODE.LE.0) THEN
  316. CALL ERREUR(5)
  317. RETURN
  318. ENDIF
  319. IFACI=IPOFAC(ITYFAC)
  320. NBELEM=0
  321. DO 611 I=IDEB,IFIN
  322. IF (NTFAC(I).NE.0) NBELEM=NBELEM+1
  323. 611 CONTINUE
  324. * write(ioimp,*) 'nbelem=',nbelem,' nnode=',nnode
  325. IF (NBELEM.EQ.0) GOTO 600
  326. NBSOU2=NBSOU2+1
  327. NBNN=NNODE
  328. SEGINI IPT1
  329. IPT1.ITYPEL=ITYEL(ITYFAC)
  330. JAUX=0
  331. DO 612 J=IDEB,IFIN
  332. IF (NTFAC(J).EQ.0) GOTO 612
  333. JAUX=JAUX+1
  334. IPT1.ICOLOR(JAUX)=IFACI(NNODE+1,KFAK(J))
  335. DO 613 I=1,NBNN
  336. IPT1.NUM(I,JAUX)=IFACI(I,KFAK(J))
  337. 613 CONTINUE
  338. 612 CONTINUE
  339. IPTFAC(ITYFAC)=IPT1
  340. * write(ioimp,*) 'ipt1=',ipt1
  341. 600 CONTINUE
  342. * on rajoute les points et les segments qui pouvaient etre dans le
  343. * maillage initial
  344. ipt5=0
  345. segact meleme
  346. ipt6=meleme
  347. do 710 io=1,max(1,lisous(/1))
  348. if (lisous(/1).ne.0) ipt6=lisous(io)
  349. segact ipt6
  350. if (ipt6.itypel.le.3) then
  351. nbsou2=nbsou2+1
  352. ipt5=ipt6
  353. endif
  354. 710 continue
  355. * write(ioimp,*) 'nbsou2=',nbsou2
  356. IF (NBSOU2.EQ.0.AND.NTPOL.EQ.0.AND.inoid.eq.0) CALL ERREUR(26)
  357. IF (NBSOU2.NE.1.OR.NTPOL.GT.0) THEN
  358. NBREF=0
  359. NBSOUS=NBSOU2+NTPOL
  360. NBNN=0
  361. NBELEM=0
  362. SEGINI IPT5
  363. I=0
  364. DO ITYFAC=1,NTYFAC
  365. IPT1=IPTFAC(ITYFAC)
  366. IF (IPT1.NE.0) THEN
  367. I=I+1
  368. IPT5.LISOUS(I)=IPT1
  369. ENDIF
  370. ENDDO
  371. segact meleme
  372. ipt1=meleme
  373. do 711 io=1,max(1,lisous(/1))
  374. if (lisous(/1).ne.0) ipt1=lisous(io)
  375. segact ipt1
  376. if (ipt1.itypel.le.3) then
  377. I=I+1
  378. IPT5.LISOUS(I)=IPT1
  379. endif
  380. 711 continue
  381. DO 720, IO = 1, NTPOL
  382. I = I+1
  383. IPT5.LISOUS(I)=IPPOL(IO)
  384. 720 CONTINUE
  385. if (ipt5.lisous(/1).eq.1) then
  386. ipt6=ipt5
  387. ipt5=ipt6.lisous(1)
  388. segsup ipt6
  389. endif
  390. ELSE
  391. * Cas ou il n'y a qu'un meleme
  392. IPTF=0
  393. DO ITYFAC=1,NTYFAC
  394. IPTF=IPTF+IPTFAC(ITYFAC)
  395. * write(ioimp,*) 'iptf=',iptf
  396. ENDDO
  397. IPTF=IPTF+IPT5
  398. IPT5=IPTF
  399. ENDIF
  400. CALL ACTOBJ('MAILLAGE',IPT5,1)
  401. CALL ECROBJ('MAILLAGE',IPT5)
  402.  
  403. SEGSUP IPTFAC,NAUX,NTFAC,KFAK,IDXFAC,IPPOL
  404. DO ITYFAC=1,NTYFAC
  405. IFACI=IPOFAC(ITYFAC)
  406. IF (IFACI.NE.0) THEN
  407. SEGSUP IFACI
  408. ENDIF
  409. ENDDO
  410. SEGSUP IPOFAC,NBFAC
  411. END
  412.  
  413.  
  414.  

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