Télécharger testma.eso

Retour à la liste

Numérotation des lignes :

testma
  1. C TESTMA SOURCE PV 20/03/30 21:25:08 10567
  2. SUBROUTINE TESTMA(IPCHE,IMEL,COND,CONST,IRET,IMODI)
  3. *______________________________________________________________________
  4. *
  5. * test sur les maillages
  6. *
  7. * entrees :
  8. * ---------
  9. * ipche chamelem dont il faut vérifier le maillage (type mchaml)
  10. * imel maillage du modèle: sert de référence (type meleme)
  11. * Maillage actif en Entree (?) et en Sortie (oui)
  12. * cond vrai si il faut tenir compte de const
  13. * const nom de constituant
  14. *
  15. *
  16. * sortie :
  17. * --------
  18. * iret chamelem réordonné si nécessaire
  19. * = 0 si pb
  20. * imodi =1 il a fallu decouper le maillage imel
  21. *
  22. *______________________________________________________________________
  23. *
  24. * declarations
  25. *
  26. IMPLICIT REAL*8(A-H,O-Z)
  27. IMPLICIT INTEGER(I-N)
  28. *
  29.  
  30. -INC PPARAM
  31. -INC CCOPTIO
  32. *
  33. -INC SMCHAML
  34. -INC SMELEME
  35. -INC SMCOORD
  36. -INC SMINTE
  37. *
  38. SEGMENT ICPR(NPT)
  39. SEGMENT ICPEL(ICPELD)
  40. SEGMENT ITRAV(NNNT)
  41. SEGMENT NIMEL(NNNT)
  42. SEGMENT IPERM(NNNT)
  43. SEGMENT NZONE(NSOUS)
  44. SEGMENT ICOM(NBEL)
  45. *
  46. SEGMENT ITRA
  47. INTEGER ICC (IA),IRE(IA,IMA)
  48. ENDSEGMENT
  49. *
  50. SEGMENT WTRAV
  51. REAL*8 QSIREF(NBPGAU),QSICHAM(NBPGAU),ETAREF(NBPGAU)
  52. REAL*8 ETACHAM(NBPGAU),DZEREF(NBPGAU),DZECHAM(NBPGAU)
  53. REAL*8 XECHAM(3,NP1),XEREF(3,NP1)
  54. INTEGER TABOK(NBPGAU),TAB(NBPGAU)
  55. ENDSEGMENT
  56. LOGICAL INIT,COND
  57. CHARACTER*(*) CONST
  58. *
  59. * executable
  60. *
  61. NPT = nbpts
  62. ICPR = 0
  63. ICPEL = 0
  64. *
  65. IRET = 0
  66. IMODI = 0
  67. *
  68. IPT1 = IMEL
  69. SEGACT IPT1
  70. NBSOU1=IPT1.LISOUS(/1)
  71. * La routine ne travaille qu'avec un maillage elementaire
  72. IF ( NBSOU1 .GT. 1) THEN
  73. C SEGDES,IPT1
  74. call erreur(25)
  75. RETURN
  76. ENDIF
  77. MELEME = IPT1
  78. *
  79. * champ servant à stocker provisoirement le resultat
  80. *
  81. MCHELM = IPCHE
  82. SEGACT MCHELM
  83. SEGINI,MCHEL1=MCHELM
  84. *
  85. * segment servant a stocker l'ordre de stockage dans mchel1
  86. NS01=0
  87. NSOUS=ICHAML(/1)
  88. SEGINI NZONE
  89. *
  90. * creation de wtrav : on l'ajustera ensuite si necessaire
  91. * mais on sort la creation destruction de la boucle
  92. * la definition doit etre coherente avec celle dans tabgau
  93. nbpgau=0
  94. np1=0
  95. wtrav=0
  96. segini wtrav
  97. *
  98. * boucle conditionnelle sur les sous zones du maillage imel
  99. *
  100. IBOU1 = 0
  101. *
  102. SEGINI ICPR
  103. 1 CONTINUE
  104. IBOU1 = IBOU1 + 1
  105. *
  106. * Initialisations de ICPR, ICPEL la premiere fois puis
  107. * remises a zero pour les fois suivantes faite apres 2
  108. *
  109. IF (NBSOU1.NE.0) THEN
  110. MELEME = IPT1.LISOUS(IBOU1)
  111. ENDIF
  112. SEGACT MELEME
  113. *
  114. * nnnt :nb d'elements de imel pour la sous zone consideree
  115. *
  116. NNNT=NUM(/2)
  117. *
  118. IA=0
  119. IMA=0
  120. *
  121. * fabrication de icpel qui dira combien de fois un point apparait
  122. * de icpr qui donne une numerotation locale
  123. *
  124. DO J=1,NNNT
  125. DO 3 K=1,NUM(/1)
  126. ID=NUM(K,J)
  127. IF(ICPR(ID).NE.0) GO TO 3
  128. * numerotation locale
  129. IA=IA+1
  130. ICPR(ID)=IA
  131. 3 CONTINUE
  132. enddo
  133. ICPELD=IA
  134. IF(ICPEL.EQ.0) THEN
  135. SEGINI ICPEL
  136. ELSE
  137. IF(ICPELD.GT.ICPEL(/1)) SEGADJ ICPEL
  138. call ooozmr(icpel(1),icpeld)
  139. ENDIF
  140.  
  141. DO J=1,NNNT
  142. DO K=1,NUM(/1)
  143. ID=NUM(K,J)
  144. ICPEL(ICPR(ID))=ICPEL(ICPR(ID))+1
  145. * on retient le max des occurences des points
  146. IMA=MAX(ICPEL(ICPR(ID)),IMA)
  147. enddo
  148. enddo
  149. *
  150. * fabrication de itra :
  151. * icc donne le nombre d'elements a partir du numero local
  152. * ire donne les numeros d'elements de imel
  153. * a partir du numero local et de l'occurence
  154. *
  155. SEGINI ITRA
  156. DO J=1,NNNT
  157. DO K=1,NUM(/1)
  158. ID=ICPR(NUM(K,J))
  159. ICC(ID)=ICC(ID)+1
  160. IA=ICC(ID)
  161. IRE(ID,IA)=J
  162. enddo
  163. enddo
  164. *
  165. * il faut maintenant regarder si dans un imache il existe
  166. * les elements de la sous zone du meleme imel
  167. * boucle sur les sous zone du chamelem
  168. *
  169. IFLAG=0
  170. * initialisons ipt2
  171. IPT2=-1
  172. DO 10 I=1,ICHAML(/1)
  173. IF (COND) THEN
  174. IF (CONCHE(I) .NE. CONST) GOTO 10
  175. ENDIF
  176. MCHAML=ICHAML(I)
  177. IPT2 =IMACHE(I)
  178. IF (IPT2 .EQ. MELEME) THEN
  179. IFLAG=1
  180. NS01=NS01+1
  181. IF(NS01 .GT. ICHAML(/1))THEN
  182. N1=NS01
  183. N3=INFCHE(/2)
  184. L1=TITCHE(/1)
  185. SEGADJ,MCHEL1
  186. DO N33=1,N3
  187. MCHEL1.INFCHE(NS01,N33)=MCHEL1.INFCHE(I,N33)
  188. ENDDO
  189. ENDIF
  190. NSOUS=NSOUS+1
  191. SEGADJ,NZONE
  192. NZONE(NS01)=I
  193. MCHEL1.ICHAML(NS01)=MCHEL1.ICHAML(I)
  194. MCHEL1.IMACHE(NS01)=MELEME
  195. MCHEL1.CONCHE(NS01)=MCHEL1.CONCHE(I)
  196. GOTO 10
  197. ENDIF
  198. SEGACT,MELEME,IPT2
  199. IF(ITYPEL.NE.IPT2.ITYPEL) GO TO 11
  200. *
  201. NBEL = NUM(/2)
  202. SEGINI ICOM
  203. *
  204. INIT=.FALSE.
  205. *
  206. ICO=0
  207. NP1=IPT2.NUM(/1)
  208. DO 12 K=1,IPT2.NUM(/2)
  209. DO 13 L=1,NP1
  210. ID=IPT2.NUM(L,K)
  211. IDD = ICPR(ID)
  212. IF(IDD.EQ.0) GO TO 12
  213. 13 CONTINUE
  214. *
  215. * ok l'element k possede ses noeuds dans imel
  216. *
  217. IF(ITYPEL.EQ.1) THEN
  218. ID=IPT2.NUM(1,K)
  219. IDD=ICPR(ID)
  220. IRE1=IRE(IDD,1)
  221. GOTO 20
  222. ENDIF
  223. *
  224. * ces noeuds correspondent-t-ils a un meme element ire1
  225. *
  226. ID = IPT2.NUM(1,K)
  227. IDD=ICPR(ID)
  228. IDE=ICC(IDD)
  229. c boucle sur les occurences de idd dans les elements
  230. DO 14 L=1,IDE
  231. IRE1=IRE(IDD,L)
  232. DO 15 M=2,NP1
  233. IDD2=ICPR(IPT2.NUM(M,K))
  234. IF(IDD2.EQ.0) GO TO 12
  235. IDE2=ICC(IDD2)
  236. DO 16 N=1,IDE2
  237. IF(IRE(IDD2,N).EQ.IRE1) GO TO 15
  238. 16 CONTINUE
  239. GO TO 14
  240. 15 CONTINUE
  241. GOTO 20
  242. 14 CONTINUE
  243. *
  244. * l'element n'appartient pas a meleme
  245. GOTO 12
  246. *
  247. * itrav donne le numero de l'element dans le chamelem
  248. * nimel donne le numero de l'element dans la sous zone
  249. * de imel
  250. cbp IPERM = 0 si pas de permutation des noeuds
  251. *
  252. 20 CONTINUE
  253. c
  254. IF(.NOT.INIT) THEN
  255. INIT=.TRUE.
  256. SEGINI ITRAV
  257. SEGINI NIMEL
  258. SEGINI IPERM
  259. ENDIF
  260. ICO=ICO+1
  261. IF (ICO .GT. NNNT) THEN
  262. CALL ERREUR(472)
  263. RETURN
  264. ENDIF
  265. ITRAV(ICO)=K
  266. NIMEL(ICO)=IRE1
  267. ICOM(IRE1)=1
  268. IPERM(ICO)=0
  269. do iii=1,NP1
  270. if(IPT2.NUM(iii,K).ne.NUM(iii,IRE1)) then
  271. IPERM(ICO)=iii
  272. goto 12
  273. endif
  274. enddo
  275. 12 CONTINUE
  276. *
  277. IF(INIT) THEN
  278. *
  279. * le nb d'elements est-il egal pour la sous zone imel et
  280. * la sous zone mchaml
  281. *
  282. IF (ICO.EQ.NNNT) THEN
  283. IFLAG=1
  284. NS01=NS01+1
  285. IF(NS01 .GT. ICHAML(/1))THEN
  286. N1=NS01
  287. N3=INFCHE(/2)
  288. L1=TITCHE(/1)
  289. SEGADJ,MCHEL1
  290. DO N33=1,N3
  291. MCHEL1.INFCHE(NS01,N33)=MCHEL1.INFCHE(I,N33)
  292. ENDDO
  293. ENDIF
  294. NSOUS=NSOUS+1
  295. SEGADJ,NZONE
  296. NZONE(NS01)=I
  297. MINTE=INFCHE(I,4)
  298. IPMINT=MINTE
  299. IF (IPMINT.NE.0) SEGACT,MINTE
  300. MCHAML=ICHAML(I)
  301. SEGINI,MCHAM1=MCHAML
  302. MCHEL1.ICHAML(NS01)=MCHAM1
  303. MCHEL1.IMACHE(NS01)=MELEME
  304. MCHEL1.CONCHE(NS01)=MCHEL1.CONCHE(I)
  305. *
  306. DO 24 ICOMP=1,MCHAM1.IELVAL(/1)
  307. MELVAL=MCHAM1.IELVAL(ICOMP)
  308. SEGACT MELVAL
  309. IMEM=0
  310. IF (MCHAM1.TYPCHE(ICOMP).EQ.'REAL*8') THEN
  311. N1PTEL=VELCHE(/1)
  312. N1EL=ICO
  313. IF (VELCHE(/2).EQ.1) THEN
  314. N1EL=1
  315. IMEM=1
  316. ENDIF
  317. NEL=N1EL
  318. N2PTEL=0
  319. N2EL =0
  320. ELSE
  321. N2PTEL=IELCHE(/1)
  322. N2EL=ICO
  323. IF (IELCHE(/2).EQ.1) THEN
  324. N2EL=1
  325. IMEM=1
  326. ENDIF
  327. NEL=N2EL
  328. N1PTEL=0
  329. N1EL =0
  330. ENDIF
  331. *
  332. SEGINI MELVA1
  333. *
  334. DO 21 J=1,NEL
  335. K =ITRAV(J)
  336. IRE1=NIMEL(J)
  337. KTEST=K
  338. ITEST=IRE1
  339. IF (NEL.EQ.1) THEN
  340. IRE1=1
  341. IF (IMEM.EQ.1) THEN
  342. K=1
  343. ENDIF
  344. ENDIF
  345. *
  346. * dans le cas des champs qui ne sont pas définis aux noeuds
  347. * vérification de la bonne orientation des maillages
  348. * sinon création d'un tableau de correspondance
  349. *
  350. NBPGAU=MAX(N1PTEL,N2PTEL)
  351. cbp -cas ou l'on doit identifier la position des pts d'integration
  352. cbp (ne marche pas avec des coques epaisses MFR=5,
  353. cbp ni avec DKT dont nombre de points >1 dans l'épaisseur)
  354. IF((IPMINT.NE.0).and.(IPERM(J).ne.0)) THEN
  355.  
  356. if (qsiref(/1).ne.nbpgau .or.
  357. & xeref(/2).ne.np1) then
  358. segadj,wtrav
  359. endif
  360. CALL TABGAU(IPMINT,ITEST,KTEST,NP1,MELEME,
  361. + IPT2,NBPGAU,IRET1,wtrav)
  362. IF (IRET1.EQ.0) THEN
  363. CALL ERREUR(757)
  364. SEGSUP MCHEL1,ITRA
  365. C if (imel.ne.ipt2) SEGDES IPT2
  366. C if (imel.ne.meleme) segdes meleme
  367. GOTO 900
  368. ENDIF
  369. IF (N1PTEL.EQ.0) THEN
  370. DO 22 IGAU=1,N2PTEL
  371. IGAU1=TAB(IGAU)
  372. MELVA1.IELCHE(IGAU,IRE1)=IELCHE(IGAU1,K)
  373. 22 CONTINUE
  374. ELSE
  375. DO 23 IGAU=1,N1PTEL
  376. IGAU1=TAB(IGAU)
  377. MELVA1.VELCHE(IGAU,IRE1)=VELCHE(IGAU1,K)
  378. 23 CONTINUE
  379. ENDIF
  380.  
  381. cbp -cas ou on recopie directement sans chercher la permutation
  382. ELSE
  383. IF (N1PTEL.EQ.0) THEN
  384. DO 25 IGAU=1,N2PTEL
  385. MELVA1.IELCHE(IGAU,IRE1)=IELCHE(IGAU,K)
  386. 25 CONTINUE
  387. ELSE
  388. DO 26 IGAU=1,N1PTEL
  389. MELVA1.VELCHE(IGAU,IRE1)=VELCHE(IGAU,K)
  390. 26 CONTINUE
  391. ENDIF
  392.  
  393. ENDIF
  394. *
  395. 21 CONTINUE
  396. MCHAM1.IELVAL(ICOMP)=MELVA1
  397. C SEGDES MELVA1,MELVAL
  398. 24 CONTINUE
  399. C SEGDES MCHAM1,MCHAML
  400. SEGSUP ITRAV,NIMEL,IPERM
  401. ELSE
  402. *
  403. * le maillage est a cheval sur plusieurs sous zones
  404. * du mchaml
  405. * il va falloir scinder le maillage en deux
  406. * et recommencer le processus
  407. *
  408. IMODI = 1
  409. * creation des deux maillages dont l'union est le maillage meleme
  410. NBNN = NUM(/1)
  411. NBELEM = ICO
  412. NBREF = 0
  413. NBSOUS = 0
  414. SEGINI IPT3
  415. NBELEM = NUM(/2) - ICO
  416. IF (NBELEM.LE.0) THEN
  417. CALL ERREUR(26)
  418. GOTO 900
  419. ENDIF
  420. SEGINI IPT4
  421. IPT3.ITYPEL = ITYPEL
  422. IPT4.ITYPEL = ITYPEL
  423. I3 = 0
  424. I4 = 0
  425. DO 45 II=1,NUM(/2)
  426. IF (ICOM(II) .EQ. 1) THEN
  427. * l'element est commun aux deux maillages
  428. I3=I3+1
  429. IPT3.ICOLOR(I3)=ICOLOR(II)
  430. DO 43 JJ=1,NUM(/1)
  431. IPT3.NUM(JJ,I3)=NUM(JJ,II)
  432. 43 CONTINUE
  433. ELSE
  434. * l'element appartient seulement a meleme
  435. I4=I4+1
  436. IPT4.ICOLOR(I4)=ICOLOR(II)
  437. DO 44 JJ=1,NUM(/1)
  438. IPT4.NUM(JJ,I4)=NUM(JJ,II)
  439. 44 CONTINUE
  440. ENDIF
  441. 45 CONTINUE
  442. *
  443. * modification de ipt1
  444. *
  445. IF (IPT1.LISOUS(/1) .LE. 1) THEN
  446. NBSOUS=2
  447. NBREF =0
  448. NBELEM=0
  449. NBNN =0
  450. SEGINI IPT5
  451. IPT5.LISOUS(1)=IPT3
  452. IPT5.LISOUS(2)=IPT4
  453. C if (imel.ne.ipt1) SEGDES IPT1
  454. IPT1=IPT5
  455. ELSE
  456. NBSOUS=IPT1.LISOUS(/1)+1
  457. NBREF =0
  458. NBELEM=0
  459. NBNN =0
  460. SEGADJ,IPT1
  461. IPT1.LISOUS(IBOU1) =IPT3
  462. IPT1.LISOUS(NBSOUS)=IPT4
  463. C if (imel.ne.meleme) SEGDES MELEME
  464. ENDIF
  465. IBOU1 =IBOU1-1
  466. NBSOU1=NBSOUS
  467. SEGSUP ITRAV,NIMEL,IPERM
  468. GOTO 2
  469. ENDIF
  470. ENDIF
  471. SEGSUP ICOM
  472. 11 CONTINUE
  473. C if (imel.ne.ipt2) SEGDES IPT2
  474. 10 CONTINUE
  475. *
  476. * on n'a pas trouve de sous zone dans le mchaml qui
  477. * correspondent a une sous zone de imel
  478. *
  479. IF (IFLAG.EQ.0) THEN
  480. C if (imel.ne.meleme) segdes meleme
  481. SEGSUP MCHEL1
  482. IF(IPT2.LE.0) THEN
  483. SEGSUP,ITRA
  484. ELSE
  485. C if (imel.ne.ipt2) SEGDES IPT2
  486. MOTERR(1:8)='MAILLAGE'
  487. MOTERR(9:16)='MCHAML'
  488. CALL ERREUR(135)
  489. ENDIF
  490. GOTO 900
  491. ELSE
  492. SEGSUP,ITRA
  493. C if (imel.ne.meleme) segdes meleme
  494. C if (imel.ne.ipt2) SEGDES IPT2
  495. ENDIF
  496.  
  497. 2 continue
  498. if (ibou1.lt.nbsou1) then
  499. * remise a zero de icpr
  500. DO J=1,NNNT
  501. DO K=1,NUM(/1)
  502. ID=NUM(K,J)
  503. ICPR(ID)=0
  504. enddo
  505. enddo
  506. goto 1
  507. endif
  508. *
  509. * fin de la boucle sur les sous zones du maillage
  510. *
  511. IF(NS01.NE.NSOUS) THEN
  512. N1=NS01
  513. N3=MCHEL1.INFCHE(/2)
  514. L1=MCHEL1.TITCHE(/1)
  515. SEGINI MCHEL2
  516. MCHEL2.TITCHE=TITCHE
  517. MCHEL2.IFOCHE=IFOCHE
  518. DO 31 ISOUS=1,NS01
  519. MCHEL2.ICHAML(ISOUS) = MCHEL1.ICHAML(ISOUS)
  520. MCHEL2.IMACHE(ISOUS) = MCHEL1.IMACHE(ISOUS)
  521. MCHEL2.CONCHE(ISOUS) = MCHEL1.CONCHE(ISOUS)
  522. NSS=NZONE(ISOUS)
  523. DO 33 N33=1,N3
  524. MCHEL2.INFCHE(ISOUS,N33) = MCHEL1.INFCHE(NSS,N33)
  525. 33 CONTINUE
  526. 31 CONTINUE
  527. SEGSUP MCHEL1
  528. MCHEL1=MCHEL2
  529. ENDIF
  530. C SEGDES,MCHEL1
  531. IRET=MCHEL1
  532.  
  533. 900 CONTINUE
  534. C SEGDES MCHELM
  535. SEGSUP,ICPR,ICPEL,NZONE,WTRAV
  536.  
  537. END
  538.  
  539.  
  540.  
  541.  
  542.  
  543.  

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