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