Télécharger opto2.eso

Retour à la liste

Numérotation des lignes :

opto2
  1. C OPTO2 SOURCE SP204843 26/02/03 21:15:32 12461
  2. SUBROUTINE OPTO2(TRAVJ,JELEM,LTOPA)
  3. IMPLICIT REAL*8 (A-H,O-Z)
  4. IMPLICIT INTEGER (I-N)
  5. C***********************************************************************
  6. C NOM : OPTO2 (anciennement optt2c)
  7. C DESCRIPTION : Une implémentation de l'amélioration d'une topologie
  8. C autour d'un élément. On reprend OPTITOPO pour le corps
  9. C du programme. On reprend l'extraction et la topologie inverse de
  10. C EXTO. Le point crucial sera d'implémenter la modification de la
  11. C topologie : enlever les anciens éléments et mettre les nouveaux.
  12. C
  13. C
  14. C Ici, on est en numérotation locale et on fait l'extraction de la
  15. C topologie proprement dite, son optimisation puis sa mise à jour.
  16. C Les segments transmis sont supposés activés en *MOD
  17. C
  18. C Le point important est de construire la topologie inverse.
  19. C
  20. C Le début est identique à exto2.eso
  21. C
  22. C Repris de optt2b : on raccourcit la subroutine en externalisant
  23. C des opérations en subroutines.
  24. C
  25. C LANGAGE : ESOPE
  26. C AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SEMT/LTA)
  27. C mél : gounand@semt2.smts.cea.fr
  28. C***********************************************************************
  29. C APPELES : EXTO4C, OPTO3
  30. C APPELES (E/S) :
  31. C APPELES (BLAS) :
  32. C APPELES (CALCUL) :
  33. C APPELE PAR : OPTO1
  34. C***********************************************************************
  35. C SYNTAXE GIBIANE :
  36. C ENTREES : JELEM
  37. C ENTREES/SORTIES : JCOORD, JTOPO
  38. C SORTIES :
  39. C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
  40. C***********************************************************************
  41. C VERSION : v1, 17/10/2017, version initiale
  42. C HISTORIQUE : v1, 17/10/2017, création
  43. C HISTORIQUE :
  44. C HISTORIQUE :
  45. C***********************************************************************
  46. -INC PPARAM
  47. -INC CCOPTIO
  48. -INC TMATOP2
  49. -INC SMLENTI
  50. POINTEUR JNBL.MLENTI
  51. POINTEUR JNNO.MLENTI,KNNO.MLENTI
  52. POINTEUR NEXTO.MLENTI
  53. -INC SMELEME
  54. *
  55. * Le nombre d'éléments de JTOPO et le nombre de points de JCOORD
  56. * vont être variables. Pour ne pas avoir à ajuster ces segments en
  57. * permanence, on va dimensionner plus large, mais du coup, il faut
  58. * aussi maintenir à la main le nombre de noeuds et d'éléments
  59. * courants.
  60. *
  61. * Le nombre d'éléments courants est NVCOU et le nombre d'éléments
  62. * max est NVMAX. Idem pour le nombre de noeuds courants et max :
  63. * NPCOU et NPMAX.
  64. *
  65. * Numerotation locale des éléments JTOPO.NUM(NBNN,NBELEM)
  66. * INTEGER NVCOU,NVMAX
  67. POINTEUR JTOPO.MELEME
  68. *del POINTEUR KTOPO.MELEME
  69. POINTEUR JELEM.MELEME
  70. POINTEUR JELEM1.MELEME
  71. POINTEUR JEXTO.MELEME,KEXTO.MELEME
  72. POINTEUR JTBES.MELEME
  73. -INC SMLOBJE
  74. POINTEUR LTOPA.MLOBJE
  75. -INC SMCOORD
  76. * Numerotation locale (de 1 à NBPTS)
  77. * INTEGER NPCOU,NPMAX
  78. *del POINTEUR JCOORD.MCOORD
  79. POINTEUR KCOORD.MCOORD
  80. -INC TMATOP1
  81. *-INC STOPINV
  82. *-INC SMETRIQ
  83. POINTEUR JCMETR.METRIQ
  84. POINTEUR KCMETR.METRIQ
  85. *-INC STRAVJ
  86. POINTEUR TRAVK.TRAVJ
  87. *-INC STRAVL
  88. -INC SMLMOTS
  89. POINTEUR JNMETR.MLMOTS
  90. POINTEUR KNMETR.MLMOTS
  91. *
  92. POINTEUR TRAVX.TRAVJ
  93. *
  94. *-INC SMLENTX
  95. POINTEUR ICPRX.MLENTX
  96. POINTEUR IDCPX.MLENTX
  97. POINTEUR ICPR2.MLENTX
  98. POINTEUR IDCP2.MLENTX
  99. *-INC SMELEMX
  100. POINTEUR KELEMX.MELEMX
  101. POINTEUR JELEM2.MELEMX
  102. *
  103. logical lchang
  104. *
  105. * Executable statements
  106. *
  107. if (impr.ge.4) WRITE(IOIMP,*) 'Entrée dans opto2.eso'
  108. *
  109. * Initialisation et extension des segments JTOPO et JCOORD
  110. *
  111. IDIMP=IDIM+1
  112.  
  113. JTOPO=TRAVJ.TOPO
  114. *
  115. * Initialisation de la topologie inverse
  116. *
  117. * CALL INTOPI(NVMAX,NPMAX,TOPINV,IMPR)
  118. CALL INTOP2(TRAVJ,IMPR)
  119. IF (IERR.NE.0) RETURN
  120. *
  121. * Remplissage de la topologie inverse avec JTOPO
  122. *
  123. * CALL RETOPI(JTOPO,NVCOU,TOPINV,IMPR)
  124. CALL RETOP2(TRAVJ,IMPR)
  125. IF (IERR.NE.0) RETURN
  126. *
  127. if (.false.) then
  128. TOPINV=TRAVJ.TOPI
  129. call ectopi(TOPINV,1)
  130. call ectopi(TOPINV,2)
  131. endif
  132. *
  133. * Initialisation LTOPA
  134. *
  135. NOBJ=0
  136. NREE=0
  137. SEGINI LTOPA
  138. LTOPA.TYPOBJ='MAILLAGE'
  139. *
  140. * Segment de travail pour l'extraction des éléments
  141. *
  142. JG=NVMAX
  143. SEGINI JNBL
  144. TRAVJ.NBL=JNBL
  145. JG=NPMAX-NPINI
  146. SEGINI JNNO
  147. TRAVJ.NNO=JNNO
  148. *
  149. * Extraction de la topologie à optimiser
  150. *
  151. * NELMOY=40
  152. IF (IDIM.EQ.2) THEN
  153. NELMOY=15
  154. NPOMOY=10
  155. ELSEIF (IDIM.EQ.3) THEN
  156. NELMOY=40
  157. NPOMOY=12
  158. * NELMOY=40
  159. * NPOMOY=25
  160. ELSE
  161. write(ioimp,*) 'idim=',idim
  162. goto 9999
  163. ENDIF
  164. *
  165. *!!! A changer plus tard
  166. *
  167. * NVXMAX=0
  168. SEGINI TRAVX
  169. *old if (isgadj.gt.0) write(ioimp,185) 'TRAVJ,TRAVX=',TRAVJ,TRAVX
  170. TRAVX.NVINI=0
  171. TRAVX.NVCOU=0
  172. TRAVX.NVMAX=NELMOY
  173. * TRAVX.NVMAX=0
  174. JG=TRAVX.NVMAX
  175. SEGINI NEXTO
  176. TRAVX.NBL=NEXTO
  177. *
  178.  
  179. NBNN=IDIMP
  180. NBELEM=TRAVX.NVMAX
  181. NBSOUS=0
  182. NBREF=0
  183. SEGINI JEXTO
  184. JEXTO.ITYPEL=JTOPO.ITYPEL
  185. TRAVX.TOPO=JEXTO
  186. * Boucle sur les éléments
  187. * write(ioimp,*) 'opto2 jelem(nno,nbnn=)',jelem.num(/1)
  188. * $ ,jelem.num(/2)
  189. NBNN1=JELEM.NUM(/1)
  190. NBNN=NBNN1
  191. NBELEM=1
  192. NBSOUS=0
  193. NBREF=0
  194. SEGINI JELEM1
  195. JELEM1.ITYPEL=JELEM.ITYPEL
  196.  
  197. * Segment de travail TRAVK pour opto3 numérotation locale à
  198. * l'élément extrait.
  199. SEGINI TRAVK
  200. TRAVK.NVINI=0
  201. TRAVK.NVCOU=0
  202. TRAVK.NVMAX=NELMOY
  203. * Important pour le segment NNO après
  204. TRAVK.NPINI=0
  205. TRAVK.NPCOU=0
  206. TRAVK.NPMAX=NPOMOY
  207. * Topologie de TRAVK (KEXTO)
  208. NBELEM=TRAVK.NVMAX
  209. NBNN=IDIMP
  210. NBSOUS=0
  211. NBREF=0
  212. SEGINI,KEXTO
  213. KEXTO.ITYPEL=JEXTO.ITYPEL
  214. TRAVK.TOPO=KEXTO
  215. * Coordonnées de TRAVK (KCOORD)
  216. NBPTS=TRAVK.NPMAX
  217. SEGINI,KCOORD
  218. TRAVK.COORD=KCOORD
  219. JNMETR=TRAVJ.NMETR
  220. IF (JNMETR.NE.0) THEN
  221. SEGINI,KNMETR=JNMETR
  222. TRAVK.NMETR=KNMETR
  223. ENDIF
  224. JCMETR=TRAVJ.CMETR
  225. IF (JCMETR.NE.0) THEN
  226. NNIN=JCMETR.XIN(/1)
  227. NNNOE=TRAVK.NPMAX
  228. SEGINI,KCMETR
  229. TRAVK.CMETR=KCMETR
  230. ENDIF
  231. *
  232. * Segment de travail pour trouver les noeuds du contour ou de
  233. * l'enveloppe pour étoiler dans topv2
  234. *
  235. JG=TRAVK.NPMAX-TRAVK.NPINI
  236. SEGINI KNNO
  237. TRAVK.NNO=KNNO
  238. * Segment de travail TRAVL pour topv2
  239. NNM=JEXTO.NUM(/1)
  240. ITYP=JEXTO.ITYPEL
  241. *del CALL TRLINI(NELMOY,JEXTO.NUM(/1),JEXTO.ITYPEL,TRAVL)
  242. CALL TRLINI(NELMOY,NNM,ITYP,TRAVL)
  243. if (iveri.ge.2) then
  244. call trlver(travl,'opto2 : Apres initialisation TRAVL')
  245. if (ierr.ne.0) return
  246. endif
  247. *
  248. * Segment de travail pour le changement de numérotation dans opto3
  249. *
  250. JGMAX=NPOMOY
  251. SEGINI ICPRX
  252. CALL mtxadj(ICPRX,0,lchang,'opto2 : ICPRX_INI')
  253. if (ierr.ne.0) return
  254. SEGINI IDCPX
  255. CALL mtxadj(IDCPX,0,lchang,'opto2 : IDCPX_INI')
  256. if (ierr.ne.0) return
  257. *
  258. * Segment de travail pour jelem en numérotation très locale dans
  259. * opto3. Ce segment a un élément et peut-être moins de noeuds que JELEM1
  260. *
  261. NNMAX=JELEM.NUM(/1)
  262. NLMAX=1
  263. SEGINI KELEMX
  264. KELEMX.ITYPEX=JELEM.ITYPEL
  265. KELEMX.NLCOU=1
  266. *
  267. * Segments de travail pour le voisinage des elements
  268. *
  269. JGMAX=NPOMOY
  270. SEGINI ICPR2
  271. CALL mtxadj(ICPR2,0,lchang,'opto2 : ICPR2_INI')
  272. if (ierr.ne.0) return
  273. SEGINI IDCP2
  274. CALL mtxadj(IDCP2,0,lchang,'opto2 : IDCP2_INI')
  275. if (ierr.ne.0) return
  276. NNMAX=1
  277. NLMAX=NPOMOY
  278. SEGINI JELEM2
  279. JELEM2.ITYPEX=JELEM.ITYPEL
  280. CALL mlxadl(JELEM2,0,lchang,'opto2 : JELEM2_INI')
  281. if (ierr.ne.0) return
  282.  
  283. NBNNT=JTOPO.NUM(/1)
  284. TOPINV=TRAVJ.TOPI
  285. *
  286. DO IAPARC=1,JELEM.NUM(/2)
  287. DO IBNN=1,NBNN1
  288. JELEM1.NUM(IBNN,1)=JELEM.NUM(IBNN,IAPARC)
  289. ENDDO
  290. JPARCO=JPARCO+1
  291. IF (IMPR.GE.4) THEN
  292. write(ioimp,*) ' opto2 : Autour de l''element ',iaparc
  293. call ecmai1(jelem1,0)
  294. segact jelem1*mod
  295. ENDIF
  296. * On fait le voisinage de Gruau uniquement pour les noeuds en 2D
  297. * et les noeuds et les aretes un 3D.
  298. * Sinon, on fait l'ancienne manière. C'est un test parce que la
  299. * méthode Gruau n'est pas terrible dans les coins.
  300. * IF (NBNN1.LT.IDIM+2) THEN
  301. IF (.true.) THEN
  302. * Voisinages des noeuds de JELEM1
  303. NCOU=TRAVJ.NPCOU
  304. CALL mtxadj(ICPR2,NCOU,lchang,'opto2 : ICPR2_dim')
  305. if (ierr.ne.0) return
  306. IIDCP2=0
  307. * Voisinage pour chaque noeud
  308. DO IBNN=1,NBNN1
  309. INOD=JELEM1.NUM(IBNN,1)
  310. IVAL=2**(IBNN-1)
  311. * Parcours de la topologie inverse et cochage dans ICPR2
  312. LAST=TIC(INOD)
  313. LDG=TDC(INOD)
  314. DO IDG=1,LDG
  315. IL=((LAST-1)/IDIMP)+1
  316. DO IBNNT=1,NBNNT
  317. IN=JTOPO.NUM(IBNNT,IL)
  318. IF (IN.NE.INOD) THEN
  319. ITOUCH=ICPR2.LECTX(IN)
  320. IF (ITOUCH.LT.IVAL) THEN
  321. ICPR2.LECTX(IN)=ITOUCH+IVAL
  322. IF (ITOUCH.EQ.0) THEN
  323. IIDCP2=IIDCP2+1
  324. CALL mtxadj(IDCP2,IIDCP2,lchang,'opto2 :
  325. $ IDCP2_dim')
  326. if (ierr.ne.0) return
  327. IDCP2.LECTX(IIDCP2)=IN
  328. ENDIF
  329. ENDIF
  330. ENDIF
  331. ENDDO
  332. LAST=TLC(LAST)
  333. ENDDO
  334. ENDDO
  335. NIDCP2=IIDCP2
  336. * JELEM2 a au moins les noeuds de JELEM1
  337. NBNN=NBNN1
  338. * segprt,icpr2
  339. * Intersection des voisinages
  340. IVAL=(2**NBNN1)-1
  341. NIDCP2=IIDCP2
  342. DO IIDCP2=1,NIDCP2
  343. IN=IDCP2.LECTX(IIDCP2)
  344. IF (ICPR2.LECTX(IN).eq.IVAL) THEN
  345. NBNN=NBNN+1
  346. ENDIF
  347. ENDDO
  348. *
  349. call mlxadl(JELEM2,NBNN,lchang,'opto2 : JELEM2_NBNN')
  350. if (ierr.ne.0) return
  351. if (iveri.ge.2) then
  352. call vemelx(jelem2,'opto2 : Apres requisition jelem2')
  353. if (ierr.ne.0) return
  354. endif
  355.  
  356. ibnn=0
  357. do i=1,nbnn1
  358. ibnn=ibnn+1
  359. jelem2.numx(1,ibnn)=jelem1.num(ibnn,1)
  360. enddo
  361. DO IIDCP2=1,NIDCP2
  362. IN=IDCP2.LECTX(IIDCP2)
  363. IF (ICPR2.LECTX(IN).eq.IVAL) THEN
  364. ibnn=ibnn+1
  365. jelem2.numx(1,ibnn)=IN
  366. ENDIF
  367. icpr2.lectx(IN)=0
  368. ENDDO
  369. IF (IMPR.GE.4) THEN
  370. write(ioimp,*)
  371. $ ' opto2 : Noeuds du voisinage de l''element '
  372. $ ,iaparc
  373. call ecmelx(jelem2,0)
  374. ENDIF
  375. *
  376. if (iveri.ge.2) call vetopi(travj,'Avant exto5c')
  377. if (ierr.ne.0) return
  378. *
  379. CALL EXTO5c(JELEM2,TRAVJ,
  380. $ TRAVX)
  381. * Menage JELEM2
  382. if (iveri.ge.2) then
  383. DO IZER=1,JELEM2.NLCOU
  384. JELEM2.NUMX(1,IZER)=0
  385. ENDDO
  386. endif
  387. CALL mlxadl(JELEM2,0,lchang,'opto2 : JELEM2_0')
  388. if (ierr.ne.0) return
  389. if (iveri.ge.2) then
  390. call vemelx(jelem2,'opto2 : Apres nettoyage jelem2')
  391. if (ierr.ne.0) return
  392. endif
  393. * verif que NBL est bien nettoyé
  394. if (iveri.ge.2) call vetopi(travj,'Apres exto5c')
  395. if (ierr.ne.0) return
  396. ELSE
  397. if (iveri.ge.2) call vetopi(travj,'Avant exto4')
  398. if (ierr.ne.0) return
  399. *
  400. CALL EXTO4c(JELEM1,TRAVJ,
  401. $ TRAVX)
  402. * verif que NBL est bien nettoyé
  403. if (iveri.ge.2) call vetopi(travj,'Apres exto4')
  404. if (ierr.ne.0) return
  405. ENDIF
  406. *tst write(ioimp,*) 'Elements de la topologie extraits :'
  407. *tst write(ioimp,187) (nexto(I),I=1,travx.nvcou)
  408. * Mise à jour de jexto
  409. nexto=travx.nbl
  410. jexto=travx.topo
  411. do iel=1,travx.nvcou
  412. do ino=1,IDIMP
  413. JEXTO.NUM(ino,iel)=JTOPO.NUM(INO,nexto.lect(iel))
  414. enddo
  415. enddo
  416. *
  417. * Optimisation de la topologie extraite
  418. *
  419. IF (IMPR.GE.4) THEN
  420. write(ioimp,*) 'opto2.eso : on a extrait la topologie : '
  421. call ecmai1(jexto,0)
  422. segact jexto
  423. ENDIF
  424. *
  425. * Init
  426. *
  427. CALL OPTO3(TRAVJ,TRAVX,JELEM1,TRAVK,TRAVL,ICPRX,IDCPX,
  428. $ KELEMX,
  429. $ JTBES,JCAND)
  430. IF (IERR.NE.0) RETURN
  431. if (iveri.ge.2) call vetopi(travj,'Apres opto3')
  432. IF (IERR.NE.0) RETURN
  433.  
  434. JEXPLO=JEXPLO+ABS(JCAND)
  435. IF (IMPR.GE.4) THEN
  436. IF (JEXTO.EQ.JTBES) THEN
  437. WRITE(IOIMP,*) 'Pas damelioration JTBES=',JTBES
  438. ELSE
  439. WRITE(IOIMP,*) 'Topologie amelioree JTBES=',JTBES
  440. CALL ECMAI1(JTBES,0)
  441. segact jtbes
  442. ENDIF
  443. ENDIF
  444. *
  445. * Si la topologie locale a été améliorée, on change la topologie
  446. * globale en conséquence
  447. *
  448. IF (JEXTO.NE.JTBES) THEN
  449. JCHANG=JCHANG+1
  450. * CALL TOPDIF(TRAVJ,TRAVX)
  451. CALL TOPDI2(TRAVJ,TRAVX)
  452. if (ierr.ne.0) return
  453. if (iveri.ge.2) call vetopi(travj,'Apres DIFF')
  454. if (ierr.ne.0) return
  455. *
  456. * On ajoute les éléments de JTBES dans JTOPO
  457. *
  458. CALL TOPFUS(TRAVJ,JTBES)
  459. if (ierr.ne.0) return
  460. if (iveri.ge.2) call vetopi(travj,'Apres ET')
  461. if (ierr.ne.0) return
  462. if (iseqm.ne.0) then
  463. NOBJ=NOBJ+3
  464. SEGADJ LTOPA
  465. * WRITE(IOIMP,*) 'NOBJ=',NOBJ
  466. * Ajouter JELEM1
  467. NBNN=JELEM1.NUM(/1)
  468. NBELEM=1
  469. NBSOUS=0
  470. NBREF=0
  471. SEGINI MELEME
  472. ITYPEL=JELEM1.ITYPEL
  473. DO IBNN=1,NBNN
  474. NUM(IBNN,1)=JELEM1.NUM(IBNN,1)
  475. ENDDO
  476. LTOPA.LISOBJ(NOBJ-2)=MELEME
  477. * write(ioimp,*) 'jelem1'
  478. * call ecmai1(meleme,0)
  479. * segact meleme*mod
  480. *
  481. * Ajouter JEXTO reduit
  482. NBNN=JEXTO.NUM(/1)
  483. NBELEM=TRAVX.NVCOU
  484. NBSOUS=0
  485. NBREF=0
  486. SEGINI MELEME
  487. ITYPEL=JEXTO.ITYPEL
  488. DO IBELEM=1,NBELEM
  489. DO IBNN=1,NBNN
  490. NUM(IBNN,IBELEM)=JEXTO.NUM(IBNN,IBELEM)
  491. ENDDO
  492. ENDDO
  493. LTOPA.LISOBJ(NOBJ-1)=MELEME
  494. * write(ioimp,*) 'jexto'
  495. * call ecmai1(meleme,0)
  496. * segact meleme*mod
  497. * Ajouter JTBES
  498. LTOPA.LISOBJ(NOBJ)=JTBES
  499. * write(ioimp,*) 'jtbes'
  500. * call ecmai1(jtbes,0)
  501. * segact jtbes*mod
  502. endif
  503. ENDIF
  504. *
  505. * Nettoyage de NEXTO et JEXTO (normalement inutile mais utilisé pour
  506. * vetopi)
  507. *
  508. if (iveri.ge.1) then
  509. nexto=travx.nbl
  510. jexto=travx.topo
  511. do iel=1,travx.nvcou
  512. do ino=1,IDIMP
  513. JEXTO.NUM(ino,iel)=0
  514. nexto.lect(iel)=0
  515. enddo
  516. enddo
  517. travx.nvcou=0
  518. endif
  519.  
  520. * Fin boucle sur les éléments
  521. ENDDO
  522. *
  523. * Il faut appeler le nettoyage avant de sortir
  524. *
  525. SEGSUP JELEM2
  526. SEGSUP IDCP2
  527. SEGSUP ICPR2
  528. SEGSUP KELEMX
  529. SEGSUP ICPRX
  530. SEGSUP IDCPX
  531. CALL TRLSUP(TRAVL)
  532. * SEGSUP IPBTL
  533. CALL TOPSUP(TRAVK)
  534. *tst topinv=travj.topi
  535. *tst write(ioimp,*) 'TOPINV Avant nettoyage elem TOPINV'
  536. *tst call ectopi(topinv,1)
  537. *tst call ectopi(topinv,2)
  538. segsup jelem1
  539. if (jtbes.ne.travx.topo.and.iseqm.eq.0) segsup jtbes
  540. call topsup(travx)
  541.  
  542. * Nettoyage des éléments vides
  543. * impr=8
  544. call topclv(travj,lchang)
  545. if (ierr.ne.0) return
  546. if (iveri.ge.2.and.lchang) call vetopi(travj
  547. $ ,'Apres nettoyage elem')
  548. if (ierr.ne.0) return
  549. *
  550. * Nettoyage des noeuds non référencés dans la topologie mais
  551. * seulement ceux ajoutés par nous, pas les autres !
  552. *
  553. IF (iseqm.eq.0) call topclp(travj,lchang)
  554. * verif
  555. if (iveri.ge.2.and.lchang) call vetopi(travj
  556. $ ,'Apres nettoyage noeuds')
  557. if (ierr.ne.0) return
  558.  
  559. *
  560. * Normal termination
  561. *
  562. RETURN
  563. *
  564. * Error handling
  565. *
  566. 9999 CONTINUE
  567. MOTERR(1:8)='OPTO2 '
  568. * 349 2
  569. *Problème non prévu dans le s.p. %m1:8 contactez votre assistance
  570. CALL ERREUR(349)
  571. RETURN
  572. *
  573. * End of subroutine OPTO2
  574. *
  575. END
  576.  
  577.  
  578.  

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