Télécharger raff.eso

Retour à la liste

Numérotation des lignes :

  1. C RAFF SOURCE GG250959 17/03/15 21:15:02 9344
  2. SUBROUTINE RAFF
  3. IMPLICIT REAL*8(A-H,O-Z)
  4. IMPLICIT INTEGER(I-N)
  5. -INC CCOPTIO
  6. -INC SMCOORD
  7. -INC CCGEOME
  8. -INC SMELEME
  9. -INC SMCHPOI
  10. -INC SMMODEL
  11. -INC SMCHAML
  12. c-INC CCREEL
  13. C
  14. C--------------------------------
  15. C B) VARIABLES RELATIVES A LA 2D !
  16. C--------------------------------
  17. SEGMENT KARETE(NBNDS,NCOL)
  18. SEGMENT KARET2(NBNDS,NCOL)
  19. SEGMENT KMILIE(NBNDS,NCOL)
  20. SEGMENT KARPOS(NBNDS)
  21. C
  22. C--------------------------------
  23. C C) VARIABLES RELATIVES A LA 3D !
  24. C--------------------------------
  25. SEGMENT IRELA(NBNDS)
  26. SEGMENT JPLANS(JPLA1,JPLA2)
  27. SEGMENT JPLAN2(JPLA1,JPLA2)
  28. SEGMENT JPLAN3(JPLA1,JPLA2)
  29. SEGMENT JPLCOM(JPLA1)
  30. SEGMENT JNOEFA(JNBFA,5)
  31. SEGMENT IWORK2(JNBSOM)
  32. SEGMENT JFARAF(JNBFA,LLLL)
  33. SEGMENT IWORK1(IJKLMN)
  34. C--------------------------
  35. C A) VARIABLES 'GENERALES' !
  36. C--------------------------
  37. SEGMENT ICPR((XCOOR(/1)/(IDIM+1)),2)
  38. SEGMENT LTYPE(NBLIS,2)
  39. SEGMENT LTYPE2(NBLIS2,4)
  40. SEGMENT LTYPE3(NBLIS3,2)
  41. SEGMENT LTYPE4(NBLIS3,2)
  42. c SEGMENT ICHPO
  43. c REAL*8 XR(NBPTCH),YR(NBPTCH),ZR(NBPTCH),DEN(NBPTCH)
  44. c ENDSEGMENT
  45. SEGMENT IDENSI
  46. REAL*8 DENSI(NBNDS)
  47. ENDSEGMENT
  48. SEGMENT XDEN((XCOOR(/1)/(IDIM+1)))
  49. XPETI= 1E-30
  50. C
  51. C---------------------------------------------------
  52. C D) VARIABLES RELATIVES AUX MAILLAGES DE RELATIONS !
  53. C---------------------------------------------------
  54. POINTEUR IPT10.MELEME
  55. POINTEUR IPT11.MELEME
  56. POINTEUR IPT12.MELEME
  57. POINTEUR IPT13.MELEME
  58.  
  59. POINTEUR IPT14.MELEME
  60. POINTEUR IPT15.MELEME
  61. POINTEUR IPT16.MELEME
  62. SEGMENT MAIREL(NBMAIL)
  63. C=======================================================================
  64. C Lecture des donnees utilisateurs
  65. C=======================================================================
  66. C On lit le maillage a raffiner et le champoint de densite a affecter
  67. c ipt1 : maillage initial à raffiner
  68. c mchpo1 : chapoint de densité voulu(que devrai être définit sur ipt1)
  69. CALL LIROBJ('MAILLAGE',IPT1,1,IRETOU)
  70. IF (IERR.NE.0) RETURN
  71. CALL LIROBJ('CHPOINT',MCHPO1,1,IRETOU)
  72. IF (IERR.NE.0) RETURN
  73. C
  74. C=======================================================================
  75. C Initialisations - Declarations
  76. C=======================================================================
  77.  
  78. NBOUCL=0
  79. NPTINI=(XCOOR(/1)/(IDIM+1))
  80. SEGACT IPT1
  81.  
  82. IJKLMN=4
  83. SEGINI IWORK1
  84. IWORK1(1)=4
  85. IWORK1(2)=8
  86. IWORK1(3)=6
  87. IWORK1(4)=10
  88.  
  89. SEGINI,XDEN
  90.  
  91. C
  92. C=======================================================================
  93. C Tests sur les donnees
  94. C=======================================================================
  95. C On teste pour savoir si le maillage IPT1 (maillage 1) a un type
  96. C d'elements implemante.
  97. c TRI3 TRI6 QUA4 QUA8 CUB8 CUB26 PRY6 PR15 TET4 TET10 SURE
  98. C les pyramides 5 et 13 ont été débranchés
  99. c création du tableau
  100. NBLIS=MAX(1,IPT1.LISOUS(/1))
  101. NSURE=0
  102. SEGINI LTYPE
  103. DO 1 INB=1,MAX(1,IPT1.LISOUS(/1))
  104. IF (IPT1.LISOUS(/1).EQ.0) THEN
  105. IPT2=IPT1
  106. ELSE
  107. IPT2=IPT1.LISOUS(INB)
  108. SEGACT IPT2
  109. ENDIF
  110. NT=IPT2.ITYPEL
  111. IF (NT.EQ.48) NSURE=NSURE+1
  112. LTYPE(INB,1)=INB
  113. LTYPE(INB,2)=NT
  114. IF (NT.EQ.16) THEN
  115. write (*,*) 'Attention le raffinement PRI6 a debrancher !'
  116. ENDIF
  117. IF (NT.EQ.23) THEN
  118. write (*,*) 'Attention le raffinement TET4 a debrancher !'
  119. ENDIF
  120. IF ((NT.NE.4).AND.(NT.NE.6).AND.(NT.NE.8).AND.(NT.NE.10).AND.
  121. + (NT.NE.14).AND.(NT.NE.15).AND.(NT.NE.16).AND.(NT.NE.17).AND.
  122. + (NT.NE.23).AND.(NT.NE.24).AND.(NT.NE.48)) THEN
  123. WRITE(*,*) 'Erreur, type d element non pris en compte'
  124. RETURN
  125. ENDIF
  126. c write(*,*) 'type d elem du maillage de depart', LTYPE(INB,1),
  127. c + LTYPE(INB,2)
  128.  
  129. C IF ((NT.NE.4).AND.(NT.NE.6).AND.(NT.NE.8).AND.(NT.NE.10).AND.
  130. C + (NT.NE.14).AND.(NT.NE.15).AND.(NT.NE.16).AND.(NT.NE.17).AND.
  131. C + (NT.NE.23).AND.(NT.NE.24).AND.(NT.NE.48).AND.(NT.NE.25).AND.
  132. C + (NT.NE.26)) RETURN
  133. IF ((IPT1.LISOUS(/1).EQ.0).AND.(NT.EQ.48)) RETURN
  134. IF (IPT1.LISOUS(/1).NE.0) SEGDES IPT2
  135. 1 CONTINUE
  136.  
  137. C
  138. C=======================================================================
  139. C Separation du maillage initial
  140. C=======================================================================
  141. C Si le maillage initial comporte des lisous de type 48, on le separe
  142. C en deux maillages dont un ne sera compose que des lisous de type 48.
  143. C Maillage d'elements de type 48 : IPT6
  144. C Maillage sans elements de type 48 : IPT5
  145. C IPT1 reste le maillage initial complet
  146. LCONF2=0
  147. IF (NSURE.EQ.0) GOTO 7
  148. LCONF2=1
  149. 7 CONTINUE
  150. C-------------------------------
  151. C A) Creation du maillage IPT6 !
  152. C-------------------------------
  153. C A.1/** Cas ou il n'y a qu'un lisous de type 48
  154. IF (NSURE.EQ.1) THEN
  155. DO 2 K=1,LTYPE(/1)
  156. IF (LTYPE(K,2).EQ.48) IPT6=IPT1.LISOUS(LTYPE(K,1))
  157. 2 CONTINUE
  158. ENDIF
  159. C
  160. C A.2/** Cas ou il y a plusieurs lisous de type 48
  161. IF (NSURE.GT.1) THEN
  162. NBREF=0
  163. NBSOUS=NSURE
  164. NBNN=0
  165. NBELEM=0
  166. SEGINI IPT6
  167. KLISOU=0
  168. DO 3 K=1,LTYPE(/1)
  169. IF (LTYPE(K,2).EQ.48) THEN
  170. KLISOU=KLISOU+1
  171. IPT6.LISOUS(KLISOU)=IPT1.LISOUS(LTYPE(K,1))
  172. ENDIF
  173. 3 CONTINUE
  174. ENDIF
  175. IF (NSURE.EQ.0) THEN
  176. NBREF=0
  177. NBSOUS=NSURE
  178. NBNN=0
  179. NBELEM=0
  180. SEGINI IPT6
  181. ENDIF
  182. C
  183. C------------------------------
  184. C B) Creation de IPT5 !
  185. C------------------------------
  186. C
  187. NBSOUS=LTYPE(/1)-NSURE
  188. NBREF=0
  189. NBNN=0
  190. NBELEM=0
  191. C B.0/** Cas ou il n'y a qu'un lisous en tout
  192. IF ((NBSOUS.EQ.1).AND.(NSURE.EQ.0)) THEN
  193. NBSOUS=0
  194. SEGINI IPT5
  195. IPT5=IPT1
  196. c write (*,*) 'un seul sous domaine en tout', IPT5
  197. ENDIF
  198. C B.1/** Cas ou il n'y a qu'un lisous de type different de 48
  199. IF ((NBSOUS.EQ.1).AND.(NSURE.NE.0)) THEN
  200. SEGINI IPT5
  201. DO 4 K=1,LTYPE(/1)
  202. IF (LTYPE(K,2).NE.48) IPT5=IPT1.LISOUS(LTYPE(K,1))
  203. 4 CONTINUE
  204. c write (*,*) 'un seul sous domaine non sure', IPT5
  205. ENDIF
  206. C
  207. C B.2/** Cas ou il y a plusieurs lisous de types differents de 48
  208. IF (NBSOUS.GT.1) THEN
  209. SEGINI IPT5
  210. c write (*,*) 'plusieurs sous domaines non sure', IPT5
  211. KLISOU=0
  212. DO 5 K=1,LTYPE(/1)
  213. IF (LTYPE(K,2).NE.48) THEN
  214. KLISOU=KLISOU+1
  215. IPT5.LISOUS(KLISOU)=IPT1.LISOUS(LTYPE(K,1))
  216. ENDIF
  217. 5 CONTINUE
  218. ENDIF
  219.  
  220. c compteur de relaton globale
  221. LCONF3=LCONF2
  222. C=======================================================================
  223. C Debut de la boucle sur les itérations de Raffinement
  224. C=======================================================================
  225. 9 CONTINUE
  226. C
  227. C=======================================================================
  228. C Renumerotation locale des noeuds du maillage
  229. C=======================================================================
  230. C ICPR(I,1)=J : Le Ieme noeud global est le Jeme local
  231. C ICPR(I,2)=J : Le Ieme noeud global appartient a J elements
  232.  
  233. SEGACT IPT5
  234.  
  235. NBOUCL=NBOUCL+1
  236. c write(*,*) 'RAFF: ITERATION',NBOUCL
  237. c write(*,*) 'ipt5 : ', ipt5
  238. c write(*,*) 'ipt6 : ', ipt6
  239. c write(*,*) 'lconf3 : ', lconf3
  240.  
  241. SEGINI ICPR
  242. NBNDS=0
  243. NCOL=0
  244. DO 11 INB=1,MAX(1,IPT5.LISOUS(/1))
  245. IF (IPT5.LISOUS(/1).EQ.0) THEN
  246. IPT2=IPT5
  247. ELSE
  248. IPT2=IPT5.LISOUS(INB)
  249. SEGACT IPT2
  250. ENDIF
  251. IF (IPT2.ITYPEL.EQ.48) GOTO 11
  252. DO 10 J=1,IPT2.NUM(/2)
  253. DO 10 I=1,IPT2.NUM(/1)
  254. ITMP=IPT2.NUM(I,J)
  255. IF (ICPR(ITMP,1).NE.0) THEN
  256. ICPR(ITMP,2)=ICPR(ITMP,2)+1
  257. NCOL=MAX(NCOL,ICPR(ITMP,2))
  258. ELSE
  259. NBNDS=NBNDS+1
  260. ICPR(ITMP,1)=NBNDS
  261. ICPR(ITMP,2)=1
  262. NCOL=MAX(NCOL,ICPR(ITMP,2))
  263. ENDIF
  264. 10 CONTINUE
  265. IF (IPT5.LISOUS(/1).NE.0) SEGDES IPT2
  266. 11 CONTINUE
  267. IF (IDIM.EQ.2) NCOL=NCOL*2
  268. IF (IDIM.EQ.3) NCOL=NCOL*4
  269. C
  270. C
  271. C
  272. C
  273. C=======================================================================
  274. C Prise en compte des noeuds de relations
  275. C=======================================================================
  276. C IRELA(I,1) = 1 si le noeud I (en numerotation locale) est un noeud
  277. C de relation ("hanging node")
  278. C IRELA(I,1) = 0 sinon
  279. SEGINI IRELA
  280. IF (LCONF3.EQ.0) GOTO 1200
  281. SEGACT IPT6
  282. DO 1218 LAA=1,MAX(1,IPT6.LISOUS(/1))
  283. IF (IPT6.LISOUS(/1).EQ.0) THEN
  284. IPT2=IPT6
  285. ELSE
  286. IPT2=IPT6.LISOUS(LAA)
  287. ENDIF
  288. SEGACT IPT2
  289. IF (IPT2.ITYPEL.NE.48) GOTO 1218
  290. DO 1220 J=1,IPT2.NUM(/2)
  291. IRELA(ICPR(IPT2.NUM(1,J),1))=1
  292. 1220 CONTINUE
  293. IF (IPT6.LISOUS(/1).NE.0) SEGDES IPT2
  294. 1218 CONTINUE
  295. 1200 CONTINUE
  296. c write (*,*) 'jusquelatoutvasbien'
  297.  
  298. C
  299. C=======================================================================
  300. C Prise en compte du champoint de densite
  301. C=======================================================================
  302. C On affecte une densite a chacun des noeuds du maillage. On calcule
  303. C cette densite a partir des densites existantes dans le champoint, a
  304. C l'aide de la fonction RAFDEN.
  305. C Remarque : on ne prend pas en compte les densites qui pourraient
  306. C exister dans le tableau XCOOR.
  307. C
  308.  
  309. C------------------------------------------
  310. C A) Recuperation des valeurs du champoint !
  311. C------------------------------------------
  312. C XDEN(I,1) = densite affectee au noeud I (en numerotation globale)
  313. IF (NBOUCL.EQ.1) THEN
  314. SEGACT MCHPO1
  315. IF (MCHPO1.IPCHP(/1).NE.1) RETURN
  316. MSOUPO=MCHPO1.IPCHP(1)
  317. SEGACT MSOUPO
  318. IPT3=IGEOC
  319. c write(*,*) 'RAFF: #',MCHPO1,IPT3
  320. MPOVA1=IPOVAL
  321. SEGACT IPT3
  322. SEGACT MPOVA1
  323. IF (MPOVA1.VPOCHA(/2).NE.1) RETURN
  324. NBPTCH=IPT3.NUM(/2)
  325. c SEGINI ICHPO
  326. DO 13 I=1,NBPTCH
  327. J=IPT3.NUM(1,I)
  328. c XR(I)=XCOOR((J-1)*(IDIM+1)+1)
  329. c YR(I)=XCOOR((J-1)*(IDIM+1)+2)
  330. c IF (IDIM.EQ.3) ZR(I)=XCOOR((J-1)*(IDIM+1)+3)
  331. c DEN(I)=MPOVA1.VPOCHA(I,1)
  332. c if(J.eq.115.or.J.eq.321.or.
  333. c . J.eq.6261.or.J.eq.5815) then
  334. C write(*,*) 'RAFF: #',J,'DEN(',I,'/',NBPTCH,')=',DEN(I)
  335. c endif
  336. XDEN(J)=MPOVA1.VPOCHA(I,1)
  337. 13 CONTINUE
  338. C Rem ;seulement pour les noeud "peres",
  339. C pour les "fils" raff2d remplit XDEN
  340. ENDIF
  341. C
  342. C--------------------------------------------------------------
  343. C B) Affectation d'une densite a chacun des noeuds du maillage !
  344. C--------------------------------------------------------------
  345. C IDENSI(I,1) = densite affectee au noeud I (en numerotation locale)
  346. SEGINI IDENSI
  347. DO 15 INB=1,MAX(1,IPT5.LISOUS(/1))
  348. IF (IPT5.LISOUS(/1).EQ.0) THEN
  349. IPT2=IPT5
  350. ELSE
  351. IPT2=IPT5.LISOUS(INB)
  352. SEGACT IPT2
  353. ENDIF
  354. c write(*,*) 'ipt2 : ', ipt2, 'type', IPT2.ITYPEL
  355. c write(*,*) IPT2.NUM(/2), IPT2.NUM(/1)
  356. c IF (IPT2.ITYPEL.EQ.48) GOTO 15
  357. DO 14 J=1,IPT2.NUM(/2)
  358. DO 14 I=1,IPT2.NUM(/1)
  359. NGLOB=IPT2.NUM(I,J)
  360. NLOC=ICPR(NGLOB,1)
  361. c XPT=XCOOR((NGLOB-1)*(IDIM+1)+1)
  362. c YPT=XCOOR((NGLOB-1)*(IDIM+1)+2)
  363. c ZPT=0
  364. c IF (IDIM.EQ.3) ZPT=XCOOR((NGLOB-1)*(IDIM+1)+3)
  365. c CALL RAFDEN(ICHPO,XPT,YPT,ZPT,IDIM,DENPT)
  366. c DENSI(NLOC)=DENPT
  367. DENSI(NLOC)=XDEN(NGLOB)
  368. c WRITE (*,*) 'densite au pt', NGLOB, XDEN(NGLOB), XPETI
  369. IF (XDEN(NGLOB).LT.XPETI) THEN
  370. c WRITE (*,*) i,j
  371. WRITE (*,*) 'Champ de densité en entrée nul au point',
  372. + NGLOB, 'impossible de raffiner'
  373. RETURN
  374. ENDIF
  375. 14 CONTINUE
  376. IF (IPT5.LISOUS(/1).NE.0) SEGDES IPT2
  377. 15 CONTINUE
  378. C
  379. C=======================================================================
  380. C Remplissage du tableau KARETE
  381. C=======================================================================
  382. C On remplit le tableau KARETE avec les donnees provenant de l'analyse
  383. C du maillage.
  384. C KARETE(i,j)=k : le ieme noeud local forme une arete
  385. C avec le keme noeud local (k>i)
  386. C KARPOS(i)=k : compteur de position (k=nb d'aretes touchant i)
  387. C KARET2(i,j)=k : l'arete i-j appartient a k elements
  388. JPLA1=NBNDS
  389. JPLA2=NCOL
  390. SEGINI KARETE,KARPOS,KARET2
  391. NBARET=0
  392. NCOL=1
  393. DO 18 INB=1,MAX(1,IPT5.LISOUS(/1))
  394. IF (IPT5.LISOUS(/1).EQ.0) THEN
  395. IPT2=IPT5
  396. ELSE
  397. IPT2=IPT5.LISOUS(INB)
  398. SEGACT IPT2
  399. ENDIF
  400. c IF (IPT2.ITYPEL.EQ.48) GOTO 18
  401. ILPL=LPL(IPT2.ITYPEL)
  402. ILPT=LPT(IPT2.ITYPEL)
  403. DO 17 J=1,IPT2.NUM(/2)
  404. DO 17 K=1,ILPL*2-1,2
  405. NPTA=IPT2.NUM(KSEGM(ILPT+K-1),J)
  406. NPTB=IPT2.NUM(KSEGM(ILPT+K),J)
  407. NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  408. NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  409. NEXIST=0
  410. DO 16 I=1,MAX(1,KARPOS(NMIN))
  411. IF (KARETE(NMIN,I).EQ.NMAX) NEXIST=I
  412. 16 CONTINUE
  413. c création d'une nouvelle arrète
  414. IF (NEXIST.EQ.0) THEN
  415. IF ((KARPOS(NMIN)+1).GT.NCOL) THEN
  416. NCOL=KARPOS(NMIN)+1
  417. SEGADJ KARETE, KARET2
  418. ENDIF
  419. KARPOS(NMIN)=KARPOS(NMIN)+1
  420. KARETE(NMIN,KARPOS(NMIN))=NMAX
  421. KARET2(NMIN,KARPOS(NMIN))=1
  422. NBARET=NBARET+1
  423. c NCOL=MAX(NCOL,KARPOS(NMIN))
  424. ENDIF
  425. IF (NEXIST.NE.0) THEN
  426. KARET2(NMIN,NEXIST)=KARET2(NMIN,NEXIST)+1
  427. ENDIF
  428. 17 CONTINUE
  429. IF (IPT5.LISOUS(/1).NE.0) SEGDES IPT2
  430. 18 CONTINUE
  431. c SEGADJ KARETE,KARET2
  432. SEGDES KARETE,KARPOS,KARET2
  433. c write (*,*) 'KARETE REMPLI'
  434. C
  435. C=======================================================================
  436. C Remplissage des tableaux JPLANS et JNOEFA
  437. C=======================================================================
  438. C On remplit les tableaux de description des faces dans le cas 3D.
  439. C JPLANS(i,j)=k : la jeme face dont le plus petit sommet en
  440. C numerotation locale est le sommet i est la face k
  441. C JNOEFA(i,1-3)=j-k-l : Les 3 plus petits sommets en numerotation
  442. C locale de la face i sont les noeuds j,k et l
  443. C JNOEFA(i,4-5)=j-k : La face identifiee par les 3 sommets precedents
  444. C appartient aux lisous de type j et k
  445. C JPLCOM(i)=k : compteur de position (k=nb de faces touchant i)
  446. IF (IDIM.NE.3) GOTO 40
  447. C----------------------------------------------------
  448. C A) Comptage du nombre d'elements total du maillage !
  449. C----------------------------------------------------
  450. NELTOT=0
  451. DO 50 INB=1,MAX(1,IPT5.LISOUS(/1))
  452. IF (IPT5.LISOUS(/1).EQ.0) THEN
  453. IPT2=IPT5
  454. ELSE
  455. IPT2=IPT5.LISOUS(INB)
  456. SEGACT IPT2
  457. ENDIF
  458. NELTOT=NELTOT+IPT2.NUM(/2)
  459. 50 CONTINUE
  460. C
  461. C-------------------------------------------------
  462. C B) Initialisation des tableaux JPLANS et JNOEFA !
  463. C-------------------------------------------------
  464. JNBFA=NELTOT*6
  465. JMAXPL=0
  466. JNUMFA=0
  467. SEGINI JPLANS,JPLCOM,JNOEFA
  468. C
  469. C-----------------------------------------------
  470. C C) Remplissage des tableaux JPLANS et JNOEFA !
  471. C-----------------------------------------------
  472. DO 51 INB=1,MAX(1,IPT5.LISOUS(/1))
  473. IF (IPT5.LISOUS(/1).EQ.0) THEN
  474. IPT2=IPT5
  475. ELSE
  476. IPT2=IPT5.LISOUS(INB)
  477. SEGACT IPT2
  478. ENDIF
  479. C IF (IPT2.ITYPEL.EQ.48) GOTO 51
  480. DO 52 KELEM=1,IPT2.NUM(/2)
  481. JTYPEL=IPT2.ITYPEL
  482. JLTEL1=LTEL(1,JTYPEL)
  483. JLTEL2=LTEL(2,JTYPEL)
  484. DO 53 IBB=1,JLTEL1
  485. JLDEL1=LDEL(1,JLTEL2-1+IBB)
  486. JTYFAC=IWORK1(JLDEL1)
  487. JLDEL2=LDEL(2,JLTEL2-1+IBB)
  488. JNBSOM=NBSOM(JTYFAC)
  489. JSPOS=NSPOS(JTYFAC)
  490. SEGINI IWORK2
  491. DO 54 IAA=1,JNBSOM
  492. NGLOBA=IPT2.NUM(LFAC(JLDEL2-1+IBSOM(JSPOS-1+IAA)),KELEM)
  493. IWORK2(IAA)=NGLOBA
  494. 54 CONTINUE
  495. C
  496. C----------------------------------------------------
  497. C D) Tri des sommets de la face par ordre croissant !
  498. C----------------------------------------------------
  499. C On va ranger les 3 plus petits sommets de la face (en numerotation
  500. C globale) dans les variables NPTA, NPTB, NPTC. On les classe par ordre
  501. C croissant.
  502. C
  503. C 1/ Classement des 3 sommets par ordre croissant (num globale)
  504. NPTA=(XCOOR(/1)/(IDIM+1))+1
  505. NPTB=NPTA+1
  506. NPTC=NPTB+1
  507. DO 55 ICC=1,JNBSOM
  508. IF (IWORK2(ICC).LT.NPTA) THEN
  509. NPTC=NPTB
  510. NPTB=NPTA
  511. NPTA=IWORK2(ICC)
  512. ELSEIF (IWORK2(ICC).LT.NPTB) THEN
  513. NPTC=NPTB
  514. NPTB=IWORK2(ICC)
  515. ELSEIF (IWORK2(ICC).LT.NPTC) THEN
  516. NPTC=IWORK2(ICC)
  517. ENDIF
  518. 55 CONTINUE
  519. C
  520. C 2/ Classement des 3 sommets par ordre croissant (num locale)
  521. NPTA2=ICPR(NPTA,1)
  522. NPTB2=ICPR(NPTB,1)
  523. NPTC2=ICPR(NPTC,1)
  524. IF ((NPTA2.LT.NPTB2).AND.(NPTA2.LT.NPTC2)) THEN
  525. NPTA=NPTA2
  526. NPTB=MIN(NPTB2,NPTC2)
  527. NPTC=MAX(NPTB2,NPTC2)
  528. ENDIF
  529. IF ((NPTB2.LT.NPTA2).AND.(NPTB2.LT.NPTC2)) THEN
  530. NPTA=NPTB2
  531. NPTB=MIN(NPTA2,NPTC2)
  532. NPTC=MAX(NPTA2,NPTC2)
  533. ENDIF
  534. IF ((NPTC2.LT.NPTA2).AND.(NPTC2.LT.NPTB2)) THEN
  535. NPTA=NPTC2
  536. NPTB=MIN(NPTA2,NPTB2)
  537. NPTC=MAX(NPTA2,NPTB2)
  538. ENDIF
  539. C
  540. C----------------------------------------
  541. C E) Renseignement des nouvelles faces !
  542. C----------------------------------------
  543. SEGACT LTYPE
  544. c IHERED=0
  545. c DO 1022 MOI=1,IPT2.NUM(/1)
  546. c IHERED=IHERED+IRELA(ICPR(IPT2.NUM(MOI,KELEM),1))
  547. c 1022 CONTINUE
  548.  
  549. C Si NPTA Ne correspond a aucune face : On cree une face
  550. IF (JPLCOM(NPTA).EQ.0) THEN
  551. JNUMFA=JNUMFA+1
  552. JPLCOM(NPTA)=JPLCOM(NPTA)+1
  553. JPLANS(NPTA,JPLCOM(NPTA))=JNUMFA
  554. JNOEFA(JNUMFA,1)=NPTA
  555. JNOEFA(JNUMFA,2)=NPTB
  556. JNOEFA(JNUMFA,3)=NPTC
  557. JNOEFA(JNUMFA,4)=LTYPE(INB,2)
  558. c IF (IHERED.NE.0) JNOEFA(JNUMFA,5)=LTYPE(INB,2)
  559. C On signale les sous face de faces supportant un sure
  560. C JNOEFA(JNUMFA,5)= 48
  561. IF (LCONF3.NE.0) THEN
  562. SEGACT IPT6
  563. DO 501 LAA1=1,MAX(1,IPT6.LISOUS(/1))
  564. IF (IPT6.LISOUS(/1).EQ.0) THEN
  565. IPT4=IPT6
  566. ELSE
  567. IPT4=IPT6.LISOUS(LAA1)
  568. ENDIF
  569. SEGACT IPT4
  570. IF (IPT4.ITYPEL.NE.48) GOTO 501
  571. IF (IPT4.ICOLOR(1).NE.3) GOTO 501
  572. DO 502,MAA1=1,IPT4.NUM(/2)
  573. NEXI1=0
  574. NEXI2=0
  575. DO 503 ICC1=1,JNBSOM
  576. IF (IWORK2(ICC1).EQ.IPT4.NUM(1, MAA1)) THEN
  577. NEXI1=NEXI1+1
  578. ENDIF
  579. 503 CONTINUE
  580. IF (NEXI1.GT.0) THEN
  581. DO 504 ICC2=1,JNBSOM
  582. DO 504 INN2=2,IPT4.NUM(/1)
  583. IF (IWORK2(ICC2).EQ.IPT4.NUM(INN2,MAA1)) THEN
  584. NEXI2=NEXI2+1
  585. ENDIF
  586. 504 CONTINUE
  587. IF (NEXI2.GT.0) JNOEFA(JNUMFA,5)= 48
  588. ENDIF
  589. 502 CONTINUE
  590. 501 CONTINUE
  591. ENDIF
  592. c IF (npta.eq.226) THEN
  593. c write (*,*) 'construct face :', JNUMFA
  594. c write (*,*) 'JNOEFA(JNUMFA,4)',JNOEFA(JNUMFA,4)
  595. c write (*,*) 'JNOEFA(JNUMFA,5)',JNOEFA(JNUMFA,5)
  596. c
  597. C endif
  598.  
  599. C Si NPTA correspond a aumoins une face
  600. ELSE
  601. NEXIST=0
  602. DO 56 IDD=1,JPLCOM(NPTA)
  603. ITEMP1=JPLANS(NPTA,IDD)
  604. NAA=JNOEFA(ITEMP1,1)
  605. NBB=JNOEFA(ITEMP1,2)
  606. NCC=JNOEFA(ITEMP1,3)
  607. C Si la face que actuelle existe déja : on renseigne JNOEFA(NEXIST,5)
  608. IF (NAA.EQ.NPTA.AND.NBB.EQ.NPTB.AND.NCC.EQ.NPTC) THEN
  609. NEXIST=ITEMP1
  610. JNOEFA(NEXIST,5)=LTYPE(INB,2)
  611. c IF (NPTA.eq.226) THEN
  612. c write (*,*) 'on repasse sur la face ;', NEXIST
  613. c write (*,*) 'JNOEFA(NEXIST,4)',JNOEFA(NEXIST,4)
  614. c write (*,*) 'JNOEFA(NEXIST,5)',JNOEFA(NEXIST,5)
  615. c endif
  616. ENDIF
  617.  
  618.  
  619. 56 CONTINUE
  620. C Si la face que actuelle n'existe pas : On cree une face
  621. IF (NEXIST.EQ.0) THEN
  622. JNUMFA=JNUMFA+1
  623. JPLCOM(NPTA)=JPLCOM(NPTA)+1
  624. JPLANS(NPTA,JPLCOM(NPTA))=JNUMFA
  625. JNOEFA(JNUMFA,1)=NPTA
  626. JNOEFA(JNUMFA,2)=NPTB
  627. JNOEFA(JNUMFA,3)=NPTC
  628. JNOEFA(JNUMFA,4)=LTYPE(INB,2)
  629. c IF (IHERED.NE.0) JNOEFA(JNUMFA,5)=LTYPE(INB,2)
  630.  
  631. C On signale les sous face de faces supportant un sure
  632. C JNOEFA(JNUMFA,5)= 48
  633. IF (LCONF3.NE.0) THEN
  634. SEGACT IPT6
  635. DO 505 LAA1=1,MAX(1,IPT6.LISOUS(/1))
  636. IF (IPT6.LISOUS(/1).EQ.0) THEN
  637. IPT4=IPT6
  638. ELSE
  639. IPT4=IPT6.LISOUS(LAA1)
  640. ENDIF
  641. SEGACT IPT4
  642. IF (IPT4.ITYPEL.NE.48) GOTO 505
  643. IF (IPT4.ICOLOR(1).NE.3) GOTO 505
  644. DO 506,MAA1=1,IPT4.NUM(/2)
  645. NEXI1=0
  646. NEXI2=0
  647. DO 507 ICC1=1,JNBSOM
  648. IF (IWORK2(ICC1).EQ.IPT4.NUM(1, MAA1)) THEN
  649. NEXI1=NEXI1+1
  650. ENDIF
  651. 507 CONTINUE
  652. IF (NEXI1.GT.0) THEN
  653. DO 508 ICC2=1,JNBSOM
  654. DO 508 INN2=2,IPT4.NUM(/1)
  655. IF (IWORK2(ICC2).EQ.IPT4.NUM(INN2,MAA1)) THEN
  656. NEXI2=NEXI2+1
  657. ENDIF
  658. 508 CONTINUE
  659. IF (NEXI2.GT.0) JNOEFA(JNUMFA,5)= 48
  660. ENDIF
  661. 506 CONTINUE
  662. 505 CONTINUE
  663. ENDIF
  664. c IF (NPTA.eq.226) THEN
  665. c write (*,*) 'const de la face ;', JNUMFA
  666. c write (*,*) 'JNOEFA(JNUMFA,4)',JNOEFA(JNUMFA,4)
  667. c write (*,*) 'JNOEFA(JNUMFA,5)',JNOEFA(JNUMFA,5)
  668. c endif
  669.  
  670. ENDIF
  671. ENDIF
  672.  
  673. JMAXPL=MAX(JMAXPL,JPLCOM(NPTA))
  674. 53 CONTINUE
  675. 52 CONTINUE
  676. 51 CONTINUE
  677. JPLA2=JMAXPL
  678. JNBFA=JNUMFA
  679. SEGADJ JPLANS,JPLCOM,JNOEFA
  680. SEGDES LTYPE
  681. c WRITE (*,*) 'JPLAN et JNOEFA remplis'
  682. C
  683. C=======================================================================
  684. C Pre-remplissage du tableau KMILIE
  685. C=======================================================================
  686. C On remplit le tableau lorsqu'il y a des maillages de relations de
  687. C type 48 (IPT6, LCONF3 non nul)
  688. C KARETE(i,j)=k ET KMILIE(i,j)=n : le noeud global n est situe
  689. C au milieu de l'arete [i,k]
  690. 40 CONTINUE
  691. SEGINI KMILIE
  692. SEGACT KARETE*MOD,KARPOS*MOD,KARET2*MOD
  693. IF (LCONF3.EQ.0) THEN
  694. GOTO 28
  695. ENDIF
  696. 22 CONTINUE
  697. NBRELA=0
  698. SEGACT IPT6
  699. DO 218 LAA=1,MAX(1,IPT6.LISOUS(/1))
  700. IF (IPT6.LISOUS(/1).EQ.0) THEN
  701. IPT2=IPT6
  702. ELSE
  703. IPT2=IPT6.LISOUS(LAA)
  704. ENDIF
  705. SEGACT IPT2
  706. IF (IPT2.ITYPEL.NE.48) GOTO 218
  707. IF ((IPT2.ICOLOR(1).NE.1).AND.(IPT2.ICOLOR(1).NE.2)) GOTO 218
  708. DO 220 J=1,IPT2.NUM(/2)
  709. c write(*,*) 'j', j
  710. NPTA=IPT2.NUM(2,J)
  711. NPTB=IPT2.NUM(3,J)
  712. NPTC=IPT2.NUM(1,J)
  713. NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  714. NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  715. NEXIST=0
  716. DO 219 I=1,MAX(1,KARPOS(NMIN))
  717. IF (KARETE(NMIN,I).EQ.NMAX) NEXIST=I
  718. 219 CONTINUE
  719. c write (*,*) 'NEXIST', NEXIST
  720.  
  721. IF (NEXIST.NE.0) then
  722. KMILIE(NMIN,NEXIST)=NPTC
  723. c WRITE (*,*) 'REPORT ANCIENNE RELA AU NOEUD:', NPTC
  724. c write (*,*) 'karret2(nmin nexist)', KARET2(NMIN,NEXIST)
  725. KARET2(NMIN,NEXIST)=2
  726. ELSE
  727. IF (IDIM.EQ.2) goto 27
  728. IF ((KARPOS(NMIN)+1).GT.NCOL) THEN
  729. NCOL=KARPOS(NMIN)+1
  730. SEGADJ KARETE, KARET2, KMILIE
  731. ENDIF
  732. JPOSI = KARPOS(NMIN)
  733. KARPOS(NMIN)=JPOSI+1
  734. KARETE(NMIN,JPOSI+1)=NMAX
  735. KARET2(NMIN,JPOSI+1)=2
  736. KMILIE(NMIN,JPOSI+1)=NPTC
  737. c WRITE (*,*) 'création de l arrete:', NMIN, NMAX
  738. ENDIF
  739. C
  740. C-----------------------------------------------
  741. C Prise en compte des relations "hereditaires" !
  742. C-----------------------------------------------
  743. C Cette partie sert a traiter le cas ou un element serait raffine
  744. C plusieurs fois alors qu'a sa frontiere les autres elements ne le
  745. C seraient pas. Il y aurait donc plusieurs points sur une arete ayant
  746. C des relations de conformite et plus seulement le noeud milieu.
  747. C
  748. c dans l'exenple ci dessous où les noeuds C et D sont des hanging nodes
  749. c L'arette [C B] n'a pas encore été créée c'est ce qui est fait ici.
  750. c
  751. c | |
  752. c | |
  753. c |A C D |B
  754. c x---------X----X----X
  755. c | | | |
  756. c | | | |
  757. 27 CONTINUE
  758. IF (IDIM.EQ.2) THEN
  759. C 1/ Test sur la premiere arete : NMIN-NPTC
  760. IF (NEXIST.NE.0) THEN
  761. NEXIS2=0
  762. NMIN2=MIN(NMIN,ICPR(NPTC,1))
  763. NMAX2=MAX(NMIN,ICPR(NPTC,1))
  764. DO 24 I=1,MAX(1,KARPOS(NMIN2))
  765. IF (KARETE(NMIN2,I).EQ.NMAX2) NEXIS2=I
  766. 24 CONTINUE
  767. IF (NEXIS2.EQ.0) THEN
  768. NBRELA=NBRELA+1
  769. JNBCOL=KARETE(/2)
  770. JPOSI=KARPOS(NMIN2)
  771. IF (JNBCOL.GT.JPOSI) THEN
  772. KARPOS(NMIN2)=JPOSI+1
  773. KARETE(NMIN2,JPOSI+1)=NMAX2
  774. KARET2(NMIN2,JPOSI+1)=2
  775. ELSEIF (JNBCOL.EQ.JPOSI) THEN
  776. NCOL=JNBCOL+1
  777. SEGADJ KARETE,KMILIE,KARET2
  778. KARETE(NMIN2,JPOSI+1)=NMAX2
  779. KARPOS(NMIN2)=JPOSI+1
  780. KARET2(NMIN2,JPOSI+1)=2
  781. ENDIF
  782. c WRITE (*,*) 'création de l arrete:', NMIN2, NMAX2
  783. ELSE
  784. KARET2(NMIN2,NEXIS2)=2
  785. ENDIF
  786. C
  787. C 2/ Test sur la deuxieme arete : NMAX-NPTC
  788. NEXIS2=0
  789. NMIN2=MIN(NMAX,ICPR(NPTC,1))
  790. NMAX2=MAX(NMAX,ICPR(NPTC,1))
  791. DO 25 I=1,MAX(1,KARPOS(NMIN2))
  792. IF (KARETE(NMIN2,I).EQ.NMAX2) NEXIS2=I
  793. 25 CONTINUE
  794. IF (NEXIS2.EQ.0) THEN
  795. NBRELA=NBRELA+1
  796. JNBCOL=KARETE(/2)
  797. JPOSI=KARPOS(NMIN2)
  798. IF (JNBCOL.GT.JPOSI) THEN
  799. KARPOS(NMIN2)=JPOSI+1
  800. KARETE(NMIN2,JPOSI+1)=NMAX2
  801. KARET2(NMIN2,JPOSI+1)=2
  802. ELSEIF (JNBCOL.EQ.JPOSI) THEN
  803. NCOL=JNBCOL+1
  804. SEGADJ KARETE,KMILIE,KARET2
  805. KARETE(NMIN2,JPOSI+1)=NMAX2
  806. KARPOS(NMIN2)=JPOSI+1
  807. KARET2(NMIN2,JPOSI+1)=2
  808. ENDIF
  809. c WRITE (*,*) 'création de l arrete:', NMIN2, NMAX2
  810. ELSE
  811. KARET2(NMIN2,NEXIS2)=2
  812. ENDIF
  813. ENDIF
  814. ENDIF
  815. 220 CONTINUE
  816. SEGDES IPT2
  817. 218 CONTINUE
  818. IF (NBRELA.NE.0) GOTO 22
  819. 28 CONTINUE
  820. C write (*,*) 'kmilieu pre-rempli'
  821. C
  822. C
  823. C=======================================================================
  824. C Pre-remplissage du tableau JFARAF
  825. C=======================================================================
  826. C On remplit le tableau JFARAF et JPLAN3 lorsqu'il y a des maillages de
  827. C relations de type 48 (IPT6, LCONF3 non nul).
  828. c
  829. C Si JPLANS(i,j)=k et JFARAF(k , (l-1)*2+1)= n : le noeud n (en
  830. C numérotation globale) est le hanging node au centre de la jieme face
  831. C dont le plus petit sommet en numerotation locale est le sommet i.
  832. C
  833. C Si de plus la face en question est un TRI6 et JFARAF(k , l*2)=m :
  834. C le noeud m (en numérotation globale) est le 4ieme noeud de l'element
  835. C 48 coul 4
  836. C
  837. C Si de plus la face en question est un QUA8-1 et JFARAF(k , l*2)=m :
  838. C le noeud m (en numérotation globale) est le 5eme noeud de l'element
  839. c 48 coul 4
  840. C
  841. C si JPLAN3(i,j)=k il y a k relation de compatibilité de faces associé a
  842. C la jieme face dont le plus petit sommet en numerotation locale est le
  843. C sommet i.
  844.  
  845. IF (IDIM.NE.3) GOTO 711
  846. JNBFA=JNOEFA(/1)
  847. LLLL=10
  848. SEGINI JFARAF
  849. SEGINI JPLAN2,JPLAN3
  850. IF (LCONF3.EQ.0) GOTO 711
  851. MPOS25=0
  852. JNBSOM=0
  853. 1287 CONTINUE
  854. NBRELA=0
  855. DO 670 MIE=1,MAX(1,IPT6.LISOUS(/1))
  856. IF (IPT6.LISOUS(/1).EQ.0) THEN
  857. IPT2=IPT6
  858. ELSE
  859. IPT2=IPT6.LISOUS(MIE)
  860. ENDIF
  861. SEGACT IPT2
  862. IF (IPT2.ITYPEL.NE.48) GOTO 670
  863. IF ((IPT2.ICOLOR(1).EQ.1).OR.(IPT2.ICOLOR(1).EQ.2)) GOTO 670
  864. LIGNUM=IPT2.NUM(/1)
  865. IF (IPT2.ICOLOR(1).EQ.4) MPOS25=4
  866. IF (IPT2.ICOLOR(1).EQ.5) MPOS25=5
  867. IF (IPT2.ICOLOR(1).EQ.3) JNBSOM=4
  868. IF (IPT2.ICOLOR(1).EQ.4) JNBSOM=3
  869. IF (IPT2.ICOLOR(1).EQ.5) JNBSOM=4
  870. IF (IPT2.ICOLOR(1).EQ.6) JNBSOM=4
  871. SEGINI IWORK2
  872. C
  873. C-------------------------
  874. C A) Recherche de la face !
  875. C-------------------------
  876. C 1/ Classement des 3 sommets par ordre croissant (num globale)
  877. DO 671,MAA=1,IPT2.NUM(/2)
  878. c write (*,*) 'MAA', MAA
  879. DO 673 MBB=1,JNBSOM
  880. IWORK2(MBB)=IPT2.NUM(LIGNUM-JNBSOM+MBB,MAA)
  881. 673 CONTINUE
  882. NPTA=(XCOOR(/1)/(IDIM+1))+1
  883. NPTB=NPTA+1
  884. NPTC=NPTB+1
  885. DO 672 ICC=1,JNBSOM
  886. IF (IWORK2(ICC).LT.NPTA) THEN
  887. NPTC=NPTB
  888. NPTB=NPTA
  889. NPTA=IWORK2(ICC)
  890. ELSEIF (IWORK2(ICC).LT.NPTB) THEN
  891. NPTC=NPTB
  892. NPTB=IWORK2(ICC)
  893. ELSEIF (IWORK2(ICC).LT.NPTC) THEN
  894. NPTC=IWORK2(ICC)
  895. ENDIF
  896. 672 CONTINUE
  897. C
  898. C 2/ Classement des 3 sommets par ordre croissant (num locale)
  899. NPTA1=NPTA
  900. NPTB1=NPTB
  901. NPTC1=NPTC
  902. NPTA2=ICPR(NPTA,1)
  903. NPTB2=ICPR(NPTB,1)
  904. NPTC2=ICPR(NPTC,1)
  905. IF ((NPTA2.LT.NPTB2).AND.(NPTA2.LT.NPTC2)) THEN
  906. NPTA=NPTA2
  907. NPTB=MIN(NPTB2,NPTC2)
  908. NPTC=MAX(NPTB2,NPTC2)
  909. ENDIF
  910. IF ((NPTB2.LT.NPTA2).AND.(NPTB2.LT.NPTC2)) THEN
  911. NPTA=NPTB2
  912. NPTB=MIN(NPTA2,NPTC2)
  913. NPTC=MAX(NPTA2,NPTC2)
  914. ENDIF
  915. IF ((NPTC2.LT.NPTA2).AND.(NPTC2.LT.NPTB2)) THEN
  916. NPTA=NPTC2
  917. NPTB=MIN(NPTA2,NPTB2)
  918. NPTC=MAX(NPTA2,NPTB2)
  919. ENDIF
  920. C
  921. C 3/ Recherche du numero de la face
  922. NEXIS2=0
  923. DO 674 IEE=1,JPLCOM(NPTA)
  924. MTMP=JPLANS(NPTA,IEE)
  925. JJ1=JNOEFA(MTMP,1)
  926. JJ2=JNOEFA(MTMP,2)
  927. JJ3=JNOEFA(MTMP,3)
  928. IF(JJ1.EQ.NPTA.AND.JJ2.EQ.NPTB.AND.JJ3.EQ.NPTC) THEN
  929. NEXIS2=IEE
  930. ENDIF
  931. 674 CONTINUE
  932. c 3.1/ si la face en question n'existe pas, dans le cas général on
  933. c s'arrète
  934. IF ((NEXIS2.EQ.0).AND.(IPT2.ICOLOR(1).NE.3)) THEN
  935. c WRITE (*,*) 'IPT2.ICOLOR(1)', IPT2.ICOLOR(1)
  936. WRITE (*,*) 'Raffinement impossible, le champ de densite
  937. + presente des variations trop brusques'
  938. WRITE (*,*) 'sure concerne', MAA
  939. WRITE (*,*) 'premier nd du sure concerne', IPT2.NUM(1,MAA)
  940. WRITE (*,*) 'numero du 1er nd de la face concernee', NPTA
  941. KARAF2=0
  942. c GOTO 679
  943. ENDIF
  944. c 3.2/ si la face en question n'existe pas dans le cas d'une relation
  945. c de face QUA4 il suffie de la créer
  946.  
  947. IF ((NEXIS2.EQ.0).AND.(IPT2.ICOLOR(1).EQ.3)) THEN
  948. JPLCOM(NPTA)= JPLCOM(NPTA)+1
  949. JNUMFA=JNUMFA+1
  950. JNUMFS=JNUMFA
  951. JNUMF1=JPLCOM(NPTA)
  952. JMAXPL=MAX(JMAXPL,JPLCOM(NPTA))
  953. c 3.2.1/ ajustement des segments
  954. JPLA2=JMAXPL
  955. JNBFA=JNUMFA
  956. SEGADJ JPLANS
  957. SEGADJ JNOEFA
  958. SEGADJ JPLAN3
  959. SEGADJ JFARAF
  960. C 3.2.1 Renseignement des iformation de la nouvelle face dans JPLANS et
  961. c JNOEFA, initialisation de JPLAN3
  962. JPLANS(NPTA,JNUMF1)=JNUMFA
  963. JPLAN3(NPTA,JNUMF1)=0
  964. JNOEFA(JNUMFA,1)= NPTA
  965. JNOEFA(JNUMFA,2)= NPTB
  966. JNOEFA(JNUMFA,3)= NPTC
  967. JNOEFA(JNUMFA,4)= 48
  968. JNOEFA(JNUMFA,5)= 48
  969. ENDIF
  970.  
  971. IF (NEXIS2.NE.0) THEN
  972. JNUMFS=JPLANS(NPTA,NEXIS2)
  973. ENDIF
  974. C
  975. C 4/ Recherche de la position de la face dans JPLANS
  976. NEXIS2=0
  977. DO 675 MAC=1,JPLCOM(NPTA)
  978. IF (JPLANS(NPTA,MAC).EQ.JNUMFS) NEXIS2=MAC
  979. 675 CONTINUE
  980. C
  981. C-------------------------------------
  982. C B) Renseignement de JFARAF et JPLAN3!
  983. C-------------------------------------
  984.  
  985.  
  986. JPLAN3(NPTA,NEXIS2)=JPLAN3(NPTA,NEXIS2)+1
  987. LTEMPO=(JPLAN3(NPTA,NEXIS2)-1)*2+1
  988. JFARAF(JNUMFS,LTEMPO)=IPT2.NUM(1,MAA)
  989. IF (MPOS25.NE.0) JFARAF(JNUMFS,LTEMPO+1)=IPT2.NUM(MPOS25,MAA)
  990. IF (JNOEFA(JNUMFS,5).EQ.0) JNOEFA(JNUMFS,5)=48
  991. c IF ((MAA.eq.102).or.(MAA.eq.66)) THEN
  992. c WRITE(*,*) 'JFARAF(,)', JFARAF(JNUMFS,LTEMPO)
  993. c WRITE(*,*) 'NPTA', NPTA1, 'NPTB', NPTB1, 'NPTC', NPTC1
  994. c WRITE(*,*) 'NPTA', NPTA, 'NPTB', NPTB, 'NPTC', NPTC
  995. c WRITE(*,*) 'JNUMFS', JNUMFS, 'NEXI', NEXIS2
  996. c WRITE(*,*) 'JPLAN3(NPTA,NEXIS2)', JPLAN3(NPTA,NEXIS2)
  997. c WRITE(*,*) 'JFARAF(JNUMFS,1)', JFARAF(JNUMFS,1)
  998. c WRITE(*,*) 'JNOEFA(JNUMFS,1)', JNOEFA(JNUMFS,1)
  999. c WRITE(*,*) 'JNOEFA(JNUMFS,2)', JNOEFA(JNUMFS,2)
  1000. c WRITE(*,*) 'JNOEFA(JNUMFS,3)', JNOEFA(JNUMFS,3)
  1001. c WRITE(*,*) 'JNOEFA(JNUMFS,4)', JNOEFA(JNUMFS,4)
  1002. c WRITE(*,*) 'JNOEFA(JNUMFS,5)', JNOEFA(JNUMFS,5)
  1003. c ENDIF
  1004. C
  1005.  
  1006.  
  1007. 671 CONTINUE
  1008. 670 CONTINUE
  1009. C IF (NBRELA.NE.0) GOTO 1287
  1010. 711 CONTINUE
  1011. C
  1012. C=======================================================================
  1013. C Calcul du critere et raffinement si necessaire
  1014. C=======================================================================
  1015. C Pour calculer le critere, on a besoin d'un modele et d'un chamelem
  1016. C constant egal a 1.
  1017. C
  1018. C------------------------------------------
  1019. C A) Creation d'un modele simple : MMODE1 !
  1020. C------------------------------------------
  1021. MN3=0
  1022. NFOR=1
  1023. NMAT=2
  1024. NPARMO=0
  1025. NOBMOD=0
  1026. N1=1
  1027. SEGINI MMODE1
  1028. SEGINI IMODE1
  1029. MMODE1.KMODEL(1)=IMODE1
  1030. IMODE1.FORMOD(1)='MECANIQUE'
  1031. IMODE1.MATMOD(1)='ELASTIQUE'
  1032. IMODE1.MATMOD(2)='ISOTROPE'
  1033. SEGDES IMODE1
  1034.  
  1035. C
  1036. KARAF2=0
  1037. c WRITE (*,*) 'nb lissou ipt5', IPT5.LISOUS(/1)
  1038. C=======================================================================
  1039. C Boucle sur les eventuels sous-domaines de IPT5
  1040. C=======================================================================
  1041.  
  1042. c-----------------------------------------------------------------------
  1043. c initialisation de LTYPE3 qui va contenir les tableau LTYP2
  1044. c correspondant à chaque sous domaine de IPT5. Ces tableau LTYP2
  1045. c contiendront eux même les maillages raffinés issue des
  1046. c différents sous domaine de IPT5.
  1047. c LTYPE3(I,1)=LTYP2 pour le sous domaine I de IPT5
  1048. c LTYPE3(I,2)=IRR le permier indice d'un maillage 48 dans LTYP2
  1049. c-----------------------------------------------------------------------
  1050. NBLIS3=MAX(1,IPT5.LISOUS(/1))
  1051. SEGINI LTYPE3
  1052. MCONF=0
  1053. NSUR = 0
  1054. DO 33 INB=1,MAX(1,IPT5.LISOUS(/1))
  1055. C WRITE (*,*) 'boucle 33 ', inb
  1056. NBNN=0
  1057. NBELEM=0
  1058. NBREF=0
  1059. NBSOUS=0
  1060. SEGINI IPT2
  1061. IF (IPT5.LISOUS(/1).EQ.0) THEN
  1062. IPT2=IPT5
  1063. ELSE
  1064. IPT2=IPT5.LISOUS(INB)
  1065. SEGACT IPT2
  1066. ENDIF
  1067. C WRITE (*,*) 'TYPE SOUS MAIL', IPT2.ITYPEL
  1068.  
  1069. SEGACT IMODE1*MOD
  1070. IMODE1.IMAMOD=IPT2
  1071. IMODE1.NEFMOD=IPT2.ITYPEL
  1072. SEGDES IMODE1
  1073. C
  1074. C-------------------------------------------------------
  1075. C B) Creation d'un chamelem constant egal a 1 : MCHEL1 !
  1076. C-------------------------------------------------------
  1077. CALL ZEROP(MMODE1,'NOEUD ',MCHEL1)
  1078. SEGACT MCHEL1*MOD
  1079. MCHEL1.IMACHE(1)=IPT2
  1080. IF (IDIM.EQ.2) MCHEL1.IFOCHE=-1
  1081. IF (IDIM.EQ.3) MCHEL1.IFOCHE=2
  1082. MCHAM1=MCHEL1.ICHAML(1)
  1083. SEGACT MCHAM1
  1084. MELVA1=MCHAM1.IELVAL(1)
  1085. SEGACT MELVA1*MOD
  1086. c DO 399 J=1,IPT2.NUM(/2)
  1087. c MELVA1.VELCHE(1,J)=1.0
  1088. c 399 CONTINUE
  1089. MELVA1.VELCHE(1,1)=1.
  1090. c write (*,*) 'support Mchel1', MCHEL1.IMACHE(1)
  1091. c write (*,*) 'nb elem', IPT2.NUM(/2)
  1092. c write (*,*) 'nb pts par elem', IPT2.NUM(/1)
  1093. c write (*,*) 'taille M1vel', Melva1.Velche(/1), Melva1.Velche(/2)
  1094. c write (*,*) 'MCHEL1.ICHAML(/1)',MCHEL1.ICHAML(/1)
  1095. c write (*,*) 'MCHEL1.INFCHE(1,6)',MCHEL1.INFCHE(1,6)
  1096. c write (*,*) 'MCHAM1.IELVAL(/1)',MCHAM1.IELVAL(/1)
  1097. c write (*,*) 'MCHAM1.NOMCHE(1)',MCHAM1.NOMCHE(1)
  1098. c write (*,*) 'MCHAM1.TYPCHE(1)',MCHAM1.TYPCHE(1)
  1099. SEGDES MELVA1
  1100. SEGDES MCHAM1
  1101. C
  1102. C-------------------------------------------------------------
  1103. C C) Calcul de la surface ou du volume des elements : MCHEL2 !
  1104. C-------------------------------------------------------------
  1105. CALL INTGCA(MMODE1,MCHEL1,0,1,IRET,XRET,MCHEL2)
  1106. IF (XRET.LT.XPETI) THEN
  1107. WRITE (*,*) 'Volume(resp surf) de la zone a raffiner nul'
  1108. RETURN
  1109. ENDIF
  1110. SEGACT IPT2
  1111. c write (*,*) 'iret ', iret, 'xret', xret
  1112. c write (*,*) 'ipt2', ipt2
  1113.  
  1114. c write (*,*) 'support Mchel2', MCHEL2.IMACHE(1)
  1115.  
  1116. C
  1117. C---------------------------------------------------------------------
  1118. C D) Creation d'un chamelem indiquant les elements a raffiner MCHAM2 !
  1119. C---------------------------------------------------------------------
  1120. SEGACT KMILIE*MOD
  1121. SEGACT MCHEL2
  1122. MCHAM2=MCHEL2.ICHAML(1)
  1123. SEGACT MCHAM2
  1124. MELVA2=MCHAM2.IELVAL(1)
  1125. SEGACT MELVA2*MOD
  1126. NACREE=0
  1127. KARAF=0
  1128.  
  1129.  
  1130. SEGACT MELVA1
  1131.  
  1132. c write (*,*) 'taille M2vel', Melva2.Velche(/1), Melva2.Velche(/2)
  1133. c write (*,*) 'MCHEL2.ICHAML(/1)',MCHEL2.ICHAML(/1)
  1134. c write (*,*) 'MCHAM2.IELVAL(/1)',MCHAM2.IELVAL(/1)
  1135. c write (*,*) 'MCHAM2.NOMCHE(1)',MCHAM2.NOMCHE(1)
  1136. c write (*,*) 'MCHEL2.INFCHE(1,6)',MCHEL2.INFCHE(1,6)
  1137. c write (*,*) 'MCHAM2.TYPCHE(1)',MCHAM2.TYPCHE(1)
  1138. c WRITE (*,*) 'MELVA1.VELCHE(1,',1,')',MELVA1.VELCHE(1,1)
  1139. c DO 39 J=1,IPT2.NUM(/2)
  1140. c WRITE (*,*) 'MELVA2.VELCHE(1,',J,')',MELVA2.VELCHE(1,J)
  1141. c WRITE (*,*) 'MELVA2.VELCHE(1,',J,')',MELVA2.VELCHE(1,J)
  1142. c 39 CONTINUE
  1143. c WRITE (*,*) 'Karaf ', KARAF, 'Karaf2 ', KARAF2
  1144.  
  1145. c boucles sur les éléments de ipt2
  1146. DO 32 J=1,IPT2.NUM(/2)
  1147. C write (*,*) J
  1148. XMOY=0
  1149. DO 29 I=1,IPT2.NUM(/1)
  1150. c write (*,*) I, J
  1151. NGLOB=IPT2.NUM(I,J)
  1152. NLOC=ICPR(NGLOB,1)
  1153. DENPT=DENSI(NLOC)
  1154. XMOY=XMOY+DENPT/(IPT2.NUM(/1))
  1155. 29 CONTINUE
  1156.  
  1157. XINT=MELVA2.VELCHE(1,J)
  1158. XDIM=IDIM*1.0
  1159. C Un element est a raffiner si la moyenne des densites affectees a
  1160. C l'ensemble de ses noeuds est inférieure a la racine carree de sa
  1161. C surface (en 2D) ou a la racine cubique de son volume (en 3D)
  1162. C WRITE (*,*) 'XMOY', XMOY, '(XINT**(1/XDIM))',(XINT**(1/XDIM))
  1163. IF ((XMOY-(XINT**(1/XDIM))).LT.0) THEN
  1164. c c y a t'il des +epsilon...?
  1165. c XDIFF = XMOY-(XINT**(1/XDIM))
  1166. c if(XDIFF.ge.0.and.XDIFF.lt.0.001*XMOY) then
  1167. c write(*,*) 'RAFF : XMOY,XDIFF=',XMOY,XDIFF,
  1168. c & 'pour l EF',J,(IPT2.NUM(iou,J),iou=1,IPT2.NUM(/1))
  1169. c endif
  1170. c IF (XDIFF.LT.0) THEN
  1171. MELVA2.VELCHE(1,J)=1
  1172. c KARAF nombre d'éléments à rafiner dans le sous domaine
  1173. KARAF=KARAF+1
  1174. ILPT=LPT(IPT2.ITYPEL)
  1175. ILPL=LPL(IPT2.ITYPEL)
  1176.  
  1177. c boucle sur les segments de l'élément à raffiner
  1178. DO 31 K=1,ILPL*2-1,2
  1179. NPTA=IPT2.NUM(KSEGM(ILPT+K-1),J)
  1180. NPTB=IPT2.NUM(KSEGM(ILPT+K),J)
  1181. NMIN=MIN(ICPR(NPTA,1),ICPR(NPTB,1))
  1182. NMAX=MAX(ICPR(NPTA,1),ICPR(NPTB,1))
  1183. NEXIST=0
  1184. DO 30 I=1,MAX(1,KARPOS(NMIN))
  1185. IF (KARETE(NMIN,I).EQ.NMAX) NEXIST=I
  1186. 30 CONTINUE
  1187. c On met dans kmilie les arêtes où il faudra mettre un nouveau noeud
  1188. c On leur associe la valeur -1
  1189. IF (KMILIE(NMIN,NEXIST).EQ.0) THEN
  1190. NACREE=NACREE+1
  1191. KMILIE(NMIN,NEXIST)=-1
  1192. ENDIF
  1193. 31 CONTINUE
  1194. C
  1195. C-------------------------------------------
  1196. C E) Comptage du nombre de faces a raffiner !
  1197. C-------------------------------------------
  1198. C 1/ Initialisations
  1199. IF (IDIM.NE.3) GOTO 32
  1200. JTYPEL=IPT2.ITYPEL
  1201. JLTEL1=LTEL(1,JTYPEL)
  1202. JLTEL2=LTEL(2,JTYPEL)
  1203. DO 73 IBB=1,JLTEL1
  1204. JLDEL1=LDEL(1,JLTEL2-1+IBB)
  1205. JTYFAC=IWORK1(JLDEL1)
  1206. JLDEL2=LDEL(2,JLTEL2-1+IBB)
  1207. C JNBPTS=NBNNE(JTYFAC)
  1208. JNBSOM=NBSOM(JTYFAC)
  1209. JSPOS=NSPOS(JTYFAC)
  1210. SEGINI IWORK2
  1211. DO 74 IAA=1,JNBSOM
  1212. NGLOBA=IPT2.NUM(LFAC(JLDEL2-1+IBSOM(JSPOS-1+IAA)),J)
  1213. IWORK2(IAA)=NGLOBA
  1214. 74 CONTINUE
  1215. C
  1216. C 2/ Classement des 3 sommets par ordre croissant (num globale)
  1217. NPTA=(XCOOR(/1)/(IDIM+1))+1
  1218. NPTB=NPTA+1
  1219. NPTC=NPTB+1
  1220. DO 75 ICC=1,JNBSOM
  1221. IF (IWORK2(ICC).LT.NPTA) THEN
  1222. NPTC=NPTB
  1223. NPTB=NPTA
  1224. NPTA=IWORK2(ICC)
  1225. ELSEIF (IWORK2(ICC).LT.NPTB) THEN
  1226. NPTC=NPTB
  1227. NPTB=IWORK2(ICC)
  1228. ELSEIF (IWORK2(ICC).LT.NPTC) THEN
  1229. NPTC=IWORK2(ICC)
  1230. ENDIF
  1231. 75 CONTINUE
  1232. C
  1233. C 3/ Classement des 3 sommets par ordre croissant (num locale)
  1234. NPTA2=ICPR(NPTA,1)
  1235. NPTB2=ICPR(NPTB,1)
  1236. NPTC2=ICPR(NPTC,1)
  1237. IF ((NPTA2.LT.NPTB2).AND.(NPTA2.LT.NPTC2)) THEN
  1238. NPTA=NPTA2
  1239. NPTB=MIN(NPTB2,NPTC2)
  1240. NPTC=MAX(NPTB2,NPTC2)
  1241. ENDIF
  1242. IF ((NPTB2.LT.NPTA2).AND.(NPTB2.LT.NPTC2)) THEN
  1243. NPTA=NPTB2
  1244. NPTB=MIN(NPTA2,NPTC2)
  1245. NPTC=MAX(NPTA2,NPTC2)
  1246. ENDIF
  1247. IF ((NPTC2.LT.NPTA2).AND.(NPTC2.LT.NPTB2)) THEN
  1248. NPTA=NPTC2
  1249. NPTB=MIN(NPTA2,NPTB2)
  1250. NPTC=MAX(NPTA2,NPTB2)
  1251. ENDIF
  1252. C
  1253. C 4/ Calcul du nombre de points a creer pour raffiner la face
  1254. NEXIS2=0
  1255. DO 76 IEE=1,JPLCOM(NPTA)
  1256. MTMP=JPLANS(NPTA,IEE)
  1257. JJ1=JNOEFA(MTMP,1)
  1258. JJ2=JNOEFA(MTMP,2)
  1259. JJ3=JNOEFA(MTMP,3)
  1260. IF(JJ1.EQ.NPTA.AND.JJ2.EQ.NPTB.AND.JJ3.EQ.NPTC) THEN
  1261. NEXIS2=IEE
  1262. ENDIF
  1263. 76 CONTINUE
  1264. IF (NEXIS2.NE.0) THEN
  1265. KFARAF=JPLAN2(NPTA,NEXIS2)
  1266. IF (KFARAF.EQ.0) THEN
  1267. NACREE=NACREE+NBINTE(JTYFAC)
  1268. JPLAN2(NPTA,NEXIS2)=JPLAN2(NPTA,NEXIS2)+1
  1269. ENDIF
  1270. ENDIF
  1271. 73 CONTINUE
  1272. ELSE
  1273.  
  1274. MELVA2.VELCHE(1,J)=0
  1275. ENDIF
  1276. 32 CONTINUE
  1277.  
  1278. c KARAF2 nombre d'éléments à rafiner sur tout les sous domaines
  1279. KARAF2=KARAF2+KARAF
  1280. c WRITE (*,*) 'Karaf ', KARAF, 'Karaf2 ', KARAF2
  1281. c si aucun élément à raffiner on stoque directement IPT2 dans LTYPE3
  1282. c on stoque également les anciennes relations et on passe au
  1283. c sous-domaine suivant.
  1284.  
  1285. c ATTENTION SI NBOUCL=1 il est possible que des relations surabondantes
  1286.  
  1287. c !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  1288. c forcage : nb max itération de raff
  1289. c IF (NBOUCL.EQ.3) THEN
  1290. c KARAF2=0
  1291. c KARAF=0
  1292. c ENDIF
  1293. c !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  1294. IF (KARAF.EQ.0) THEN
  1295. IF (NBOUCL.EQ.1) THEN
  1296. IF (LCONF3.GT.0) THEN
  1297. NBLIS2=1+MAX(1,IPT6.LISOUS(/1))
  1298. SEGINI LTYPE2
  1299.  
  1300. LTYPE2(1,1)=1
  1301. LTYPE2(1,2)=IPT2.ITYPEL
  1302. LTYPE2(1,3)=IPT2
  1303. LTYPE2(1,4)=IPT2.ICOLOR(1)
  1304. LTYPE3(INB,2)=0
  1305. DO 676 MIE=1,MAX(1,IPT6.LISOUS(/1))
  1306. IF (IPT6.LISOUS(/1).EQ.0) THEN
  1307. IPT4=IPT6
  1308. ELSE
  1309. IPT4=IPT6.LISOUS(MIE)
  1310. ENDIF
  1311. NSUR=1
  1312. SEGACT IPT4
  1313. LTYPE2(1+MIE,1)=1+MIE
  1314. LTYPE2(1+MIE,2)=IPT4.ITYPEL
  1315. LTYPE2(1+MIE,3)=IPT4
  1316. LTYPE2(1+MIE,4)=IPT4.ICOLOR(1)
  1317. LTYPE3(INB,2)=2
  1318. 676 CONTINUE
  1319. ELSE
  1320. NBLIS2=1
  1321. SEGINI LTYPE2
  1322. LTYPE2(1,1)=1
  1323. LTYPE2(1,2)=IPT2.ITYPEL
  1324. LTYPE2(1,3)=IPT2
  1325. LTYPE2(1,4)=IPT2.ICOLOR(1)
  1326. LTYPE3(INB,2)=0
  1327. ENDIF
  1328. LTYPE3(INB,1)=LTYPE2
  1329. ELSE
  1330. SEGACT LTYPE4
  1331. LTYPE3(INB,1)=LTYPE4(INB,1)
  1332. LTYPE3(INB,2)=LTYPE4(INB,2)
  1333. c write (*,*) 'LTYPE4(INB,:)',LTYPE4(INB,1), LTYPE4(INB,2)
  1334. ENDIF
  1335. c WRITE (*,*) 'KARAF=0 ','ip5 ip6 ',LTYPE2(1,3), LTYPE2(2,3)
  1336. c write (*,*) 'LTYPE3(INB,:)', LTYPE3(INB,1), LTYPE3(INB,2)
  1337. GOTO 33
  1338. ENDIF
  1339. C
  1340. C----------------------------
  1341. C F) Raffinement du maillage !
  1342. C----------------------------
  1343. NML=2
  1344. C
  1345. C 1/ Raffinement d'un maillage 2D
  1346. C On recupere dans MELEME le maillage raffine sans ses relations de
  1347. C conformite et dans IPT8 les relations de ce maillage
  1348. C (IPT8 encore incomplet en sortie de RAFF2D)
  1349. c entrees de raff2D:
  1350. c - IPT2: maillage élémentaire à raffiner
  1351. c (ici sous domaine de ipt5)
  1352. c - ICPR: tableau de passage noeuds local/global
  1353. c (ici à partir de ipt5 interne à la boucle 9)
  1354. c - KARPOS: tableau du nb d'arretes par noeuds
  1355. c (ici à partir de ipt5, mis a jour avec ipt6)
  1356. c - KARETE: tableau des d'arretes KARETE(i,j)=k
  1357. c (ici à partir de ipt5, mis a jour avec ipt6)
  1358. c - KMILIE: tableau des hanging nodes associés au arrètes
  1359. c en entrée: à partir des relations de ipt6 et
  1360. c des nouveau noeuds à construires
  1361. c si relation issue de ipt6
  1362. c KMILIE(i,j)=n : le noeud global n est situe
  1363. c au milieu de l'arete [i,k]
  1364. c si nouvelle relation
  1365. c KMILIE(i,j)=-1 : il faut construire un noeud
  1366. c au milieu de l'arete [i,k]
  1367. c - MELVA2: champ par elem valait 1 pour les élem à raffiner
  1368. c et 0 pour les autres
  1369. c (ici à partir de ipt2 et de densi)
  1370. c - NACREE: nombre de noeuds à créer dans ipt2
  1371. c - KARAF: nombre d'élément à raffiner dans ipt2
  1372. c - IDCP : ??? variable non utilisé ???
  1373. c - XDEN : densité aux noeuds en notation globale
  1374. c (ici à partir de mchpo1)
  1375. c - KARET2: tableau du nombre d'éléments touchant les arretes
  1376. c (ici à partir de ipt5, mis a jour avec ipt6)
  1377. csorties de raff2D :
  1378. c - KMILIE: tableau des hanging nodes associés au arrètes
  1379. c en sortie: à partir des relations donnée en
  1380. c entrée et de celles crées
  1381. c si anciene relation à garder
  1382. c KMILIE(i,j)=n : le noeud global n est situe
  1383. c au milieu de l'arete [i,k]
  1384. c si anciene relation à suprimer
  1385. c KMILIE(i,j)=0 :
  1386. c si nouvelle relation
  1387. c KMILIE(i,j)=n : le noeud n à été construit au
  1388. c milieu de l'arete [i,k]
  1389. c - MELEME: Maillage élémentaire raffiné à partir de ipt2
  1390. c sans relations
  1391. c - IPT8 : Maillage de relation incomplet
  1392. c si J nouvel élément 48:
  1393. c IPT8.NUM(1,j) = nouveau hanging node
  1394. c IPT8.NUM(2:4,j) = autres noeuds de la relation
  1395. c IPT8.NUM(5,j) = 1 si arête a 2 noeuds
  1396. c 2 si arête a 3 noeuds
  1397. c si J nouvel élément 48:
  1398. c IPT8.NUM(1,j) = nouveau hanging node
  1399. c IPT8.NUM(2:5,j) = 0
  1400. c - XDEN : densité aux noeuds en notation globale
  1401. c interpolé aux nouveaux noeuds.
  1402. IF (IDIM.EQ.2) CALL RAFF2D(IPT2,ICPR,KARPOS,KARETE,KMILIE,
  1403. + MELVA2,NACREE,KARAF,MELEME,IDCP,IPT8,KARET2,XDEN)
  1404. C
  1405. C 2/ Raffinement d'un maillage 3D
  1406. IF (IDIM.EQ.3) THEN
  1407. NML=0
  1408. IF ((IPT2.ITYPEL.EQ.25).OR.(IPT2.ITYPEL.EQ.26)) NML=1
  1409. ENDIF
  1410. C a) Maillage 3D contenant des pyramides
  1411. IPT9=0
  1412. IF (NML.EQ.1) CALL RAFPYR(IPT2,ICPR,KARPOS,KARETE,KMILIE,
  1413. + MELVA2,NACREE,KARAF,MELEME,JPLANS,JPLAN3,JPLCOM,JNOEFA,
  1414. + IPT8,JFARAF,KARET2,IPT9,XDEN)
  1415. C b) Maillage 3D sans pyramides
  1416. IF (NML.EQ.0) CALL RAFF3D(IPT2,IPT6,ICPR,KARPOS,KARETE,KMILIE,
  1417. + MELVA2,NACREE,KARAF,MELEME,JPLANS,JPLAN3,JPLCOM,JNOEFA,
  1418. + IPT8,JFARAF,KARET2,XDEN)
  1419. C
  1420. C=======================================================================
  1421. C Creation du maillage de relations de conformite
  1422. C=======================================================================
  1423. C--------------------
  1424. C A) Initialisations !
  1425. C--------------------
  1426. IPT10=0
  1427. IPT11=0
  1428. IPT12=0
  1429. IPT13=0
  1430. IPT14=0
  1431. IPT15=0
  1432. IPT16=0
  1433. MCONF=0
  1434. c WRITE (*,*) 'Karaf ', KARAF, 'Karaf2 ', KARAF2
  1435.  
  1436. c IF (KARAF2.EQ.0) GOTO 33
  1437. IF (IPT8.NE.0) THEN
  1438. C
  1439. C---------------------------------------
  1440. C B) Traitement des anciennes relations !
  1441. C---------------------------------------
  1442. C A ce stade, IPT8 contient les relations crees lors du dernier
  1443. C raffinement (noeud de relation + noeuds formant la relation) et les
  1444. C anciennes relations a reporter (noeud de relation seulement)
  1445.  
  1446. C 1/ Recherche de l'existence d'anciennes relations a reporter
  1447. NB10=0
  1448. MLONG=IPT8.NUM(/1)
  1449. DO 775 MM=1,IPT8.NUM(/2)
  1450. IF (IPT8.NUM(MLONG,MM).EQ.0) NB10=NB10+1
  1451. 775 CONTINUE
  1452.  
  1453. IF (NB10.NE.0) THEN
  1454. C
  1455. C 2/ Creation de IPT10
  1456. c MCOL=0
  1457. c NBNN=9
  1458. c NBELEM=NB10
  1459. c NBREF=0
  1460. c NBSOUS=0
  1461. c SEGINI IPT10
  1462. c IPT10.ITYPEL=48
  1463. C
  1464. C 3/ Recherche des anciennes relations a reporter
  1465. C Cette partie permet de compléter IPT8 avec les noeuds formant les
  1466. C relations, pour les anciennes relations à reporter
  1467. DO 628 MM11=1,IPT8.NUM(/2)
  1468. IF (IPT8.NUM(MLONG,MM11).NE.0) GOTO 628
  1469. INRELA=IPT8.NUM(1,MM11)
  1470. c On active pous le cas ou IPT6=IPT2.
  1471. SEGACT IPT6
  1472. DO 625 KOU=1,MAX(1,IPT6.LISOUS(/1))
  1473. IF (IPT6.LISOUS(/1).EQ.0) THEN
  1474. IPT2=IPT6
  1475. ELSE
  1476. IPT2=IPT6.LISOUS(KOU)
  1477. ENDIF
  1478. SEGACT IPT2
  1479. IF (IPT2.ITYPEL.NE.48) GOTO 625
  1480. DO 626 JKG=1,IPT2.NUM(/2)
  1481. INOEGL=IPT2.NUM(1,JKG)
  1482. IF (INOEGL.NE.INRELA) GOTO 626
  1483. MCOL=MCOL+1
  1484. DO 627 MTMP=1,IPT2.NUM(/1)
  1485. IPT8.NUM(MTMP,MM11)=IPT2.NUM(MTMP,JKG)
  1486. IPT8.NUM(MLONG,MM11)= IPT2.ICOLOR(JKG)
  1487. 627 CONTINUE
  1488. GOTO 628
  1489. 626 CONTINUE
  1490. 625 CONTINUE
  1491. 628 CONTINUE
  1492. ENDIF
  1493. c à ce stade le maillage de relation ipt8 est complet
  1494. C
  1495. C---------------------------------------
  1496. C C) Traitement des nouvelles relations !
  1497. C---------------------------------------
  1498. C Ici, IPT8 contient toutes les relations du nouveau maillage
  1499. C On les trie par type de relations
  1500.  
  1501. C 1/ Recherche de l'existence de nouvelles relations
  1502. NB11=0
  1503. NB12=0
  1504. NB13=0
  1505. NB14=0
  1506. NB15=0
  1507. NB16=0
  1508. MLONG=IPT8.NUM(/1)
  1509. DO 776 MM=1,IPT8.NUM(/2)
  1510. c write (*,*) 'IPT8.NUM(MLONG,MM)', IPT8.NUM(MLONG,MM)
  1511. IF (IPT8.NUM(MLONG,MM).EQ.1) NB11=NB11+1
  1512. IF (IPT8.NUM(MLONG,MM).EQ.2) NB12=NB12+1
  1513. IF (IPT8.NUM(MLONG,MM).EQ.3) NB13=NB13+1
  1514. IF (IPT8.NUM(MLONG,MM).EQ.4) NB14=NB14+1
  1515. IF (IPT8.NUM(MLONG,MM).EQ.5) NB15=NB15+1
  1516. IF (IPT8.NUM(MLONG,MM).EQ.6) NB16=NB16+1
  1517. 776 CONTINUE
  1518. C
  1519. C 2/ Cas des relations 'aretes a 2 noeuds'
  1520. IF (NB11.NE.0) THEN
  1521. MCONF=MCONF+1
  1522. c write (*,*) 'incr Mconf coul : 1'
  1523. MCOL=0
  1524. NBNN=3
  1525. NBELEM=NB11
  1526. NBREF=0
  1527. NBSOUS=0
  1528. SEGINI IPT11
  1529. IPT11.ITYPEL=48
  1530. DO 778 MM11=1,IPT8.NUM(/2)
  1531. IF (IPT8.NUM(MLONG,MM11).NE.1) GOTO 778
  1532. MCOL=MCOL+1
  1533. DO 777 MTMP=1,3
  1534. IPT11.NUM(MTMP,MCOL)=IPT8.NUM(MTMP,MM11)
  1535. IPT11.ICOLOR(MCOL)=1
  1536. 777 CONTINUE
  1537. 778 CONTINUE
  1538. ENDIF
  1539. C
  1540. C 3/ Cas des relations 'aretes a 3 noeuds'
  1541. IF (NB12.NE.0) THEN
  1542. MCONF=MCONF+1
  1543.  
  1544. c write (*,*) 'incr Mconf coul : 2'
  1545. MCOL=0
  1546. NBNN=4
  1547. NBELEM=NB12
  1548. NBREF=0
  1549. NBSOUS=0
  1550. SEGINI IPT12
  1551. IPT12.ITYPEL=48
  1552. DO 668 MM11=1,IPT8.NUM(/2)
  1553. IF (IPT8.NUM(MLONG,MM11).NE.2) GOTO 668
  1554. MCOL=MCOL+1
  1555. DO 667 MTMP=1,4
  1556. IPT12.NUM(MTMP,MCOL)=IPT8.NUM(MTMP,MM11)
  1557. IPT12.ICOLOR(MCOL)=2
  1558. 667 CONTINUE
  1559. 668 CONTINUE
  1560. ENDIF
  1561. C
  1562. C 4/ Cas des relations 'face QUA4'
  1563. IF (NB13.NE.0) THEN
  1564. MCONF=MCONF+1
  1565. c write (*,*) 'incr Mconf coul : 3'
  1566. MCOL=0
  1567. NBNN=5
  1568. NBELEM=NB13
  1569. NBREF=0
  1570. NBSOUS=0
  1571. SEGINI IPT13
  1572. IPT13.ITYPEL=48
  1573. DO 658 MM11=1,IPT8.NUM(/2)
  1574. IF (IPT8.NUM(MLONG,MM11).NE.3) GOTO 658
  1575. MCOL=MCOL+1
  1576. DO 657 MTMP=1,5
  1577. IPT13.NUM(MTMP,MCOL)=IPT8.NUM(MTMP,MM11)
  1578. IPT13.ICOLOR(MCOL)=3
  1579. 657 CONTINUE
  1580. 658 CONTINUE
  1581. ENDIF
  1582. C
  1583. C 5/ Cas des relations 'face TRI6'
  1584. IF (NB14.NE.0) THEN
  1585. MCONF=MCONF+1
  1586. c write (*,*) 'incr Mconf coul : 4'
  1587. MCOL=0
  1588. NBNN=7
  1589. NBELEM=NB14
  1590. NBREF=0
  1591. NBSOUS=0
  1592. SEGINI IPT14
  1593. IPT14.ITYPEL=48
  1594. DO 648 MM11=1,IPT8.NUM(/2)
  1595. IF (IPT8.NUM(MLONG,MM11).NE.4) GOTO 648
  1596. MCOL=MCOL+1
  1597. DO 647 MTMP=1,7
  1598. IPT14.NUM(MTMP,MCOL)=IPT8.NUM(MTMP,MM11)
  1599. IPT14.ICOLOR(MCOL)=4
  1600. 647 CONTINUE
  1601. 648 CONTINUE
  1602. ENDIF
  1603. C
  1604. C 6/ Cas des relations 'face QUA8-1'
  1605. IF (NB15.NE.0) THEN
  1606. MCONF=MCONF+1
  1607. c write (*,*) 'incr Mconf coul : 5'
  1608. MCOL=0
  1609. NBNN=9
  1610. NBELEM=NB15
  1611. NBREF=0
  1612. NBSOUS=0
  1613. SEGINI IPT15
  1614. IPT15.ITYPEL=48
  1615. DO 638 MM11=1,IPT8.NUM(/2)
  1616. IF (IPT8.NUM(MLONG,MM11).NE.5) GOTO 638
  1617. MCOL=MCOL+1
  1618. DO 637 MTMP=1,9
  1619. IPT15.NUM(MTMP,MCOL)=IPT8.NUM(MTMP,MM11)
  1620. IPT15.ICOLOR(MCOL)=5
  1621. 637 CONTINUE
  1622. 638 CONTINUE
  1623. ENDIF
  1624. C
  1625. C 7/ Cas des relations 'face QUA8-2'
  1626. IF (NB16.NE.0) THEN
  1627. MCONF=MCONF+1
  1628. c write (*,*) 'incr Mconf coul : 6'
  1629. MCOL=0
  1630. NBNN=9
  1631. NBELEM=NB16
  1632. NBREF=0
  1633. NBSOUS=0
  1634. SEGINI IPT16
  1635. IPT16.ITYPEL=48
  1636. DO 678 MM11=1,IPT8.NUM(/2)
  1637. IF (IPT8.NUM(MLONG,MM11).NE.6) GOTO 678
  1638. MCOL=MCOL+1
  1639. DO 677 MTMP=1,9
  1640. IPT16.NUM(MTMP,MCOL)=IPT8.NUM(MTMP,MM11)
  1641. IPT16.ICOLOR(MCOL)=6
  1642. 677 CONTINUE
  1643. 678 CONTINUE
  1644. ENDIF
  1645. C
  1646. ENDIF
  1647. C
  1648. C--------------------------------------
  1649. C D) Creation du maillage de relations !
  1650. C--------------------------------------
  1651. C 1/ Initialisation de MAIREL
  1652. C MAIREL contient tous les maillages de relations, certains étant vides
  1653. C car il n'y a aucune relation du type considéré dans le nveau maillage
  1654. NBMAIL=6
  1655. SEGINI MAIREL
  1656. MAIREL(1)=IPT11
  1657. MAIREL(2)=IPT12
  1658. MAIREL(3)=IPT13
  1659. MAIREL(4)=IPT14
  1660. MAIREL(5)=IPT15
  1661. MAIREL(6)=IPT16
  1662. C
  1663. C 2/ Initialisation de LTYPE2
  1664. C LTYPE2 va contenir tous les sous-objets du nouveau maillage,
  1665. C relations et non relations
  1666. NBLIS2=MCONF+MAX(1,MELEME.LISOUS(/1))
  1667. SEGINI LTYPE2
  1668. C
  1669. C 3/ Renseignement des lisous non 48 dans LTYPE2
  1670. C S'il y a N sous-objets de non-relations dans le nouveau maillage,
  1671. C les N-ièmes premières colonnes de LTYPE2 contiennent ces sous-objets
  1672. DO 685 IZZ=1,MAX(1,MELEME.LISOUS(/1))
  1673. IF (MELEME.LISOUS(/1).EQ.0) THEN
  1674. IPT2=MELEME
  1675. ELSE
  1676. IPT2=MELEME.LISOUS(IZZ)
  1677. SEGACT IPT2
  1678. ENDIF
  1679. LTYPE2(IZZ,1)=IZZ
  1680. LTYPE2(IZZ,2)=IPT2.ITYPEL
  1681. LTYPE2(IZZ,3)=IPT2
  1682. LTYPE2(IZZ,4)=0
  1683. 685 CONTINUE
  1684. C
  1685. C 4/ Renseignement des lisous 48 dans LTYPE2
  1686. C Les colonnes suivantes contiennent les sous-objets de relations,
  1687. C chaque sous-objet correspondant a un type de relations de conformité
  1688. IZZ=NBLIS2-MCONF
  1689. DO 686 IXX=1,6
  1690. IF (MAIREL(IXX).EQ.0) GOTO 686
  1691. IPT2=MAIREL(IXX)
  1692. IZZ=IZZ+1
  1693. LTYPE2(IZZ,1)=IZZ
  1694. LTYPE2(IZZ,2)=IPT2.ITYPEL
  1695. LTYPE2(IZZ,3)=IPT2
  1696. LTYPE2(IZZ,4)=IPT2.ICOLOR(1)
  1697. 686 CONTINUE
  1698. C--------------------------------------
  1699. C E) Renseignement de LTYPE2 interne a
  1700. c la boucle 33 dans LTYPE3 externe a
  1701. c la boucle 33
  1702. C--------------------------------------C
  1703. C
  1704. c write (*,*) 'LTYPE2', LTYPE2, 'LTYPE3(INB,1)', LTYPE3(INB,1)
  1705. c write (*,*) 'NBLIS2-MCONF+1', NBLIS2-MCONF+1
  1706. c write (*,*) 'stokage du tableau Ltyp2 pour l itération', INB
  1707. LTYPE3(INB,1)=LTYPE2
  1708. IF (MCONF.EQ.0) THEN
  1709. LTYPE3(INB,2)=0
  1710. ELSE
  1711. LTYPE3(INB,2)=NBLIS2-MCONF+1
  1712. ENDIF
  1713. c write (*,*) 'LTYPE3(INB,1)', LTYPE3(INB,1)
  1714. c write (*,*) 'LTYPE3(INB,2)', LTYPE3(INB,2)
  1715. c fin de boucle sur les sous domaines d'ipt5
  1716. 33 CONTINUE
  1717. c write (*,*) 'SORTIE DE LA BOUCLE 33'
  1718. c SEGSUP IPT2
  1719.  
  1720.  
  1721. C Mise a jour des maillages en entrée sortie de la boucle 9
  1722. C IPT5, IPT6, IPT7
  1723. C=======================================================================
  1724. 679 CONTINUE
  1725. C--------------------------------------C
  1726. C A) initialisation
  1727. C--------------------------------------C
  1728.  
  1729. c 1/ initialisaton de IPT7 maillage stockant tout les sous maillage
  1730. c à la fin d'une itération de la boucle 9
  1731. NBS7=0
  1732. DO 331 I3=1,LTYPE3(/1)
  1733. c write (*,*) 'I3' , I3
  1734. c write (*,*) 'LTYPE3(I3, 1)' , LTYPE3(I3, 1)
  1735. LTYPE2=LTYPE3(I3, 1)
  1736. SEGACT LTYPE2
  1737. c write (*,*) 'LTYPE2(/1)', LTYPE2(/1)
  1738. NBS7=NBS7+LTYPE2(/1)
  1739. 331 CONTINUE
  1740. IF (NBS7.EQ.1) NBS7=0
  1741. c write (*,*) 'NBS7', NBS7
  1742. NBNN=0
  1743. NBELEM=0
  1744. NBREF=0
  1745. NBSOUS=NBS7
  1746. SEGINI IPT7
  1747.  
  1748. c 2/ initialisaton de IPT6 maillage stockant tout les sous maillage
  1749. c de relation à la fin d'une itération de la boucle 9
  1750. NBS6=0
  1751.  
  1752. DO 332 I3=1,LTYPE3(/1)
  1753. c write (*,*) 'LTYPE3(I3,2)', LTYPE3(I3,2)
  1754. LTYPE2=LTYPE3(I3, 1)
  1755. SEGACT LTYPE2
  1756. IF (LTYPE3(I3,2).NE.0) NBS6=NBS6+(LTYPE2(/1)-LTYPE3(I3,2)+1)
  1757. 332 CONTINUE
  1758. IF (NBS6.EQ.1) NBS6=0
  1759. c write (*,*) 'NBS6', NBS6
  1760. NBNN=0
  1761. NBELEM=0
  1762. NBREF=0
  1763. NBSOUS=NBS6
  1764. SEGINI IPT6
  1765. c 3/ initialisaton de IPT5 maillage stockant tout les sous maillage
  1766. c non relation à la fin d'une itération de la boucle 9
  1767. NBS5=0
  1768. DO 333 I3=1,LTYPE3(/1)
  1769. LTYPE2=LTYPE3(I3, 1)
  1770. SEGACT LTYPE2
  1771. IF (LTYPE3(I3,2).EQ.0) THEN
  1772. NBS5=NBS5+LTYPE2(/1)
  1773. ELSE
  1774. NBS5=NBS5+(LTYPE3(I3,2)-1)
  1775. ENDIF
  1776. 333 CONTINUE
  1777. IF (NBS5.EQ.1) NBS5=0
  1778. c write (*,*) 'NBS5', NBS5
  1779. NBNN=0
  1780. NBELEM=0
  1781. NBREF=0
  1782. NBSOUS=NBS5
  1783. SEGINI IPT5
  1784. C--------------------------------------C
  1785. C B) assemblage des nouveaux maillage à
  1786. c partir des sous maillages stockés
  1787. c dans LTYP3
  1788. C--------------------------------------C
  1789. IND5=0
  1790. IND6=0
  1791. IND7=0
  1792. c 1/assemblage des nouveaux maillage non relation à partir des sous
  1793. c maillages stockés dans LTYP3
  1794.  
  1795. DO 35 I3=1,LTYPE3(/1)
  1796. c write (*,*) 'I3' , I3
  1797. c extraction d'un LTYPE2 de LTYPE3
  1798. LTYPE2=LTYPE3(I3, 1)
  1799. SEGACT LTYPE2
  1800. c NBNN=0
  1801. c NBELEM=0
  1802. c NBREF=0
  1803. c NBSOUS=0
  1804. DO 613 I2=1,MAX(1,(LTYPE3(I3,2)-1))
  1805. c write (*,*) 'I2' , I2
  1806. c extraction d'un MAILLAGE non relation de LTYPE2
  1807. IPT2=LTYPE2(I2,3)
  1808. c write (*,*) 'EXTRACTION D UN IPT2'
  1809. SEGACT IPT2
  1810. c write (*,*) 'IPT2: Type' , IPT2.ITYPEL,'nbelem', IPT2.NUM(/2)
  1811. IND5=IND5+1
  1812. IF (NBS5.EQ.0) THEN
  1813. IPT5=IPT2
  1814. c write (*,*) 'IPT5' , IPT5
  1815. ELSE
  1816. IPT5.LISOUS(IND5)=IPT2
  1817. c write (*,*) 'IPT5.LISOUS(',IND5,')' , IPT5.LISOUS(IND5)
  1818. ENDIF
  1819. IND7=IND7+1
  1820. IF (NBS7.EQ.0) THEN
  1821. IPT7=IPT2
  1822. c write (*,*) 'IPT7' , IPT7
  1823. ELSE
  1824. IPT7.LISOUS(IND7)=IPT2
  1825. c write (*,*) 'IPT7.LISOUS(',IND7,')' , IPT7.LISOUS(IND7)
  1826. ENDIF
  1827. 613 CONTINUE
  1828. 35 CONTINUE
  1829.  
  1830. c 2/assemblage des nouveaux maillage de relation à partir des sous
  1831. c maillages stockés dans LTYPE3
  1832.  
  1833. DO 36 I3=1,LTYPE3(/1)
  1834. c extraction d'un LTYPE2 de LTYPE3
  1835. c write (*,*) 'I3', I3
  1836. c write (*,*) 'LTYPE3(I3,2)', LTYPE3(I3,2)
  1837. IF (LTYPE3(I3,2).EQ.0) GOTO 36
  1838. SEGACT LTYPE2
  1839. LTYPE2=LTYPE3(I3, 1)
  1840. NBNN=0
  1841. NBELEM=0
  1842. NBREF=0
  1843. NBSOUS=0
  1844. DO 614 I2=LTYPE3(I3,2), LTYPE2(/1)
  1845. c extraction d'un MAILLAGE de relation de LTYPE2
  1846. c write (*,*) 'I2', I2
  1847. SEGINI IPT2
  1848. IPT2=LTYPE2(I2,3)
  1849. SEGACT IPT2
  1850. c write (*,*) 'IPT2: Type', IPT2.ITYPEL,'nbelem', IPT2.NUM(/2)
  1851. c write (*,*) 'coul IPT2' , IPT2.ICOLOR(1)
  1852. IND6=IND6+1
  1853. IF (NBS6.EQ.0) THEN
  1854. IPT6=IPT2
  1855. c write (*,*) 'IPT6' , IPT6
  1856. ELSE
  1857. IPT6.LISOUS(IND6)=IPT2
  1858. c write (*,*) 'IPT6.LISOUS(',IND6,')' , IPT6.LISOUS(IND6)
  1859. ENDIF
  1860. IND7=IND7+1
  1861. IPT7.LISOUS(IND7)=IPT2
  1862. c write (*,*) 'IPT7.LISOUS(',IND7, ')', IPT7.LISOUS(I2)
  1863. 614 CONTINUE
  1864. 36 CONTINUE
  1865. c stockage de LTYPE3 dans LTYPE4
  1866.  
  1867. SEGINI LTYPE4
  1868. LTYPE4 = LTYPE3
  1869. c WRITE (*,*) 'LTYPE4', LTYPE4, 'LTYPE3', LTYPE4
  1870. c mise a jour du nombre global de relations
  1871. c WRITE (*,*) 'Mconf', MCONF
  1872. LCONF3=MCONF
  1873.  
  1874.  
  1875. C
  1876. C=======================================================================
  1877. C Test sur la progression du raffinement
  1878. C=======================================================================
  1879. C
  1880. C
  1881. SEGDES IPT5,MCHPO1,MSOUPO,MPOVA1,IPT3,KMILIE
  1882. c SEGSUP ICPR,KARETE,KARPOS,KARET2,ICHPO,IDENSI
  1883. c SEGSUP ICPR,KARETE,KARPOS,KARET2,IDENSI
  1884. SEGSUP KARETE,KARPOS,KARET2
  1885. SEGSUP MCHAM1,MELVA1,MCHAM2,MELVA2,IMODE1
  1886. SEGSUP MMODE1,MCHEL1,MCHEL2
  1887.  
  1888. IF ((IDIM.EQ.3).AND.(NBOUCL.NE.1)) SEGSUP JPLANS,JPLAN2,JPLAN3,
  1889. +JPLCOM,JNOEFA,IWORK2,JFARAF, IPT9
  1890.  
  1891.  
  1892. C S'il y a eu a raffiner des elements lors de l'iteration courante,
  1893. C on passe a l'iteration suivante
  1894. C
  1895. c !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  1896. c forcage : nb max itération de raff
  1897. c IF (NBOUCL.EQ.3) THEN
  1898. c KARAF2=0
  1899. c KARAF=0
  1900. c ENDIF
  1901. c !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  1902. IF (KARAF2.NE.0) GOTO 9
  1903.  
  1904.  
  1905.  
  1906. c WRITE (*,*) 'sortie de la boucle 9 '
  1907. C
  1908. C ### SORTIE TEMPORAIRE DU PROGRAMME ###
  1909. C Si on est a la premiere iteration de RAFF et qu'il n'y a pas
  1910. C d'elements a raffiner, on retourne en sortie le maillage de depart,
  1911. C donné en entree
  1912. IF ((KARAF2.EQ.0).AND.(NBOUCL.EQ.1)) THEN
  1913. WRITE(*,*) 'il n est pas necessaire de raffiner le maillage'
  1914. MELEME=IPT1
  1915. GOTO 998
  1916.  
  1917. ENDIF
  1918.  
  1919.  
  1920. C
  1921. C
  1922. C=======================================================================
  1923. C Mises a jour diverses
  1924. C=======================================================================
  1925. C ** On affecte une densite a chaque noeud que l'on a cree
  1926. c WRITE (*,*) 'NPTINI ', NPTINI
  1927. c WRITE (*,*) 'NPTFIN', XCOOR(/1)/(IDIM+1)
  1928. c WRITE (*,*) 'DENSI', DENSI
  1929. c WRITE (*,*) 'DENSI(1)', DENSI(1)
  1930.  
  1931. c DO 37 I=NPTINI+1,XCOOR(/1)/(IDIM+1)
  1932. c WRITE (*,*) i
  1933. c WRITE (*,*) ICPR(I,1)
  1934. c WRITE(*,*) DENSI(ICPR(I,1))
  1935. c XCOOR(I*(IDIM+1))=DENSI(ICPR(I,1))
  1936. c 37 CONTINUE
  1937.  
  1938. C
  1939. C=======================================================================
  1940. C Fin du programme
  1941. C=======================================================================
  1942.  
  1943. MELEME=IPT7
  1944.  
  1945. SEGSUP XDEN
  1946. SEGDES IPT1
  1947. 998 CONTINUE
  1948.  
  1949. c WRITE (*,*) 'appel a ecrobj', MELEME
  1950. CALL ECROBJ('MAILLAGE',MELEME)
  1951.  
  1952. c WRITE (*,*) 'presque fin'
  1953. RETURN
  1954. c WRITE (*,*) 'fin'
  1955. END
  1956.  
  1957.  
  1958.  
  1959.  
  1960.  
  1961.  
  1962.  
  1963.  
  1964.  
  1965.  
  1966.  
  1967.  
  1968.  
  1969.  

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