Télécharger numopt.eso

Retour à la liste

Numérotation des lignes :

  1. C NUMOPT SOURCE PV 19/04/19 21:15:01 10203
  2. C RACINE DE LA NUMEROTATION POUR LA SORTIE SUR FAC
  3. C
  4. SUBROUTINE NUMOPT(MELEME,ICPR,NODES)
  5. IMPLICIT INTEGER(I-N)
  6. -INC SMELEME
  7. -INC SMCOORD
  8. -INC CCOPTIO
  9. -INC CCPRECO
  10. CHARACTER*4 MVAL
  11. DATA MVAL/'NOOP'/
  12. SEGMENT MEMJT1(NKON)
  13. SEGMENT MEMJT2(NKON)
  14. SEGMENT IPOME(NODES+1)
  15. SEGMENT JNT(NODES)
  16. SEGMENT JMEM(NODES)
  17. SEGMENT ICPR(XCOOR(/1)/(IDIM+1))
  18. SEGMENT ICPREN(XCOOR(/1)/(IDIM+1))
  19. SEGMENT IDCP(XCOOR(/1)/(IDIM+1))
  20. SEGMENT JOINT(NODES)
  21. SEGMENT NEWJT(NODES)
  22. icpren=0
  23. segact mcoord
  24. segact meleme
  25. * cas trivial ou le maillage est fait de points isoles
  26. CALL LIRMOT(MVAL,1,IRET,0)
  27. if (itypel.eq.1.or.iret.eq.1.or.nucrou.eq.-1) then
  28. segact icpr*MOD
  29. ia=0
  30. ipt1=meleme
  31. segact meleme
  32. do k=1,max(1,lisous(/1))
  33. if (lisous(/1).ne.0) ipt1=lisous(k)
  34. segact ipt1
  35. do j= ipt1.num(/2),1,-1
  36. do i=1,ipt1.num(/1)
  37. if (icpr(ipt1.num(i,j)).eq.0) then
  38. ia=ia+1
  39. icpr(ipt1.num(i,j))=ia
  40. endif
  41. enddo
  42. enddo
  43. enddo
  44. nodes=ia
  45. return
  46. endif
  47. nucro=0
  48. * si plus de 4000 noeuds dans le meleme numérotation NS
  49. segact icpr*mod
  50. do i=1,icpr(/1)
  51. icpr(i)=0
  52. enddo
  53. segact meleme
  54. ia=0
  55. ib=0
  56. ipt1=meleme
  57. ifrot=0
  58. do 60 isous=1,max(1,lisous(/1))
  59. if (lisous(/1).ne.0) ipt1=lisous(isous)
  60. segact ipt1
  61. if (ipt1.itypel.eq.-22) then
  62. * write (6,*) ' frottement '
  63. ifrot=1
  64. endif
  65. ib=max(ib,ipt1.num(/1))
  66. do 50 i=1,ipt1.num(/2)
  67. do 50 j=1,ipt1.num(/1)
  68. if (icpr(ipt1.num(j,i)).eq.0) then
  69. ia=ia+1
  70. icpr(ipt1.num(j,i))=ia
  71. endif
  72. 50 continue
  73. 60 continue
  74. if (ia.lt.ib*5 ) then
  75. * write (6,*) ' superelement detecte ',ia,ib
  76. ia=1
  77. endif
  78. do i=1,icpr(/1)
  79. icpr(i)=0
  80. enddo
  81. * choix de la numerotation en fonction du nombre de noeuds
  82. if (ia.gt.4000) nucro=2
  83. if (ia.le.4000) nucro=0
  84. * Reverse CM en iteratif
  85. IF(NUCROU.EQ.1) nucro=0
  86. IF(NUCROU.EQ.2) nucro=2
  87. * en cas de frottement forcément nested dissection pour placer correctement les mult de frot
  88. if (ifrot.eq.1) nucro=2
  89. * ajuster le nombre de zero consecutif
  90. if (nucro.eq.0) izrosf=48
  91. if (nucro.eq.2) izrosf=16
  92. CALL LIRMOT(MVAL,1,IRET,0)
  93. if (iret.eq.0) then
  94. * a t'on deja fait l'operation ???
  95. segini icpren,idcp
  96. incren=0
  97. do isous=2,lisous(/1)
  98. ipt3=lisous(isous)
  99. segact ipt3
  100. * SG 2011/10/07 on a échangé l'ordre des boucles pour le remettre dans
  101. * l'ordre "usuel"
  102. do j=1,ipt3.num(/2)
  103. do i=1,ipt3.num(/1)
  104. ip=ipt3.num(i,j)
  105. if (icpren(ip).eq.0) then
  106. incren=incren+1
  107. icpren(ip)=incren
  108. endif
  109. enddo
  110. enddo
  111. enddo
  112. do i=1,icpren(/1)
  113. if (icpren(i).ne.0) idcp(icpren(i))=i
  114. enddo
  115.  
  116. do 150 ip=1,30
  117. if (prenum(ip).ne.0) then
  118. ipt1=prenum(ip)
  119. segact ipt1
  120. * write (6,*) ' prenum ',ip,lisous(/1),ipt1.lisous(/1)
  121. if (lisous(/1)+1.eq.ipt1.lisous(/1)) then
  122. *SG Attention ! On a osé stocker le type de renumérotation qui avait été
  123. *SG calculée dans le itypel
  124. ipt2=ipt1.lisous(ipt1.lisous(/1))
  125. segact ipt2
  126. nucroa=ipt2.itypel
  127. segdes ipt2
  128. if (nucro.ne.nucroa) goto 11
  129. * verif maillage identique
  130. * write (6,*) ' maillage 1 ',lisous(/1)
  131. * write (6,*) ' maillage 2 ',ipt1.lisous(/1)
  132. do isous=2,lisous(/1)
  133. ipt3=lisous(isous)
  134. ipt4=ipt1.lisous(isous)
  135. segact ipt3,ipt4
  136. * write (6,*) ipt3.num(/1),ipt4.num(/1),ipt3.num(/2),ipt4.num(/2)
  137. * write (6,*) ' itypel ',ipt3.itypel,ipt4.itypel
  138. if (ipt3.itypel.ne.ipt4.itypel) goto 10
  139. if (ipt3.num(/1).ne.ipt4.num(/1)) goto 10
  140. if (ipt3.num(/2).ne.ipt4.num(/2)) goto 10
  141. * SG 2011/10/07 on a échangé l'ordre des boucles pour le remettre dans
  142. * l'ordre "usuel"
  143. do j=1,ipt3.num(/2)
  144. do i=1,ipt3.num(/1)
  145. * write (6,*) i,j,ipt3.num(i,j),
  146. * > icpren(ipt3.num(i,j)),ipt4.num(i,j)
  147. if (icpren(ipt3.num(i,j)).ne.ipt4.num(i,j)) goto 10
  148. enddo
  149. enddo
  150. segdes ipt3,ipt4
  151. enddo
  152. if (iimpi.ne.0) write (6,*)
  153. > ' preconditionnement de la numerotation ',ip,
  154. > (prenum(i),i=1,30)
  155. ipt2=ipt1.lisous(ipt1.lisous(/1))
  156. segact ipt2
  157. segact icpr*MOD
  158. nodes=ipt2.num(/2)
  159. * write (6,*) ' icpr(/1),nodes ',icpr(/1),nodes
  160. do i=1,nodes
  161. if (idcp(ipt2.num(1,i)).ne.0) icpr(idcp(ipt2.num(1,i)))=i
  162. enddo
  163. ia=nodes
  164. segdes ipt1,ipt2
  165. * write (6,*) ' apres precond dans numopt '
  166. * write (6,*) (icpr(i),i=1,icpr(/1))
  167. return
  168.  
  169. 10 continue
  170. segdes ipt3,ipt4
  171. 11 continue
  172. endif
  173. endif
  174. 150 continue
  175. if (nucro.eq.1) then
  176. call numop1(MELEME,ICPR,NODES)
  177. goto 300
  178. endif
  179. if (nucro.eq.2) then
  180. call numop2(MELEME,ICPR,NODES)
  181. segact icpr*mod
  182. goto 300
  183. endif
  184. endif
  185. if (iret.eq.1) call refus
  186. IENORM=2000000000
  187. NODES=XCOOR(/1)/(IDIM+1)
  188. SEGACT ICPR*MOD
  189. SEGACT MELEME
  190. DO 200 I=1,ICPR(/1)
  191. 200 ICPR(I)=0
  192. IPT1=MELEME
  193. IKOU=0
  194. DO 202 IO=1,MAX(1,LISOUS(/1))
  195. IF (LISOUS(/1).NE.0) THEN
  196. IPT1=LISOUS(IO)
  197. SEGACT IPT1
  198. ENDIF
  199. DO 203 I=1,IPT1.NUM(/1)
  200. DO 203 J=1,IPT1.NUM(/2)
  201. IJ=IPT1.NUM(I,J)
  202. IF (ICPR(IJ).NE.0) GOTO 203
  203. IKOU=IKOU+1
  204. ICPR(IJ)=IKOU
  205. 203 CONTINUE
  206. 202 CONTINUE
  207. NODES=IKOU
  208. SEGINI JNT,JMEM,NEWJT,JOINT,IPOME
  209. DO 1 I=1,NODES
  210. 1 JMEM(I)=0
  211. IPT1=MELEME
  212. DO 3 IO=1,MAX(1,LISOUS(/1))
  213. IF (LISOUS(/1).NE.0) IPT1=LISOUS(IO)
  214. DO 4 I=1,IPT1.NUM(/1)
  215. DO 4 J=1,IPT1.NUM(/2)
  216. JMEM(ICPR(IPT1.NUM(I,J)))=JMEM(ICPR(IPT1.NUM(I,J)))+1
  217. 4 CONTINUE
  218. 3 CONTINUE
  219. IPOME(1)=0
  220. DO 6 I=1,NODES
  221. IPOME(I+1)=IPOME (I) + JMEM(I)
  222. 6 CONTINUE
  223. NKON=IPOME(NODES+1)
  224. SEGINI MEMJT1,memjt2
  225. DO 90 J=1,NODES
  226. 90 JMEM(J)=0
  227. IPT1=MELEME
  228. DO 101 IO=1,MAX(1,LISOUS(/1))
  229. IF (LISOUS(/1).NE.0) IPT1=LISOUS(IO)
  230. DO 100 I=1,IPT1.NUM(/2)
  231. DO 100 J=1,IPT1.NUM(/1)
  232. IND=ICPR(IPT1.NUM(J,I))
  233. JMEM(IND)=JMEM(IND)+1
  234. MEMJT1(IPOME(IND)+JMEM(IND))=I
  235. MEMJT2(IPOME(IND)+JMEM(IND))=IO
  236. 100 CONTINUE
  237. 101 CONTINUE
  238. IDIFF=1
  239. CALL OPTNUM(MELEME,MEMJT1,memjt2,JMEM,JNT,NEWJT,JOINT,IDIFF,IPOME,
  240. # ICPR,NODES,NOOPTI)
  241. DO 110 I=1,NODES
  242. JNT(I)=NODES-JNT(I)+1
  243. 110 CONTINUE
  244. C PERMUTER LES COORDONNEES ET CORRIGER IPAD
  245. DO 111 I=1,icpr(/1)
  246. IF (ICPR(I).EQ.0) GOTO 111
  247. ICPR(I)=JNT(ICPR(I))
  248. 111 CONTINUE
  249. IF(NOOPTI .EQ.0) GO TO 116
  250. IA=0
  251. DO 115 I=1,icpr(/1)
  252. IF(ICPR(I).EQ.0) GO TO 115
  253. IA=IA+1
  254. ICPR(I)=IA
  255. 115 CONTINUE
  256. 116 CONTINUE
  257. SEGSUP MEMJT1,MEMJT2,JMEM,NEWJT,JOINT,IPOME
  258. SEGSUP JNT
  259. * Mise en commentaire du retuen par SG car comme cela on préconditionne
  260. * aussi le Cuthill-McKee
  261. * RETURN
  262. 300 continue
  263. * sauver le resultat au cas ou on veuille recommencer
  264. segini,ipt1=meleme
  265. nbsous=ipt1.lisous(/1)+1
  266. nbref=ipt1.lisref(/1)
  267. nbnn=ipt1.num(/1)
  268. nbelem=ipt1.num(/2)
  269. segadj ipt1
  270. nbnn=1
  271. nbsous=0
  272. nbref=0
  273. nbelem=nodes
  274. segini ipt2
  275. *SG Attention ! On ose stocker le type de renumérotation qui a été
  276. *SG calculée dans le itypel
  277. ipt2.itypel=nucro
  278. ia=0
  279. do 120 i=1,icpr(/1)
  280. if (icpr(i).gt.nodes.or.icpr(i).eq.0) goto 120
  281. ia=ia+1
  282. ipt2.num(1,icpr(i))=icpren(i)
  283. 120 continue
  284. * call ecmail(ipt2)
  285. ipt1.lisous(ipt1.lisous(/1))=ipt2
  286. segdes ipt2
  287. do 125 is=1,ipt1.lisous(/1)-1
  288. ipt2=ipt1.lisous(is)
  289. segini,ipt3=ipt2
  290. segdes ipt2
  291. * SG 2011/10/07 on a échangé l'ordre des boucles pour le remettre dans
  292. * l'ordre "usuel"
  293. do il=1,ipt3.num(/2)
  294. do ip=1,ipt3.num(/1)
  295. ipt3.num(ip,il)=icpren(ipt3.num(ip,il))
  296. enddo
  297. enddo
  298. segdes ipt3
  299. ipt1.lisous(is)=ipt3
  300. 125 continue
  301. if (ipt1.lisous(/1).eq.1) goto 131
  302. segdes ipt1
  303. do 130 ip=29,1,-1
  304. prenum(ip+1)=prenum(ip)
  305. 130 continue
  306. prenum(1)=ipt1
  307. 131 continue
  308. if (icpren.ne.0) segsup icpren,idcp
  309. segact icpr*mod
  310. * write (6,*) ' fin de numopt '
  311. * write (6,*) (icpr(i),i=1,icpr(/1))
  312. *
  313. RETURN
  314. END
  315.  
  316.  
  317.  
  318.  
  319.  
  320.  
  321.  
  322.  
  323.  
  324.  
  325.  
  326.  
  327.  
  328.  
  329.  
  330.  
  331.  

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