Télécharger fuchpo.eso

Retour à la liste

Numérotation des lignes :

  1. C FUCHPO SOURCE PV 13/04/12 21:15:40 7756
  2. SUBROUTINE FUCHPO(IP1,IP2,IRET)
  3. C======================================================================
  4. C fonction:
  5. C sous routine pour fusionner deux champs par points diffus
  6. C
  7. C arguments:
  8. C ip1 (E) pointeur sur le premier des deux champ par point
  9. C ip2 (E) pointeur sur le second des deux champ par point
  10. C iret (S) pointeur sur le champ par point resultat
  11. C
  12. C variables:
  13. C
  14. C * mtrav : - bb(i,j) est la valeur de la ieme inconnue de champ pour
  15. C le jieme noeud du tableau igeo .
  16. C - inco(nnin) contient le nom des nnin inconnues differentes
  17. C - ibin(i,j)=1 ou 0 indique que la ieme inconnue du champ
  18. C existe pour le jieme noeud du tableau igeo .
  19. C - igeo(i) est le numero a mettre dans un objet meleme pour
  20. C referencer le ieme noeud .
  21. C
  22. C * mtra : - nopoin(i) adresse de colonne dans bb et ibin des valeurs
  23. C correspondant au noeud i .
  24. C
  25. C * mtr1 : - ipcom liste des noms des inconnues permet de creer inco .
  26. C
  27. C * mtr2 : - iicom adresse dans ipcom des inconnues correspondant au
  28. C 2ieme ch point .
  29. C
  30. C * mtr3 : - index tableau de correspondance entre les supports geome-
  31. C triques du 1er chpoint et du 2ieme chpoint .
  32. C
  33. C auteur: A de Gayffier 13/06/94
  34. C======================================================================
  35. IMPLICIT INTEGER(I-N)
  36. IMPLICIT REAL*8(A-H,O-Z)
  37. -INC SMCHPOI
  38. -INC CCOPTIO
  39. -INC SMELEME
  40. -INC SMCOORD
  41. -INC TMTRAV
  42. SEGMENT/MTRA/(NOPOIN(XCOOR(/1)/(IDIM+1)))
  43. SEGMENT MTR1
  44. CHARACTER*4 IPCOM(0)
  45. ENDSEGMENT
  46. SEGMENT/MTR2/(IICOM(0))
  47. SEGMENT/MTR3/(INDEX(0))
  48. SEGMENT/MTR4/(IPHAR(0))
  49. C ordre de grandeur des composantes
  50. SEGMENT/MTR5/(DMOY(NNIN))
  51. C
  52. DIMENSION IPO(2)
  53. CHARACTER*8 MOT
  54. DATA UN,ZERO/1.D0,0.D0/
  55. C
  56. segact mcoord
  57. IRET = 0
  58. MCHPO1 = IP1
  59. MCHPO2 = IP2
  60. SEGACT MCHPO1,MCHPO2
  61. C
  62. C verification de la compatibilité des natures
  63. C
  64. NAT1 = MCHPO1.JATTRI(1)
  65. NAT2 = MCHPO2.JATTRI(1)
  66. C
  67. IF ( (NAT1*NAT2) .EQ. 0) THEN
  68. C une des deux nature est indeterminée
  69. CALL ERREUR(650)
  70. SEGDES MCHPO1,MCHPO2
  71. RETURN
  72. ELSE
  73. NANA = (NAT1 + NAT2) / NAT1 / NAT2
  74. IF ((NAT1 .EQ. 2) .AND. (NAT2 .EQ. 2)) THEN
  75. C les deux champ sont discrets: on somme les composantes communes
  76. CALL ADCHPO(IP1,IP2,IRET,1D0,1D0)
  77. RETURN
  78. ELSE IF ((NAT1 .NE. 1) .OR. (NAT2 .NE. 1)) THEN
  79. C les natures ne sont pas compatibles
  80. CALL ERREUR(649)
  81. RETURN
  82. ENDIF
  83. ENDIF
  84. C
  85. C les deux champ sont de nature diffuse
  86. C on moyenne les composantes communes
  87. C
  88. IF ( IP1 .NE. IP2) GOTO 60
  89. C
  90. C *** cas ou les 2 pointeurs ipo1 et ipo2 sont egaux
  91. C
  92. MCHPO1=IP1
  93. SEGACT MCHPO1
  94. NSOUPO=MCHPO1.IPCHP(/1)
  95. NAT=MCHPO1.JATTRI(/1)
  96. SEGINI MCHPOI
  97. DO 10 I=1,NAT
  98. JATTRI(I)=MCHPO1.JATTRI(I)
  99. 10 CONTINUE
  100. IRET=MCHPOI
  101. MOCHDE=MCHPO1.MOCHDE
  102. MTYPOI=MCHPO1.MTYPOI
  103. IFOPOI=IFOUR
  104. DO 50 IA=1,NSOUPO
  105. MSOUP1=MCHPO1.IPCHP(IA)
  106. SEGACT MSOUP1
  107. NC=MSOUP1.NOCOMP(/2)
  108. SEGINI MSOUPO
  109. IPCHP(IA)=MSOUPO
  110. DO 20 IB=1,NC
  111. NOCOMP(IB)=MSOUP1.NOCOMP(IB)
  112. NOHARM(IB)=MSOUP1.NOHARM(IB)
  113. 20 CONTINUE
  114. IGEOC=MSOUP1.IGEOC
  115. MPOVA1=MSOUP1.IPOVAL
  116. SEGACT MPOVA1
  117. N =MPOVA1.VPOCHA(/1)
  118. NC1=MPOVA1.VPOCHA(/2)
  119. IF(NC1.EQ.NC) GOTO 30
  120. C
  121. C erreur pb dimension tableau voir routine adchpo
  122. C
  123. IRET=0
  124. SEGDES MPOVA1,MSOUP1,MCHPO1
  125. SEGSUP MSOUPO,MCHPOI
  126. CALL ERREUR(114)
  127. C on sort de la sous routine
  128. GOTO 666
  129. C
  130. 30 CONTINUE
  131. SEGINI MPOVAL
  132. IPOVAL=MPOVAL
  133. DO 40 IC=1,NC
  134. DO 40 IB=1,N
  135. VPOCHA(IB,IC)=MPOVA1.VPOCHA(IB,IC)
  136. 40 CONTINUE
  137. SEGDES MPOVAL,MPOVA1,MSOUPO,MSOUP1
  138. 50 CONTINUE
  139. C
  140. SEGDES MCHPOI,MCHPO1
  141. C on sort de la sous routine
  142. GOTO 666
  143. C
  144. C *** cas ou les pointeurs ip1 et ip2 sont differents
  145. C
  146. 60 CONTINUE
  147. MCHPO1=IP1
  148. MCHPO2=IP2
  149. IPO(1)=IP1
  150. IPO(2)=IP2
  151. SEGACT MCHPO1,MCHPO2
  152. MOT=MCHPO1.MTYPOI
  153. IF(MOT.NE.MCHPO2.MTYPOI) THEN
  154. MOT='ET OU +'
  155. ENDIF
  156. SEGDES MCHPO1,MCHPO2
  157. C
  158. C on verifie que les nbres de sous paquets sont egaux
  159. C
  160. MCHPO1=IP1
  161. MCHPO2=IP2
  162. SEGACT MCHPO1,MCHPO2
  163. NSOUP1=MCHPO1.IPCHP(/1)
  164. NSOUP2=MCHPO2.IPCHP(/1)
  165. NAT1 = MCHPO1.JATTRI(/1)
  166. NAT2 = MCHPO2.JATTRI(/1)
  167. IF(NSOUP1.EQ.NSOUP2) GO TO 75
  168. SEGDES MCHPO1,MCHPO2
  169. C traitement par la methode générale
  170. GO TO 300
  171. C
  172. C on regarde si les supports geometriques sont identiques a une
  173. C permutation pres
  174. C
  175. 75 SEGINI MTR3
  176. DO 100 IA=1,NSOUP1
  177. MSOUP1=MCHPO1.IPCHP(IA)
  178. SEGACT MSOUP1
  179. DO 80 IB=1,NSOUP2
  180. MSOUP2=MCHPO2.IPCHP(IB)
  181. SEGACT MSOUP2
  182. IF(MSOUP1.IGEOC.EQ.MSOUP2.IGEOC) GO TO 90
  183. SEGDES MSOUP2
  184. 80 CONTINUE
  185. C
  186. C il n y a pas egalite des supports geometriques a une permutation
  187. C pres
  188. C
  189. SEGDES MSOUP1,MCHPO1,MCHPO2
  190. SEGSUP MTR3
  191. C traitement par la methode générale
  192. GO TO 300
  193. C
  194. 90 CONTINUE
  195. C la permutation est rangée dans index
  196. INDEX(**)=IB
  197. SEGDES MSOUP1,MSOUP2
  198. 100 CONTINUE
  199. C
  200. C *** cas ou il y a egalite des supports geometriques a une permutation
  201. C pres
  202. C
  203. NSOUPO=NSOUP1
  204. NAT=MAX(NAT1,NAT2,1)
  205. SEGINI MCHPOI
  206. JATTRI(1) = 1
  207. IRET=MCHPOI
  208. MTYPOI=MOT
  209. MOCHDE=MCHPO1.MOCHDE
  210. IFOPOI=IFOUR
  211. C
  212. DO 250 IA=1,NSOUP1
  213. MSOUP1=MCHPO1.IPCHP(IA)
  214. MSOUP2=MCHPO2.IPCHP(INDEX(IA))
  215. SEGACT MSOUP1,MSOUP2
  216. C
  217. C comparaison des noms des composantes
  218. C
  219. SEGINI MTR1,MTR4
  220. NC1=MSOUP1.NOCOMP(/2)
  221. NC2=MSOUP2.NOCOMP(/2)
  222. DO 130 IB=1,NC1
  223. IPCOM(**)=MSOUP1.NOCOMP(IB)
  224. IPHAR(**)=MSOUP1.NOHARM(IB)
  225. 130 CONTINUE
  226. SEGINI MTR2
  227. DO 160 IB=1,NC2
  228. DO 140 IC=1,NC1
  229. IF(MSOUP2.NOCOMP(IB).NE.MSOUP1.NOCOMP(IC)) GOTO 140
  230. IF(MSOUP2.NOHARM(IB).EQ.MSOUP1.NOHARM(IC)) GOTO 150
  231. 140 CONTINUE
  232. C la composante du IB n'est pas commune
  233. IPCOM(**)=MSOUP2.NOCOMP(IB)
  234. IPHAR(**)=MSOUP2.NOHARM(IB)
  235. IICOM(**)=IPCOM(/2)
  236. GO TO 160
  237. 150 CONTINUE
  238. C la composante est commune
  239. IICOM(**)=IC
  240. 160 CONTINUE
  241. C
  242. MPOVA1=MSOUP1.IPOVAL
  243. MPOVA2=MSOUP2.IPOVAL
  244. SEGACT MPOVA1,MPOVA2
  245. N1=MPOVA1.VPOCHA(/1)
  246. N2=MPOVA2.VPOCHA(/1)
  247. NCX1=MPOVA1.VPOCHA(/2)
  248. NCX2=MPOVA2.VPOCHA(/2)
  249. IF(NCX1.NE.NC1) GOTO 170
  250. IF(NCX2.NE.NC2) GOTO 170
  251. IF(N1.NE.N2) GOTO 170
  252. GOTO 180
  253. 170 CONTINUE
  254. SEGDES MSOUP1,MSOUP2,MPOVA1,MPOVA2,MCHPO1,MCHPO2
  255. SEGSUP MTR1,MTR2,MTR3,MCHPOI,MTR4
  256. C
  257. C pb avec les dimensions des chpoints
  258. C
  259. CALL ERREUR(114)
  260. IRET=0
  261. C on sort de la sous routine
  262. GOTO 666
  263. C
  264. 180 CONTINUE
  265. N=N1
  266. NC=IPCOM(/2)
  267. SEGINI MPOVAL
  268. NNIN = NC
  269. SEGINI MTR5
  270. C
  271. C mise a 0 de vpocha
  272. C
  273. * DO 190 IB=1,N
  274. * DO 190 IC=1,NC
  275. * VPOCHA(IB,IC)=ZERO
  276. * 190 CONTINUE
  277. C
  278. C addition des chpoints
  279. C
  280. C on place les valeurs de MCHPO1
  281. DO 210 IC=1,NC1
  282. DO 200 IB=1,N
  283. VPOCHA(IB,IC) = MPOVA1.VPOCHA(IB,IC)+VPOCHA(IB,IC)
  284. DMOY(IC) = DMOY(IC) + ABS(VPOCHA(IB,IC)/N)
  285. 200 CONTINUE
  286. IF (IIMPI.EQ.123)
  287. & write (IOIMP,*) ' ic dmoy(ic) ',ic,dmoy(ic)
  288. 210 CONTINUE
  289. C
  290. DO 230 IC=1,NC2
  291. IIC=IICOM(IC)
  292. DO 220 IB=1,N
  293. IF ( IIC .LE. NC1 ) THEN
  294. C il s'agit d'ne composante commune
  295. XX1 = MPOVA2.VPOCHA(IB,IC)
  296. XX2 = VPOCHA(IB,IIC)
  297. DXX = ABS ( XX2 - XX1)
  298. * SXX = MIN(ABS ( XX1 + XX2 ) / 2.D0,1.D-50)
  299. SXX = DMOY(IIC)
  300. IF ( (DXX / SXX) .LT. 1.D-3) THEN
  301. VPOCHA(IB,IIC)= ( XX1 + XX2 ) / 2.D0
  302. ELSE
  303. C les valeurs des champ diffus au meme point sont différentes
  304. IF (IIMPI.EQ.123)
  305. & write (IOIMP,*) xx1,xx2,SXX
  306. CALL ERREUR(651)
  307. SEGDES MSOUP1,MSOUP2,MPOVA1,MPOVA2,MCHPO1,MCHPO2
  308. SEGSUP MTR1,MTR2,MTR3,MCHPOI,MTR4
  309. C on sort
  310. GOTO 666
  311. ENDIF
  312. ELSE
  313. C composantes non communes
  314. VPOCHA(IB,IIC) = MPOVA2.VPOCHA(IB,IC)+VPOCHA(IB,IIC)
  315. ENDIF
  316. 220 CONTINUE
  317. 230 CONTINUE
  318. SEGDES MPOVA2
  319. SEGDES MPOVA1
  320. C
  321. SEGINI MSOUPO
  322. DO 240 IB=1,NC
  323. NOCOMP(IB)=IPCOM(IB)
  324. NOHARM(IB)=IPHAR(IB)
  325. 240 CONTINUE
  326. SEGSUP MTR1,MTR2,MTR4
  327. IPOVAL=MPOVAL
  328. SEGDES MPOVAL
  329. IPT2=MSOUP1.IGEOC
  330. **** SEGINI,IPT1=IPT2
  331. IPT1 = IPT2
  332. IGEOC=IPT1
  333. SEGDES IPT1
  334. IPCHP(IA)=MSOUPO
  335. SEGDES MSOUPO,MSOUP1,MSOUP2
  336. SEGSUP MTR5
  337. 250 CONTINUE
  338. C
  339. SEGSUP MTR3
  340. SEGDES MCHPOI,MCHPO1,MCHPO2
  341. C on sort
  342. GOTO 666
  343. C
  344. C *** cas ou les supports geometriques ne sont pas identiques
  345. C a une permutation pres
  346. C
  347. 300 CONTINUE
  348. C
  349. C **** a-t-on affaires a des champoints vides?
  350. C
  351. MCHPOI=IP1
  352. SEGACT MCHPOI
  353. NS1=IPCHP(/1)
  354. MCHPO2=IP2
  355. SEGACT MCHPO2
  356. NS2=MCHPO2.IPCHP(/1)
  357. IF(NS1*NS2.NE.0) GO TO 3001
  358. SEGDES MCHPOI,MCHPO2
  359. IF(NS1+NS2.NE.0) THEN
  360. C un seul des chpoints est vide
  361. IF(NS1.EQ.0) IP1=IP2
  362. CALL ECRCHA('GEOM')
  363. CALL ECROBJ('CHPOINT ',IP1)
  364. CALL COPIER
  365. CALL LIROBJ('CHPOINT',IRET,1,IRETOU)
  366. ELSE
  367. C les deux chpoints sont vides
  368. NSOUPO=0
  369. NAT=1
  370. SEGINI MCHPOI
  371. IRET = MCHPOI
  372. ENDIF
  373. RETURN
  374. 3001 CONTINUE
  375. C
  376. SEGINI MTRA,MTR1,MTR4
  377. C
  378. C mise a zero de nopoin
  379. C
  380. * DO 305 IA=1,XCOOR(/1)/(IDIM+1)
  381. * NOPOIN(IA)=0
  382. * 305 CONTINUE
  383. C
  384. MCHPOI=IP1
  385. SEGACT MCHPOI
  386. MSOUPO=IPCHP(1)
  387. SEGACT MSOUPO
  388. IPCOM(**)=NOCOMP(1)
  389. IPHAR(**)=NOHARM(1)
  390. SEGDES MSOUPO,MCHPOI
  391. NC=1
  392. IK=0
  393. C
  394. C creation de nopoin et de ipcom
  395. C
  396. DO 360 ICH=1,2
  397. MCHPOI=IPO(ICH)
  398. SEGACT MCHPOI
  399. NSOUPO=IPCHP(/1)
  400. C
  401. C boucle sur les sous references d 1 chpoint
  402. C
  403. DO 350 IA=1,NSOUPO
  404. MSOUPO=IPCHP(IA)
  405. SEGACT MSOUPO
  406. MELEME=IGEOC
  407. SEGACT MELEME
  408. NBNN =NUM(/1)
  409. NBELEM=NUM(/2)
  410. DO 310 IB=1,NBELEM
  411. K=NUM(1,IB)
  412. IF(NOPOIN(K).NE.0) GO TO 310
  413. IK=IK+1
  414. NOPOIN(K)=IK
  415. 310 CONTINUE
  416. SEGDES MELEME
  417. NCBBB=NOCOMP(/2)
  418. DO 330 IB=1,NCBBB
  419. NC=IPCOM(/2)
  420. DO 320 IC=1,NC
  421. IF(IPCOM(IC).NE.NOCOMP(IB)) GO TO 320
  422. IF(IPHAR(IC).EQ.NOHARM(IB)) GO TO 330
  423. 320 CONTINUE
  424. IPCOM(**)=NOCOMP(IB)
  425. IPHAR(**)=NOHARM(IB)
  426. NC=NC+1
  427. 330 CONTINUE
  428. SEGDES MSOUPO
  429. 350 CONTINUE
  430. SEGDES MCHPOI
  431. 360 CONTINUE
  432. C
  433. NNIN=NC
  434. NNNOE=IK
  435. SEGINI MTRAV
  436. C
  437. C initialisation a zero des tableaux
  438. C
  439. SEGINI MTR5
  440. C
  441. * DO 370 IB=1,NNNOE
  442. * DO 370 IA=1,NNIN
  443. * BB(IA,IB)=ZERO
  444. * IBIN(IA,IB)=0
  445. * IMOY(IA,IB) = 0
  446. * 370 CONTINUE
  447. C
  448. C creation de inco
  449. C
  450. DO 380 IA=1,NNIN
  451. INCO(IA)=IPCOM(IA)
  452. NHAR(IA)=IPHAR(IA)
  453. 380 CONTINUE
  454. C
  455. C creation de bb,ibin,igeo
  456. C
  457. DO 450 ICH=1,2
  458. MCHPOI=IPO(ICH)
  459. SEGACT MCHPOI
  460. NSOUPO=IPCHP(/1)
  461. DO 430 IA=1,NSOUPO
  462. MSOUPO=IPCHP(IA)
  463. SEGACT MSOUPO
  464. MELEME=IGEOC
  465. SEGACT MELEME
  466. MPOVAL=IPOVAL
  467. SEGACT MPOVAL
  468. NBELEM=NUM(/2)
  469. N=VPOCHA(/1)
  470. NC=VPOCHA(/2)
  471. NC1=NOCOMP(/2)
  472. C
  473. DO 420 IB=1,NC1
  474. DO 390 IC=1,NNIN
  475. IF(NOCOMP(IB).NE.IPCOM(IC)) GO TO 390
  476. IF(NOHARM(IB).EQ.IPHAR(IC)) GO TO 400
  477. 390 CONTINUE
  478. 400 CONTINUE
  479. DO 411 ID=1,NBELEM
  480. DMOY(IB)=DMOY(IB)+ABS(VPOCHA(ID,IB)/NBELEM)
  481. 411 CONTINUE
  482. DO 410 ID=1,NBELEM
  483. KI=NOPOIN(NUM(1,ID))
  484. IGEO(KI)=NUM(1,ID)
  485. IF ( IBIN(IC,KI) .EQ. 1) THEN
  486. C la valeur au point est defini dans les deux champs
  487. XX1 = BB(IC,KI)
  488. XX2 = VPOCHA(ID,IB)
  489. DXX = ABS ( XX2 - XX1 )
  490. SXX = DMOY(IB)
  491. IF ( DXX .LE. (1.D-4*SXX)) THEN
  492. BB(IC,KI) = ( XX1 + XX2 ) / 2.D0
  493. ELSE
  494. C les valeurs des champs au meme point sont différentes
  495. IF (IIMPI.EQ.123)
  496. & write (IOIMP,*) xx1,xx2,sxx
  497. CALL ERREUR (651)
  498. SEGSUP MTRAV,MTRA,MTR1,MTR4,MTR5
  499. GOTO 666
  500. ENDIF
  501. ELSE
  502. BB(IC,KI)=BB(IC,KI)+VPOCHA(ID,IB)
  503. * DMOY(IC) = DMOY(IC) +ABS(BB(IC,KI)/NNNOE)
  504. ENDIF
  505. IBIN(IC,KI)=1
  506. 410 CONTINUE
  507. 420 CONTINUE
  508. SEGDES MELEME,MSOUPO,MPOVAL
  509. 430 CONTINUE
  510. SEGDES MCHPOI
  511. 450 CONTINUE
  512. ITRAV=MTRAV
  513. SEGDES MTRAV
  514. C
  515. C reconstuction de la partition
  516. C
  517. CALL CRECHP(ITRAV,ICHPOI)
  518. C
  519. SEGSUP MTRAV,MTRA,MTR1,MTR4,MTR5
  520. IRET=ICHPOI
  521. MCHPOI=ICHPOI
  522. MCHPO1 = IP1
  523. MCHPO2 = IP2
  524. SEGACT MCHPOI,MCHPO1,MCHPO2
  525. NAT1 = MCHPO1.JATTRI(/1)
  526. NAT2 = MCHPO2.JATTRI(/1)
  527. NAT=MAX(NAT1,NAT2,1)
  528. NSOUPO = IPCHP(/1)
  529. SEGADJ MCHPOI
  530. IRET=MCHPOI
  531. MTYPOI=MOT
  532. JATTRI(1) = 1
  533. SEGDES MCHPOI,MCHPO1,MCHPO2
  534. C
  535. 666 CONTINUE
  536. RETURN
  537. END
  538.  
  539.  
  540.  
  541.  
  542.  
  543.  
  544.  
  545.  
  546.  
  547.  

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