Télécharger raff2d.eso

Retour à la liste

Numérotation des lignes :

  1. C RAFF2D SOURCE GG250959 17/03/15 21:15:09 9344
  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 CCOPTIO
  52. -INC SMCOORD
  53. -INC CCGEOME
  54. -INC SMELEME
  55. -INC SMCHPOI
  56. -INC SMMODEL
  57. -INC SMCHAML
  58. C
  59. C=======================================================================
  60. C Initialisations - Declarations
  61. C=======================================================================
  62. SEGMENT ICPR((XCOOR(/1)/(IDIM+1)),2)
  63. SEGMENT XDEN((XCOOR(/1)/(IDIM+1)))
  64. SEGMENT KARETE(NBNDS,NCOL)
  65. SEGMENT KARET2(NBNDS,NCOL)
  66. SEGMENT KMILIE(NBNDS,NCOL)
  67. SEGMENT KARPOS(NBNDS)
  68. SEGMENT NUMNOE(INUMNO)
  69. SEGACT IPT2,ICPR,KARPOS,KARETE,KMILIE*MOD,MELVA2,KARET2
  70. NBNN=NBNNE(LNELM(2,(IPT2.ITYPEL-1)*2+1))
  71. INELM=LNELM(1,(IPT2.ITYPEL-1)*2+1)
  72. C
  73. C=======================================================================
  74. C Phase de raffinement 2D
  75. C=======================================================================
  76. C
  77. C 1/ Creation du squelette du maillage resultat
  78. NBELEM=IPT2.NUM(/2)-KARAF+INELM*KARAF
  79. NBSOUS=0
  80. NBREF=0
  81. SEGINI IPT4
  82. IPT4.ITYPEL=LNELM(2,(IPT2.ITYPEL-1)*2+1)
  83. LPO2=LPOS2(IPT2.ITYPEL)
  84. NBPT1=XCOOR(/1)/(IDIM+1)
  85. NBPTS=NBPT1+NACREE+KARAF*NBINTE(IPT2.ITYPEL)
  86. SEGADJ MCOORD
  87. SEGADJ XDEN
  88. INUMNO=NBRAF(IPT2.ITYPEL)
  89. SEGINI NUMNOE
  90. LCOMP=1
  91. C
  92. C 2/ Recherche des elements a raffiner
  93. DO 6 IARAF=1,MELVA2.VELCHE(/2)
  94. IF (MELVA2.VELCHE(1,IARAF).NE.1) THEN
  95. DO 1 IJKL=1,IPT2.NUM(/1)
  96. IPT4.NUM(IJKL,NBELEM)=IPT2.NUM(IJKL,IARAF)
  97. 1 CONTINUE
  98. NBELEM=NBELEM-1
  99. GOTO 6
  100. ENDIF
  101. C
  102. C 3/ Raffinement de l'element
  103. JCOMPT=0
  104. JPOS5=LPOS5(IPT2.ITYPEL)
  105. DO 4 I=1,NBRAF(IPT2.ITYPEL)
  106. JPOS1=LPOS1(1,I+LPOS2(IPT2.ITYPEL)-1)
  107. JLONG=LPOS1(2,I+LPOS2(IPT2.ITYPEL)-1)
  108. JLISN=LPOS3(IPT2.ITYPEL)+JCOMPT
  109. LTYPNO=JTYPNO(JPOS5-1+I)
  110. IF (LTYPNO.EQ.0) THEN
  111. NPTA=IPT2.NUM(LISNOE(JLISN),IARAF)
  112. NPTB=IPT2.NUM(LISNOE(JLISN+1),IARAF)
  113. NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  114. NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  115. DO 2 K=1,MAX(1,KARPOS(NMIN))
  116. IF (KARETE(NMIN,K).EQ.NMAX) NEXIST=K
  117. 2 CONTINUE
  118. IF (KMILIE(NMIN,NEXIST).GT.0) THEN
  119. NUMNOE(I)=KMILIE(NMIN,NEXIST)
  120. KMILIE(NMIN,NEXIST)=0
  121. JCOMPT=JCOMPT+JLONG
  122. GOTO 4
  123. ELSE
  124. NBPT1=NBPT1+1
  125. NUMNOE(I)=NBPT1
  126. KMILIE(NMIN,NEXIST)=NBPT1
  127. ENDIF
  128. ELSEIF (LTYPNO.EQ.7) THEN
  129. NBPT1=NBPT1+1
  130. NUMNOE(I)=NBPT1
  131. ELSE
  132. WRITE(*,*) 'ERREUR, on ne devrait jamais passer ici !'
  133. ENDIF
  134. C
  135. C 4/ Creation des nouveaux points
  136. XPT=0.
  137. YPT=0.
  138. XDEN1=0.D0
  139. XDENMIN = 1000.
  140. DO 3 J=1,JLONG
  141. NGLOB=IPT2.NUM(LISNOE(JLISN-1+J),IARAF)
  142. XINI=XCOOR((NGLOB-1)*(IDIM+1)+1)
  143. YINI=XCOOR((NGLOB-1)*(IDIM+1)+2)
  144. XPT=XPT+XINI*XCOEFF(JPOS1-1+J)
  145.  
  146. YPT=YPT+YINI*XCOEFF(JPOS1-1+J)
  147.  
  148. XDEN1=XDEN1+XDEN(NGLOB)*XCOEFF(JPOS1-1+J)
  149. IF (XDEN(NGLOB).LT.XDENMIN) XDENMIN = XDEN(NGLOB)
  150. c if (iaraf .eq. 10) then
  151. c write (*,*) 'nglob', nglob ,'nbpt1', nbpt1
  152. c write (*,*) 'xpt ', xpt
  153. c write (*,*) 'ypt ', ypt
  154. c write (*,*) 'xden1 ', xden1
  155. c write (*,*) 'XDEN(NGLOB)', XDEN(NGLOB)
  156. c write (*,*) 'xcoeff1 ', XCOEFF(JPOS1-1+J)
  157. c endif
  158. 3 CONTINUE
  159. c seuil de densité a xdenmin
  160. IF (XDEN1.LT.XDENMIN) XDEN1 = XDENMIN
  161.  
  162. XCOOR((NBPT1-1)*(IDIM+1)+1)=XPT
  163. XCOOR((NBPT1-1)*(IDIM+1)+2)=YPT
  164. XCOOR((NBPT1-1)*(IDIM+1)+3)=XDEN1
  165. XDEN(NBPT1)=XDEN1
  166. JCOMPT=JCOMPT+JLONG
  167. 4 CONTINUE
  168. JPOS4=LPOS4(IPT2.ITYPEL)
  169. C
  170. C 5/ Creation des nouveaux elements
  171. C On remplit la portion de IPT4 relative aux elements crees a partir
  172. C de la division de l'element IARAF (indice de boucle 1).
  173. C Cette portion de IPT4 contient les colonnes dont la valeur s'etend
  174. C de INELM*(LCOMP-1)+1 a INELM*LCOMP.
  175. DO 5 J=1,INELM
  176. DO 5 I=1,IPT4.NUM(/1)
  177. NTEMP=LIELM(JPOS4-1+NBNN*(J-1)+I)
  178. IF (NTEMP.GT.NBNN) THEN
  179. IPT4.NUM(I,INELM*(LCOMP-1)+J)=NUMNOE(NTEMP-NBNN)
  180. ELSE
  181. IPT4.NUM(I,INELM*(LCOMP-1)+J)=IPT2.NUM(NTEMP,IARAF)
  182. ENDIF
  183. 5 CONTINUE
  184. LCOMP=LCOMP+1
  185. 6 CONTINUE
  186. C
  187. C=======================================================================
  188. C Preparation du maillage de relations
  189. C=======================================================================
  190. C On cree un maillage IPT5 contenant tous les noeuds soumis a des
  191. C relations
  192. C
  193. C 1/ Comptage du nombre de noeuds soumis a des relations
  194. NBELEM=0
  195. DO 7 J=1,KMILIE(/2)
  196. DO 7 I=1,KMILIE(/1)
  197. IF (KARET2(I,J).NE.2) GOTO 7
  198. IF (KMILIE(I,J).GT.0) NBELEM=NBELEM+1
  199. 7 CONTINUE
  200. C
  201. C 2/ Creation de IPT5
  202. IPT5=0
  203. IF (NBELEM.EQ.0) GOTO 999
  204. NBNN=5
  205. NBREF=0
  206. NBSOUS=0
  207. SEGINI IPT5
  208. IPT5.ITYPEL=48
  209. C
  210. C 3/ Renseignement des noeuds support des relations
  211. DO 8 J=1,KMILIE(/2)
  212. DO 8 I=1,KMILIE(/1)
  213. IF (KARET2(I,J).NE.2) GOTO 8
  214. IF (KMILIE(I,J).GT.0) THEN
  215. NBREF=NBREF+1
  216. IPT5.NUM(1,NBREF)=KMILIE(I,J)
  217. ENDIF
  218. 8 CONTINUE
  219. C
  220. C 4/ Recherche des noeuds formant les relations
  221. DO 14 IARAF=1,MELVA2.VELCHE(/2)
  222. IF (MELVA2.VELCHE(1,IARAF).NE.1) GOTO 14
  223. JCOMPT=0
  224. JPOS5=LPOS5(IPT2.ITYPEL)
  225. DO 13 I=1,NBRAF(IPT2.ITYPEL)
  226. JPOS1=LPOS1(1,I+LPOS2(IPT2.ITYPEL)-1)
  227. JLONG=LPOS1(2,I+LPOS2(IPT2.ITYPEL)-1)
  228. JLISN=LPOS3(IPT2.ITYPEL)+JCOMPT
  229. LTYPNO=JTYPNO(JPOS5-1+I)
  230. IF (LTYPNO.NE.0) GOTO 12
  231. NPTA=IPT2.NUM(LISNOE(JLISN),IARAF)
  232. NPTB=IPT2.NUM(LISNOE(JLISN+1),IARAF)
  233. NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  234. NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  235. DO 9 K=1,MAX(1,KARPOS(NMIN))
  236. IF (KARETE(NMIN,K).EQ.NMAX) NEXIST=K
  237. 9 CONTINUE
  238. IF (KARET2(NMIN,NEXIST).NE.2) GOTO 12
  239. IF (KMILIE(NMIN,NEXIST).EQ.0) GOTO 12
  240. NEXIS5=0
  241. DO 10 MM=1,IPT5.NUM(/2)
  242. INOEGL=KMILIE(NMIN,NEXIST)
  243. INRELA=IPT5.NUM(1,MM)
  244. IF (INOEGL.EQ.INRELA) NEXIS5=MM
  245. 10 CONTINUE
  246. C
  247. C 5/ Renseignement des noeuds formant les relations
  248. DO 11 IHH=1,JLONG
  249. IPT5.NUM(1+IHH,NEXIS5)=IPT2.NUM(LISNOE(JLISN-1+IHH),IARAF)
  250. 11 CONTINUE
  251. IF (JLONG.EQ.2) IPT5.NUM(5,NEXIS5)=1
  252. IF (JLONG.EQ.3) IPT5.NUM(5,NEXIS5)=2
  253. 12 CONTINUE
  254. JCOMPT=JCOMPT+JLONG
  255. 13 CONTINUE
  256. 14 CONTINUE
  257. C
  258. C=======================================================================
  259. C Fin du programme
  260. C=======================================================================
  261. 999 CONTINUE
  262. RETURN
  263. END
  264.  
  265.  
  266.  
  267.  
  268.  
  269.  
  270.  
  271.  
  272.  

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