Télécharger la plaquette Cast3M


fr en

Télécharger raff2d.eso

Retour à la liste

Numérotation des lignes : 

raff2d

  1. C RAFF2D    SOURCE    CB215821  24/04/12    21:17:01     11897          
  2.       SUBROUTINE RAFF2D(IPT2,ICPR,KARPOS,KARETE,KMILIE,MELVA2,NACREE,
  3.      .KARAF,IPT4,IDCP,IPT5,KARET2,XDEN)
  4. c entrees de raff2D:
  5. c           - IPT2:   maillage élémentaire à raffiner 
  6. c           - ICPR:   tableau de passage noeuds local/global
  7. c           - KARPOS: tableau du nb d'arretes par noeuds
  8. c           - KARETE: tableau des d'arretes KARETE(i,j)=k 
  9. c           - KMILIE: tableau des hanging nodes associés au arêtes
  10. c                     si ancienne relation 
  11. c                        KMILIE(i,j)=n : le noeud global n est situe
  12. c                        au milieu de l'arete [i,k]
  13. c                     si nouvelle relation
  14. c                        KMILIE(i,j)=-1 : il faut construire un noeud
  15. c                        au milieu de l'arete [i,k]
  16. c           - MELVA2: champ par elem valait 1 pour les élem à raffiner
  17. c                     et 0 pour les autres 
  18. c           - NACREE: nombre de noeuds à créer dans ipt2
  19. c           - KARAF:  nombre d'élément à raffiner dans ipt2
  20. c           - IDCP :  ??? variable non utilisé ???
  21. c           - XDEN :  densité aux noeuds en notation globale
  22. c           - KARET2: tableau du nombre d'éléments touchant les arretes
  23. c sorties: 
  24. c           - KMILIE: tableau des hanging nodes associés au arêtes 
  25. c                     (en sortie: à partir des relations donnée en 
  26. c                     entrée et de celles crées
  27. c                     si anciene relation à garder
  28. c                       KMILIE(i,j)=n : le noeud global n est situe
  29. c                       au milieu de l'arete [i,k]
  30. c                     si anciene relation à supprimer
  31. c                       KMILIE(i,j)=0 : 
  32. c                     si nouvelle relation
  33. c                       KMILIE(i,j)=n : le noeud n à été construit au
  34. c                       milieu de l'arete [i,k]
  35. c           - MELEME: Maillage élémentaire raffiné à partir de ipt2
  36. c                     sans relations
  37. c           - IPT8 :  Maillage de relation incomplet 
  38. c                     si J nouvel élément 48:
  39. c                       IPT8.NUM(1,j) = nouveau hanging node
  40. c                       IPT8.NUM(2:4,j) = autres noeuds de la relation
  41. c                       IPT8.NUM(5,j) = 1 si arête a 2 noeuds
  42. c                                       2 si arête a 3 noeuds
  43. c                     si J ancien élément 48:
  44. c                       IPT8.NUM(1,j) = nouveau hanging node
  45. c                       IPT8.NUM(2:5,j) = 0
  46. c           - XDEN :  densité aux noeuds en notation globale
  47. c                     interpolé aux nouveaux noeuds.
  48.  
  49.       IMPLICIT REAL*8(A-H,O-Z)
  50.       IMPLICIT INTEGER(I-N)
  51. -INC PPARAM
  52. -INC CCOPTIO
  53. -INC SMCOORD
  54. -INC CCGEOME
  55. -INC SMELEME
  56. -INC SMCHPOI
  57. -INC SMMODEL
  58. -INC SMCHAML
  59. C
  60. C=======================================================================
  61. C                   Initialisations - Declarations
  62. C=======================================================================
  63.       SEGMENT ICPR(nbpts,2)
  64.       SEGMENT XDEN(nbpts)
  65.       SEGMENT KARETE(NBNDS,NCOL)
  66.       SEGMENT KARET2(NBNDS,NCOL)
  67.       SEGMENT KMILIE(NBNDS,NCOL)
  68.       SEGMENT KARPOS(NBNDS)
  69.       SEGMENT NUMNOE(INUMNO)
  70.       SEGACT IPT2,ICPR,KARPOS,KARETE,KMILIE*MOD,MELVA2,KARET2
  71.       NBNN=NBNNE(LNELM(2,(IPT2.ITYPEL-1)*2+1))
  72.       INELM=LNELM(1,(IPT2.ITYPEL-1)*2+1)
  73. C
  74. C=======================================================================
  75. C                   Phase de raffinement 2D
  76. C=======================================================================
  77. C
  78. C 1/ Creation du squelette du maillage resultat
  79.       NBELEM=IPT2.NUM(/2)-KARAF+INELM*KARAF
  80.       NBSOUS=0
  81.       NBREF=0
  82.       SEGINI IPT4
  83.       IPT4.ITYPEL=LNELM(2,(IPT2.ITYPEL-1)*2+1)
  84.       LPO2=LPOS2(IPT2.ITYPEL)
  85.       segact mcoord*mod
  86.       NBPT1=nbpts
  87.       NBPTS=NBPT1+NACREE+KARAF*NBINTE(IPT2.ITYPEL)
  88.       SEGADJ MCOORD
  89.       SEGADJ XDEN
  90.       INUMNO=NBRAF(IPT2.ITYPEL)
  91.       SEGINI NUMNOE
  92.       LCOMP=1
  93. C
  94. C 2/ Recherche des elements a raffiner
  95.       DO 6 IARAF=1,MELVA2.VELCHE(/2)
  96.         IF (MELVA2.VELCHE(1,IARAF).NE.1) THEN
  97.           DO 1 IJKL=1,IPT2.NUM(/1)
  98.             IPT4.NUM(IJKL,NBELEM)=IPT2.NUM(IJKL,IARAF)
  99.    1      CONTINUE
  100.           NBELEM=NBELEM-1
  101.           GOTO 6
  102.         ENDIF
  103. C
  104. C 3/ Raffinement de l'element
  105.         JCOMPT=0
  106.         JPOS5=LPOS5(IPT2.ITYPEL)
  107.         DO 4 I=1,NBRAF(IPT2.ITYPEL)
  108.           JPOS1=LPOS1(1,I+LPOS2(IPT2.ITYPEL)-1)
  109.           JLONG=LPOS1(2,I+LPOS2(IPT2.ITYPEL)-1)
  110.           JLISN=LPOS3(IPT2.ITYPEL)+JCOMPT
  111.           LTYPNO=JTYPNO(JPOS5-1+I)
  112.           IF (LTYPNO.EQ.0) THEN
  113.             NPTA=IPT2.NUM(LISNOE(JLISN),IARAF)
  114.             NPTB=IPT2.NUM(LISNOE(JLISN+1),IARAF)
  115.             NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  116.             NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  117.             DO 2 K=1,MAX(1,KARPOS(NMIN))
  118.               IF (KARETE(NMIN,K).EQ.NMAX) NEXIST=K
  119.    2        CONTINUE
  120.             IF (KMILIE(NMIN,NEXIST).GT.0) THEN
  121.               NUMNOE(I)=KMILIE(NMIN,NEXIST)
  122.               KMILIE(NMIN,NEXIST)=0
  123.               JCOMPT=JCOMPT+JLONG
  124.               GOTO 4
  125.             ELSE
  126.               NBPT1=NBPT1+1
  127.               NUMNOE(I)=NBPT1
  128.               KMILIE(NMIN,NEXIST)=NBPT1
  129.             ENDIF
  130.           ELSEIF (LTYPNO.EQ.7) THEN
  131.             NBPT1=NBPT1+1
  132.             NUMNOE(I)=NBPT1
  133.           ELSE
  134.             WRITE(*,*) 'ERREUR, on ne devrait jamais passer ici !'
  135.           ENDIF
  136. C
  137. C 4/ Creation des nouveaux points
  138.           XPT=0.
  139.           YPT=0.
  140.           XDEN1=0.D0
  141.           XDENMIN = 1000.
  142.           DO 3 J=1,JLONG
  143.             NGLOB=IPT2.NUM(LISNOE(JLISN-1+J),IARAF)
  144.             XINI=XCOOR((NGLOB-1)*(IDIM+1)+1)
  145.             YINI=XCOOR((NGLOB-1)*(IDIM+1)+2)
  146.             XPT=XPT+XINI*XCOEFF(JPOS1-1+J)
  147.  
  148.             YPT=YPT+YINI*XCOEFF(JPOS1-1+J)
  149.  
  150.             XDEN1=XDEN1+XDEN(NGLOB)*XCOEFF(JPOS1-1+J)
  151.             IF (XDEN(NGLOB).LT.XDENMIN) XDENMIN = XDEN(NGLOB) 
  152. c            if (iaraf .eq. 10) then
  153. c              write (*,*) 'nglob', nglob ,'nbpt1', nbpt1
  154. c              write (*,*) 'xpt ', xpt
  155. c              write (*,*) 'ypt ', ypt
  156. c              write (*,*) 'xden1 ', xden1
  157. c              write (*,*) 'XDEN(NGLOB)', XDEN(NGLOB)
  158. c              write (*,*) 'xcoeff1 ', XCOEFF(JPOS1-1+J)
  159. c            endif
  160.    3      CONTINUE
  161. c seuil de densité a xdenmin
  162.           IF (XDEN1.LT.XDENMIN) XDEN1 = XDENMIN
  163.  
  164.           XCOOR((NBPT1-1)*(IDIM+1)+1)=XPT
  165.           XCOOR((NBPT1-1)*(IDIM+1)+2)=YPT
  166.           XCOOR((NBPT1-1)*(IDIM+1)+3)=XDEN1
  167.           XDEN(NBPT1)=XDEN1
  168.           JCOMPT=JCOMPT+JLONG
  169.    4    CONTINUE
  170.         JPOS4=LPOS4(IPT2.ITYPEL)
  171. C
  172. C 5/ Creation des nouveaux elements
  173. C On remplit la portion de IPT4 relative aux elements crees a partir
  174. C de la division de l'element IARAF (indice de boucle 1).
  175. C Cette portion de IPT4 contient les colonnes dont la valeur s'etend
  176. C de INELM*(LCOMP-1)+1 a INELM*LCOMP.
  177.         DO 5 J=1,INELM
  178.         DO 5 I=1,IPT4.NUM(/1)
  179.           NTEMP=LIELM(JPOS4-1+NBNN*(J-1)+I)
  180.           IF (NTEMP.GT.NBNN) THEN
  181.             IPT4.NUM(I,INELM*(LCOMP-1)+J)=NUMNOE(NTEMP-NBNN)
  182.           ELSE
  183.             IPT4.NUM(I,INELM*(LCOMP-1)+J)=IPT2.NUM(NTEMP,IARAF)
  184.           ENDIF
  185.    5    CONTINUE
  186.         LCOMP=LCOMP+1
  187.    6  CONTINUE
  188. C
  189. C=======================================================================
  190. C              Preparation du maillage de relations
  191. C=======================================================================
  192. C On cree un maillage IPT5 contenant tous les noeuds soumis a des
  193. C relations
  194. C
  195. C 1/ Comptage du nombre de noeuds soumis a des relations
  196.       NBELEM=0
  197.       DO 7 J=1,KMILIE(/2)
  198.       DO 7 I=1,KMILIE(/1)
  199.         IF (KARET2(I,J).NE.2) GOTO 7
  200.         IF (KMILIE(I,J).GT.0) NBELEM=NBELEM+1
  201.    7  CONTINUE
  202. C
  203. C 2/ Creation de IPT5
  204.       IPT5=0
  205.       IF (NBELEM.EQ.0) GOTO 999
  206.       NBNN=5
  207.       NBREF=0
  208.       NBSOUS=0
  209.       SEGINI IPT5
  210.       IPT5.ITYPEL=48
  211. C
  212. C 3/ Renseignement des noeuds support des relations
  213.       DO 8 J=1,KMILIE(/2)
  214.       DO 8 I=1,KMILIE(/1)
  215.         IF (KARET2(I,J).NE.2) GOTO 8
  216.         IF (KMILIE(I,J).GT.0) THEN
  217.           NBREF=NBREF+1
  218.           IPT5.NUM(1,NBREF)=KMILIE(I,J)
  219.         ENDIF
  220.    8  CONTINUE
  221. C
  222. C 4/ Recherche des noeuds formant les relations
  223.       DO 14 IARAF=1,MELVA2.VELCHE(/2)
  224.         IF (MELVA2.VELCHE(1,IARAF).NE.1) GOTO 14
  225.         JCOMPT=0
  226.         JPOS5=LPOS5(IPT2.ITYPEL)
  227.         DO 13 I=1,NBRAF(IPT2.ITYPEL)
  228.           JPOS1=LPOS1(1,I+LPOS2(IPT2.ITYPEL)-1)
  229.           JLONG=LPOS1(2,I+LPOS2(IPT2.ITYPEL)-1)
  230.           JLISN=LPOS3(IPT2.ITYPEL)+JCOMPT
  231.           LTYPNO=JTYPNO(JPOS5-1+I)
  232.           IF (LTYPNO.NE.0) GOTO 12
  233.           NPTA=IPT2.NUM(LISNOE(JLISN),IARAF)
  234.           NPTB=IPT2.NUM(LISNOE(JLISN+1),IARAF)
  235.           NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  236.           NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  237.           DO 9 K=1,MAX(1,KARPOS(NMIN))
  238.             IF (KARETE(NMIN,K).EQ.NMAX) NEXIST=K
  239.    9      CONTINUE
  240.           IF (KARET2(NMIN,NEXIST).NE.2) GOTO 12
  241.           IF (KMILIE(NMIN,NEXIST).EQ.0) GOTO 12
  242.           NEXIS5=0
  243.           DO 10 MM=1,IPT5.NUM(/2)
  244.             INOEGL=KMILIE(NMIN,NEXIST)
  245.             INRELA=IPT5.NUM(1,MM)
  246.             IF (INOEGL.EQ.INRELA) NEXIS5=MM
  247.   10      CONTINUE
  248. C
  249. C 5/ Renseignement des noeuds formant les relations
  250.           DO 11 IHH=1,JLONG
  251.             IPT5.NUM(1+IHH,NEXIS5)=IPT2.NUM(LISNOE(JLISN-1+IHH),IARAF)
  252.   11      CONTINUE
  253.           IF (JLONG.EQ.2) IPT5.NUM(5,NEXIS5)=1
  254.           IF (JLONG.EQ.3) IPT5.NUM(5,NEXIS5)=2
  255.   12      CONTINUE
  256.         JCOMPT=JCOMPT+JLONG
  257.   13    CONTINUE
  258.   14  CONTINUE
  259. C
  260. C=======================================================================
  261. C                      Fin du programme
  262. C=======================================================================
  263.  999  CONTINUE
  264.       RETURN
  265.       END
  266.  
  267.  
  268.  
  269.  
  270.  
  271.  
  272.  
  273.  
  274.  
  275.  
  276.  
  277.  

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