Télécharger la plaquette Cast3M


fr en

Télécharger topv2.eso

Retour à la liste

Numérotation des lignes : 

topv2

  1. C TOPV2     SOURCE    GOUNAND   21/04/06    21:15:33     10940          
  2. *ijob      SUBROUTINE TOPV2(TRAVK,KELEM,IJOB,XVTOL,QTOL,IMET,XDENS,
  3. *kelemx      SUBROUTINE TOPV2(TRAVK,KELEM,IALGO,IAJNO,XVTOL,QTOL,IMET,XDENS,
  4.       SUBROUTINE TOPV2(TRAVK,KELEMX,IALGO,IAJNO,XVTOL,QTOL,IMET,IMOMET
  5.      $     ,XDENS,INCMA,ISTMA,KTBES,JCAND,JNASCM,IVERI,impr,TRAVL
  6.      $     ,lchtop)
  7.       IMPLICIT REAL*8 (A-H,O-Z)
  8.       IMPLICIT INTEGER (I-N)
  9. C***********************************************************************
  10. C NOM         : TOPV2
  11. C DESCRIPTION : IJOB=0
  12. C               Minimise le volume d'une topologie de maillage
  13. C               en le maintenant supérieur à 0
  14. C               IJOB=1
  15. C               Minimise le volume, mais on a le droit d'ajouter des
  16. C               noeuds internes
  17. C               IJOB=2
  18. C               La topologie de maillage est supposée être un maillage
  19. C               On essaie de l'améliorer en conservant son volume
  20. C               mais en augmentant sa qualité grace a l'adjonction
  21. C               de noeuds internes
  22. C     
  23. * 2017/11/30 : On remplace par ialgo (0 ou 1 : génération ou
  24. *              optimisation de maillage) et iajno (autorise-t-on
  25. *              l'algorithme à ajouter des noeuds.)
  26. C     
  27. C     
  28. C LANGAGE     : ESOPE
  29. C AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  30. C               mél : gounand@semt2.smts.cea.fr
  31. C***********************************************************************
  32. C VERSION    : v1, 05/05/2013, version initiale
  33. C HISTORIQUE : v1, 05/05/2013, création
  34. C HISTORIQUE :
  35. C HISTORIQUE :
  36. C***********************************************************************
  37. -INC PPARAM
  38. -INC CCOPTIO
  39. -INC CCREEL
  40. -INC SMELEME
  41. *      POINTEUR KELEM.MELEME
  42.       POINTEUR KEXTO.MELEME
  43.       POINTEUR IBTLOC.MELEME
  44.       POINTEUR IPBTL2.MELEME
  45.       POINTEUR KTBES.MELEME
  46.       POINTEUR KTBES2.MELEME
  47. * anc  POINTEUR IMCAND.MELEME
  48. -INC TMATOP1      
  49. *-INC SMELEMX
  50.       POINTEUR LMCANS.MELEMX
  51.       POINTEUR IPBTL.MELEMX
  52.       POINTEUR KELEMX.MELEMX
  53. -INC SMLENTI
  54. *anc      POINTEUR KNNO.MLENTI
  55.       POINTEUR LIDXCA.MLENTI
  56.       POINTEUR LOKVOL.MLENTI
  57.       POINTEUR LNQUAL.MLENTI
  58.       POINTEUR LINDI.MLENTI
  59.       POINTEUR LINDJ.MLENTI
  60. -INC SMLREEL
  61.       POINTEUR IQUAL.MLREEL
  62.       POINTEUR LQUALS.MLREEL
  63. -INC SMCOORD
  64.       POINTEUR KCOORD.MCOORD
  65. *-INC SMETRIQ
  66.       POINTEUR KCMETR.METRIQ
  67. *-INC STRAVJ
  68.       POINTEUR TRAVK.TRAVJ
  69. *-INC STRAVL
  70.  
  71. *     
  72.       LOGICAL LOK
  73. *anc      LOGICAL LTOIBO
  74. *anc      LOGICAL LTOIBA
  75.       INTEGER JCAND
  76.       LOGICAL LCHANG
  77.       LOGICAL LCHTOP
  78. * Liste de topologies de maillages candidates
  79. *      SEGMENT ITCAND(0)
  80. * Liste de topologies de maillages candidats de plus petit volume non nul
  81. *      SEGMENT ITVOL(JG)
  82. * Liste de topologies de maillages candidats de plus petit volume
  83. * et de meilleure qualité
  84. *      SEGMENT ILQUAL(JG)
  85. *      SEGMENT ILIND(JG)
  86. *      SEGMENT JLIND(JG)
  87. *     
  88. * Executable statements
  89. *     
  90. *      WRITE(IOIMP,*) 'coucou topv2'
  91. * Il vaudrait mieux la lire !!
  92. *      XVTOL=XZPREC*1.D2
  93.       IDIMP1=IDIM+1
  94.       KCOORD=TRAVK.COORD
  95.       KPVIRT=TRAVK.PVIRT
  96.       KEXTO=TRAVK.TOPO
  97.       KCMETR=TRAVK.CMETR
  98. *
  99.       LMCANS=TRAVL.MCANS
  100.       LIDXCA=TRAVL.IDXCA
  101.       LOKVOL=TRAVL.OKVOL
  102.       LQUALS=TRAVL.QUALS
  103.       LNQUAL=TRAVL.NQUAL
  104.       LINDI=TRAVL.INDI
  105.       LINDJ=TRAVL.INDJ
  106. *      IPBTL=TRAVL.PBTL
  107. *     
  108. * Initialisations SEGMENT TRAVL : nb candidats=0     
  109. *                                 nbel maillage candidats=0
  110. *     
  111.       CALL TRLADJ(TRAVL,0,0,lchang,'topv2 : TRAVL')
  112.       IF (IERR.NE.0) RETURN
  113. *     
  114.       IPOPL2=0
  115.       KTBES=KEXTO
  116.       JCAND=0
  117. * Tests sur le maillage initial : uniquement des éléments à 3 noeuds en 2D
  118. * et des éléments à 4 noeuds en 3D
  119. *      SEGACT KEXTO
  120. *!      NLTLOC=KEXTO.NUM(/2)
  121.       NLTLOC=TRAVK.NVCOU
  122.       IF (NLTLOC.EQ.0) GOTO 7
  123.       NBSOUS=KEXTO.LISOUS(/1)
  124.       IF (NBSOUS.NE.0) THEN
  125.          WRITE(IOIMP,*) 'We want only elementary meshes'
  126.          GOTO 9999
  127.       ENDIF
  128.       NBNN=KEXTO.NUM(/1)
  129.       IF (NBNN.NE.IDIM+1) THEN
  130.          WRITE(IOIMP,*) 'Only simplices are allowed'
  131.          GOTO 9999
  132.       ENDIF
  133. *
  134. *      IARET=KELEM.NUM(/1)
  135. *      NBELEM=KELEM.NUM(/2)
  136.       IARET=KELEMX.NNCOU
  137.       NBELEM=KELEMX.NLCOU
  138. * Une vérif qu'on aurait dû faire depuis longtemps
  139. *      IF (IARET.LT.0.OR.NBELEM.LT.0) THEN
  140.       IF (IARET.LE.0.OR.NBELEM.LE.0) THEN
  141.          WRITE(IOIMP,*) 'KELEM nul dim=',IARET,NBELEM
  142.          goto 9999
  143.       ENDIF
  144.       LOK=.TRUE.
  145.       DO IIAR=1,IARET
  146. *         INOD=KELEM.NUM(IIAR,1)
  147.          INOD=KELEMX.NUMX(IIAR,1)
  148.          IF (INOD.EQ.0) THEN
  149.             WRITE(IOIMP,*) 'KELEMX noeud nul iiar=',IIAR
  150.             LOK=.FALSE.
  151.          ENDIF
  152.       ENDDO
  153.       IF (.NOT.LOK) THEN
  154.          WRITE(IOIMP,*) 'Ce cas nest pas bon du tout.'
  155.          goto 9999
  156.       ENDIF         
  157. * On remet le noeud KPVIRT au barycentre de KELEM (modifie le MCOORD !!!)
  158.       IF (KPVIRT.NE.0) THEN
  159. *     Pas besoin a priori.         segact mcoord*mod
  160. * fait plus haut         KCOORD=TRAVK.COORD
  161. *         write(ioimp,*) 'kpvirt=',kpvirt,'nntot=',xcoor(/1)/idimp1
  162.          DO L=1,IDIM
  163. *            write(ioimp,*) 'l=',l
  164.             xig=0.D0
  165. *            npoin=KELEM.NUM(/1)
  166.             npoin=KELEMX.NNCOU
  167.             DO in=1,npoin
  168. *               ip=kelem.num(in,1)
  169.                ip=kelemx.numx(in,1)
  170.                iref=idimp1*(ip-1)
  171.                xig=xig+kcoord.xcoor(iref+l)
  172.             ENDDO
  173.             iref=idimp1*(kpvirt-1)
  174. *            xold=xcoor(iref+l)
  175. *            write(ioimp,*) 'iref=',iref
  176.             kcoord.xcoor(iref+l)=xig/npoin
  177. *            dx=xcoor(iref+l)-xold
  178. *            if (abs(dx).gt.xzprec) write(ioimp,*) 'dx=',dx
  179. *            write(ioimp,*) 'dx=',dx
  180. *            write(ioimp,*) 'x=',xcoor(iref+l)
  181.          ENDDO
  182.       ENDIF
  183. * Calcul du volume de la topologie initiale
  184. *!!Erreur pour calculer le volume, il ne faut pas utiliser la métrique !
  185. *      CALL VOMSIM(KEXTO,IMET,KPVIRT,XVTLOC)
  186. *     CALL VOMSIM(KEXTO,0,KPVIRT,XVTLOC,XVTLOV)
  187.       IELDEB=1
  188.       IELFIN=TRAVK.NVCOU
  189.       CALL VOMSI2(KEXTO,IELDEB,IELFIN,0,KPVIRT,XVTLOC,XVTLOV)
  190.       IF (IERR.NE.0) RETURN
  191. *      JCAND=0
  192. * Note : on saute les topologies locales ayant un volume nul
  193. * car on n'arrivera pas à les améliorer
  194.       IF (XVTLOC.LE.XVTOL) GOTO 7
  195. *
  196. *     Génération des topologies candidates (stockage dans TRAVL: MCANS
  197. *     indexé par IDXCA)
  198. *     ICBES est le meilleur candidat par défaut (souvent=1, la topologie initiale)
  199. *     IPOPL2 est le candidat où on a créé un nouveau noeud (souvent le
  200. *     dernier, égal à travl.nccou)
  201. *       
  202. *ijob      CALL TOPV3(TRAVK,KELEM,IJOB,TRAVL,ICBES,IPOPL2)
  203. *kelemx      CALL TOPV3(TRAVK,KELEM,IAJNO,TRAVL,INCMA,ISTMA,JNASCM,ICBES,IPOPL2
  204.       CALL TOPV3(TRAVK,KELEMX,IAJNO,TRAVL,INCMA,ISTMA,JNASCM,ICBES
  205.      $     ,IPOPL2,iveri,impr) 
  206.       IF (IERR.NE.0) RETURN
  207. *      
  208. * On a les candidats dans ITCAND
  209. * On sélectionne ceux de volume minimum non nul
  210. *      
  211.       JCAND=TRAVL.NCCOU
  212.       if (impr.ge.4) then
  213.          write(ioimp,*) 'Apres gen candidat n=',jcand
  214.       endif
  215. * Sortie anticipée s'il n'y a qu'un candidat
  216. *      WRITE(IOIMP,*) 'JCAND=',JCAND
  217.       IF (JCAND.EQ.1) GOTO 8
  218. *      
  219.       XVMIN=XVTLOC
  220.       if (impr.ge.4) then
  221.          write(ioimp,*) 'volume vise XVMIN=',XVMIN,' virtuel=',XVTLOV
  222.       endif
  223. *
  224. *ijob      CALL TOPV4(TRAVL,IJOB,XVMIN,XVTLOV,XVTOL,KPVIRT)
  225.       CALL TOPV4(TRAVL,IALGO,XVMIN,XVTLOV,XVTOL,KPVIRT,impr)
  226.       IF (IERR.NE.0) RETURN
  227. *      
  228.       if (impr.ge.4) then
  229. *anc         write(ioimp,*) 'Apres selection vol ncandidat=',ITVOL(/1)
  230. *anc         Write(ioimp,*) (itvol(iii),iii=1,itvol(/1))
  231. *         write(ioimp,*) 'Apres selection vol ncandidat=',nvocou
  232. *         Write(ioimp,*) (lokvol.lect(iii),iii=1,nvocou)
  233.       endif
  234. * S'il n'en reste qu'un, on peut sauter le calcul des qualités ??????
  235.       IF (TRAVL.NVOCOU.EQ.1) THEN
  236.          ICBES=LOKVOL.LECT(NVOCOU)
  237.          GOTO 8
  238.       ENDIF
  239. *     
  240. *     Les candidats ayant le bon volume sont dans ITVOL
  241. *     On calcule les qualités de chaque élément des candidats et on ordonne.
  242. *
  243.       CALL TOPV5(TRAVL,XVTOL,IMET,IMOMET,XDENS,KCMETR,KPVIRT)
  244.       IF (IERR.NE.0) RETURN
  245.  
  246.       if (impr.ge.4) then
  247.          write(ioimp,*) 'Apres calcul qualité candidats'
  248.          DO i=1,NVOCOU
  249.             Write(ioimp,*) 'Candidat i=',i
  250.             ICAND=LOKVOL.LECT(i)
  251.             IELDEB=LIDXCA.LECT(ICAND)
  252.             IELFIN=LIDXCA.LECT(ICAND+1)-1
  253.             Write(ioimp,*) (lquals.prog(iii),iii=ieldeb,ielfin)
  254. *     $           ,lnqual.lect(icand))
  255.          enddo
  256.       endif
  257. *
  258. * Calcul des meilleurs candidats par maximum lexical
  259. *
  260.       CALL TOPV6(TRAVL,QTOL,ICBES)
  261. *
  262.       if (impr.ge.4) then
  263.          write(ioimp,*) 'Apres tri qualité candidats, le meilleur est :'
  264.      $        ,icbes
  265.       endif
  266.  
  267.  8    CONTINUE
  268. *
  269. *anc      KTBES2=KTBES
  270.       ICAND=ICBES
  271. *dbg      write(ioimp,*) 'icand,nccou=',icand,travl.nccou
  272.       IF (ICAND.NE.1) THEN
  273.          IELDEB=LIDXCA.LECT(ICAND)
  274.          IELFIN=LIDXCA.LECT(ICAND+1)-1
  275.          NBNN=LMCANS.NNCOU
  276.          NBELEM=IELFIN-IELDEB+1
  277.          NBSOUS=0
  278.          NBREF=0
  279. *            write(ioimp,*) 'icand,ieldeb,ielfin,nbnn,nbelem=',icand
  280. *     $           ,ieldeb,ielfin,nbnn,nbelem
  281.          SEGINI KTBES
  282.          KTBES.ITYPEL=LMCANS.ITYPEX
  283.          IDX=1
  284.          DO IEL=IELDEB,IELFIN
  285.             DO INO=1,NBNN
  286.                KTBES.NUM(INO,IDX)=LMCANS.NUMX(INO,IEL)
  287.             ENDDO
  288.             IDX=IDX+1
  289.          ENDDO
  290.       ENDIF
  291. *
  292. * Si on n'a pas sélectionné le candidat avec le noeud supplémentaire
  293. * on peut enlever le noeud
  294. *ijob      IF (IJOB.EQ.1.OR.IJOB.EQ.2) THEN 
  295.       IF (IAJNO.EQ.1) THEN 
  296.          IF (IPOPL2.NE.0.AND.ICBES.NE.IPOPL2) THEN
  297.             NPCOUN=TRAVK.NPCOU-1
  298. *     Remise à zéro : nécessaire uniquement pour la vérification ??
  299.             if (iveri.ge.1) then
  300.                IREF=(TRAVK.NPCOU-1)*IDIMP1
  301.                DO 11 I=1,IDIMP1
  302.                   KCOORD.XCOOR(IREF+I)=0.D0
  303.  11            CONTINUE
  304.                IF (KCMETR.NE.0) THEN
  305.                   DO 12 ININ=1,KCMETR.XIN(/1)
  306.                      KCMETR.XIN(ININ,TRAVK.NPCOU)=0.D0
  307.  12               CONTINUE
  308.                ENDIF
  309.             endif
  310.             TRAVK.NPCOU=NPCOUN
  311. * Désactivation temporaire car pas de ménage...            
  312. *tmp            if (iveri.ge.2) then
  313. *tmp            call vetopi(travk,'topv2 : Apres reduction travk')
  314. *tmp               if (ierr.ne.0) return
  315. *tmp            endif
  316.          ELSE
  317. *            write(ioimp,*) 'Nouveau noeud cree ','IPOPL2,ICBES,NPCOU='
  318. *     $           ,IPOPL2,ICBES,TRAVK.NPCOU
  319. *            IREF=(TRAVK.NPCOU-1)*IDIMP1
  320. *            write(ioimp,*)  (kcoord.xcoor(iref+iii),iii=1,idimp1)
  321.          ENDIF
  322.       ENDIF
  323.  7    CONTINUE
  324.       IF (KEXTO.EQ.KTBES) JCAND=-JCAND
  325.       IF (KEXTO.EQ.KTBES) THEN
  326.          LCHTOP=.FALSE.
  327.       ELSE
  328.          LCHTOP=.TRUE.
  329.       ENDIF
  330.       RETURN
  331. *     
  332. *     
  333. *     
  334.  9999 CONTINUE
  335.       MOTERR(1:8)='TOPV2   '
  336. * 349 2
  337. *Problème non prévu dans le s.p. %m1:8 contactez votre assistance
  338.       CALL ERREUR(349)
  339.       RETURN
  340. *     
  341. * End of subroutine TOPV2
  342. *     
  343.       END
  344.  
  345.  
  346.  

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