Télécharger fuchpo.eso

Retour à la liste

Numérotation des lignes :

  1. C FUCHPO SOURCE CB215821 19/08/20 21:17:56 10287
  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. character*4 mcle(1)
  56. data mcle/'NOER'/
  57. C
  58. noer=0
  59. call lirmot(mcle,1,noer,0)
  60. segact mcoord
  61. IRET = 0
  62. MCHPO1 = IP1
  63. MCHPO2 = IP2
  64. SEGACT MCHPO1,MCHPO2
  65. C
  66. C verification de la compatibilité des natures
  67. C
  68. NAT1 = MCHPO1.JATTRI(1)
  69. NAT2 = MCHPO2.JATTRI(1)
  70. C
  71. IF ( (NAT1*NAT2) .EQ. 0) THEN
  72. C une des deux nature est indeterminée
  73. CALL ERREUR(650)
  74. RETURN
  75. ELSE
  76. NANA = (NAT1 + NAT2) / NAT1 / NAT2
  77. IF ((NAT1 .EQ. 2) .AND. (NAT2 .EQ. 2)) THEN
  78. C les deux champ sont discrets: on somme les composantes communes
  79. CALL ADCHPO(IP1,IP2,IRET,1D0,1D0)
  80. RETURN
  81. ELSE IF ((NAT1 .NE. 1) .OR. (NAT2 .NE. 1)) THEN
  82. C les natures ne sont pas compatibles
  83. CALL ERREUR(649)
  84. RETURN
  85. ENDIF
  86. ENDIF
  87. C
  88. C les deux champ sont de nature diffuse
  89. C on moyenne les composantes communes
  90. C
  91. IF ( IP1 .NE. IP2) GOTO 60
  92. C
  93. C *** cas ou les 2 pointeurs ipo1 et ipo2 sont egaux
  94. C
  95. MCHPO1=IP1
  96. SEGACT MCHPO1
  97. NSOUPO=MCHPO1.IPCHP(/1)
  98. NAT=MCHPO1.JATTRI(/1)
  99. SEGINI MCHPOI
  100. DO 10 I=1,NAT
  101. JATTRI(I)=MCHPO1.JATTRI(I)
  102. 10 CONTINUE
  103. IRET=MCHPOI
  104. MOCHDE=MCHPO1.MOCHDE
  105. MTYPOI=MCHPO1.MTYPOI
  106. IFOPOI=IFOUR
  107. DO 50 IA=1,NSOUPO
  108. MSOUP1=MCHPO1.IPCHP(IA)
  109. SEGACT MSOUP1
  110. NC=MSOUP1.NOCOMP(/2)
  111. SEGINI MSOUPO
  112. IPCHP(IA)=MSOUPO
  113. DO 20 IB=1,NC
  114. NOCOMP(IB)=MSOUP1.NOCOMP(IB)
  115. NOHARM(IB)=MSOUP1.NOHARM(IB)
  116. 20 CONTINUE
  117. IGEOC=MSOUP1.IGEOC
  118. MPOVA1=MSOUP1.IPOVAL
  119. SEGACT MPOVA1
  120. N =MPOVA1.VPOCHA(/1)
  121. NC1=MPOVA1.VPOCHA(/2)
  122. IF(NC1.EQ.NC) GOTO 30
  123. C
  124. C erreur pb dimension tableau voir routine adchpo
  125. C
  126. IRET=0
  127. SEGSUP MSOUPO,MCHPOI
  128. CALL ERREUR(114)
  129. C on sort de la sous routine
  130. GOTO 666
  131. C
  132. 30 CONTINUE
  133. SEGINI MPOVAL
  134. IPOVAL=MPOVAL
  135. DO 40 IC=1,NC
  136. DO 40 IB=1,N
  137. VPOCHA(IB,IC)=MPOVA1.VPOCHA(IB,IC)
  138. 40 CONTINUE
  139. 50 CONTINUE
  140. C
  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. C
  157. C on verifie que les nbres de sous paquets sont egaux
  158. C
  159. MCHPO1=IP1
  160. MCHPO2=IP2
  161. SEGACT MCHPO1,MCHPO2
  162. NSOUP1=MCHPO1.IPCHP(/1)
  163. NSOUP2=MCHPO2.IPCHP(/1)
  164. NAT1 = MCHPO1.JATTRI(/1)
  165. NAT2 = MCHPO2.JATTRI(/1)
  166. IF(NSOUP1.EQ.NSOUP2) GO TO 75
  167. C traitement par la methode générale
  168. GO TO 300
  169. C
  170. C on regarde si les supports geometriques sont identiques a une
  171. C permutation pres
  172. C
  173. 75 SEGINI MTR3
  174. DO 100 IA=1,NSOUP1
  175. MSOUP1=MCHPO1.IPCHP(IA)
  176. SEGACT MSOUP1
  177. DO 80 IB=1,NSOUP2
  178. MSOUP2=MCHPO2.IPCHP(IB)
  179. SEGACT MSOUP2
  180. IF(MSOUP1.IGEOC.EQ.MSOUP2.IGEOC) GO TO 90
  181. 80 CONTINUE
  182. C
  183. C il n y a pas egalite des supports geometriques a une permutation
  184. C pres
  185. C
  186. SEGSUP MTR3
  187. C traitement par la methode générale
  188. GO TO 300
  189. C
  190. 90 CONTINUE
  191. C la permutation est rangée dans index
  192. INDEX(**)=IB
  193. 100 CONTINUE
  194. C
  195. C *** cas ou il y a egalite des supports geometriques a une permutation
  196. C pres
  197. C
  198. NSOUPO=NSOUP1
  199. NAT=MAX(NAT1,NAT2,1)
  200. SEGINI MCHPOI
  201. JATTRI(1) = 1
  202. IRET=MCHPOI
  203. MTYPOI=MOT
  204. MOCHDE=MCHPO1.MOCHDE
  205. IFOPOI=IFOUR
  206. C
  207. DO 250 IA=1,NSOUP1
  208. MSOUP1=MCHPO1.IPCHP(IA)
  209. MSOUP2=MCHPO2.IPCHP(INDEX(IA))
  210. SEGACT MSOUP1,MSOUP2
  211. C
  212. C comparaison des noms des composantes
  213. C
  214. SEGINI MTR1,MTR4
  215. NC1=MSOUP1.NOCOMP(/2)
  216. NC2=MSOUP2.NOCOMP(/2)
  217. DO 130 IB=1,NC1
  218. IPCOM(**)=MSOUP1.NOCOMP(IB)
  219. IPHAR(**)=MSOUP1.NOHARM(IB)
  220. 130 CONTINUE
  221. SEGINI MTR2
  222. DO 160 IB=1,NC2
  223. DO 140 IC=1,NC1
  224. IF(MSOUP2.NOCOMP(IB).NE.MSOUP1.NOCOMP(IC)) GOTO 140
  225. IF(MSOUP2.NOHARM(IB).EQ.MSOUP1.NOHARM(IC)) GOTO 150
  226. 140 CONTINUE
  227. C la composante du IB n'est pas commune
  228. IPCOM(**)=MSOUP2.NOCOMP(IB)
  229. IPHAR(**)=MSOUP2.NOHARM(IB)
  230. IICOM(**)=IPCOM(/2)
  231. GO TO 160
  232. 150 CONTINUE
  233. C la composante est commune
  234. IICOM(**)=IC
  235. 160 CONTINUE
  236. C
  237. MPOVA1=MSOUP1.IPOVAL
  238. MPOVA2=MSOUP2.IPOVAL
  239. SEGACT MPOVA1,MPOVA2
  240. N1=MPOVA1.VPOCHA(/1)
  241. N2=MPOVA2.VPOCHA(/1)
  242. NCX1=MPOVA1.VPOCHA(/2)
  243. NCX2=MPOVA2.VPOCHA(/2)
  244. IF(NCX1.NE.NC1) GOTO 170
  245. IF(NCX2.NE.NC2) GOTO 170
  246. IF(N1.NE.N2) GOTO 170
  247. GOTO 180
  248. 170 CONTINUE
  249. SEGSUP MTR1,MTR2,MTR3,MCHPOI,MTR4
  250. C
  251. C pb avec les dimensions des chpoints
  252. C
  253. CALL ERREUR(114)
  254. IRET=0
  255. C on sort de la sous routine
  256. GOTO 666
  257. C
  258. 180 CONTINUE
  259. N=N1
  260. NC=IPCOM(/2)
  261. SEGINI MPOVAL
  262. NNIN = NC
  263. SEGINI MTR5
  264. C
  265. C mise a 0 de vpocha
  266. C
  267. * DO 190 IB=1,N
  268. * DO 190 IC=1,NC
  269. * VPOCHA(IB,IC)=ZERO
  270. * 190 CONTINUE
  271. C
  272. C addition des chpoints
  273. C
  274. C on place les valeurs de MCHPO1
  275. DO 210 IC=1,NC1
  276. DO 200 IB=1,N
  277. VPOCHA(IB,IC) = MPOVA1.VPOCHA(IB,IC)+VPOCHA(IB,IC)
  278. DMOY(IC) = DMOY(IC) + ABS(VPOCHA(IB,IC)/N)
  279. 200 CONTINUE
  280. IF (IIMPI.EQ.123)
  281. & write (IOIMP,*) ' ic dmoy(ic) ',ic,dmoy(ic)
  282. 210 CONTINUE
  283. C
  284. DO 230 IC=1,NC2
  285. IIC=IICOM(IC)
  286. DO 220 IB=1,N
  287. IF ( IIC .LE. NC1 ) THEN
  288. C il s'agit d'ne composante commune
  289. XX1 = MPOVA2.VPOCHA(IB,IC)
  290. XX2 = VPOCHA(IB,IIC)
  291. DXX = ABS ( XX2 - XX1)
  292. * SXX = MIN(ABS ( XX1 + XX2 ) / 2.D0,1.D-50)
  293. SXX = DMOY(IIC)
  294. IF ( (DXX / SXX) .LT. 1.D-3 .or.noer.eq.1) THEN
  295. VPOCHA(IB,IIC)= ( XX1 + XX2 ) / 2.D0
  296. ELSE
  297. C les valeurs des champ diffus au meme point sont différentes
  298. IF (IIMPI.EQ.123)
  299. & write (IOIMP,*) xx1,xx2,SXX
  300. CALL ERREUR(651)
  301. SEGSUP MTR1,MTR2,MTR3,MCHPOI,MTR4
  302. C on sort
  303. GOTO 666
  304. ENDIF
  305. ELSE
  306. C composantes non communes
  307. VPOCHA(IB,IIC) = MPOVA2.VPOCHA(IB,IC)+VPOCHA(IB,IIC)
  308. ENDIF
  309. 220 CONTINUE
  310. 230 CONTINUE
  311. C
  312. SEGINI MSOUPO
  313. DO 240 IB=1,NC
  314. NOCOMP(IB)=IPCOM(IB)
  315. NOHARM(IB)=IPHAR(IB)
  316. 240 CONTINUE
  317. SEGSUP MTR1,MTR2,MTR4
  318. IPOVAL=MPOVAL
  319. IPT2=MSOUP1.IGEOC
  320. **** SEGINI,IPT1=IPT2
  321. IPT1 = IPT2
  322. IGEOC=IPT1
  323. IPCHP(IA)=MSOUPO
  324. SEGSUP MTR5
  325. 250 CONTINUE
  326. C
  327. SEGSUP MTR3
  328. C on sort
  329. GOTO 666
  330. C
  331. C *** cas ou les supports geometriques ne sont pas identiques
  332. C a une permutation pres
  333. C
  334. 300 CONTINUE
  335. C
  336. C **** a-t-on affaires a des champoints vides?
  337. C
  338. MCHPOI=IP1
  339. SEGACT MCHPOI
  340. NS1=IPCHP(/1)
  341. MCHPO2=IP2
  342. SEGACT MCHPO2
  343. NS2=MCHPO2.IPCHP(/1)
  344. IF(NS1*NS2.NE.0) GO TO 3001
  345. IF(NS1+NS2.NE.0) THEN
  346. C un seul des chpoints est vide
  347. IF(NS1.EQ.0) IP1=IP2
  348. CALL ECRCHA('GEOM')
  349. CALL ECROBJ('CHPOINT ',IP1)
  350. CALL COPIER
  351. CALL LIROBJ('CHPOINT',IRET,1,IRETOU)
  352. ELSE
  353. C les deux chpoints sont vides
  354. NSOUPO=0
  355. NAT=1
  356. SEGINI MCHPOI
  357. IRET = MCHPOI
  358. ENDIF
  359. RETURN
  360. 3001 CONTINUE
  361. C
  362. SEGINI MTRA,MTR1,MTR4
  363. C
  364. C mise a zero de nopoin
  365. C
  366. * DO 305 IA=1,XCOOR(/1)/(IDIM+1)
  367. * NOPOIN(IA)=0
  368. * 305 CONTINUE
  369. C
  370. MCHPOI=IP1
  371. SEGACT MCHPOI
  372. MSOUPO=IPCHP(1)
  373. SEGACT MSOUPO
  374. IPCOM(**)=NOCOMP(1)
  375. IPHAR(**)=NOHARM(1)
  376. NC=1
  377. IK=0
  378. C
  379. C creation de nopoin et de ipcom
  380. C
  381. DO 360 ICH=1,2
  382. MCHPOI=IPO(ICH)
  383. SEGACT MCHPOI
  384. NSOUPO=IPCHP(/1)
  385. C
  386. C boucle sur les sous references d 1 chpoint
  387. C
  388. DO 350 IA=1,NSOUPO
  389. MSOUPO=IPCHP(IA)
  390. SEGACT MSOUPO
  391. MELEME=IGEOC
  392. SEGACT MELEME
  393. NBNN =NUM(/1)
  394. NBELEM=NUM(/2)
  395. DO 310 IB=1,NBELEM
  396. K=NUM(1,IB)
  397. IF(NOPOIN(K).NE.0) GO TO 310
  398. IK=IK+1
  399. NOPOIN(K)=IK
  400. 310 CONTINUE
  401. NCBBB=NOCOMP(/2)
  402. DO 330 IB=1,NCBBB
  403. NC=IPCOM(/2)
  404. DO 320 IC=1,NC
  405. IF(IPCOM(IC).NE.NOCOMP(IB)) GO TO 320
  406. IF(IPHAR(IC).EQ.NOHARM(IB)) GO TO 330
  407. 320 CONTINUE
  408. IPCOM(**)=NOCOMP(IB)
  409. IPHAR(**)=NOHARM(IB)
  410. NC=NC+1
  411. 330 CONTINUE
  412. 350 CONTINUE
  413. 360 CONTINUE
  414. C
  415. NNIN=NC
  416. NNNOE=IK
  417. SEGINI MTRAV
  418. C
  419. C initialisation a zero des tableaux
  420. C
  421. SEGINI MTR5
  422. C
  423. * DO 370 IB=1,NNNOE
  424. * DO 370 IA=1,NNIN
  425. * BB(IA,IB)=ZERO
  426. * IBIN(IA,IB)=0
  427. * IMOY(IA,IB) = 0
  428. * 370 CONTINUE
  429. C
  430. C creation de inco
  431. C
  432. DO 380 IA=1,NNIN
  433. INCO(IA)=IPCOM(IA)
  434. NHAR(IA)=IPHAR(IA)
  435. 380 CONTINUE
  436. C
  437. C creation de bb,ibin,igeo
  438. C
  439. DO 450 ICH=1,2
  440. MCHPOI=IPO(ICH)
  441. SEGACT MCHPOI
  442. NSOUPO=IPCHP(/1)
  443. DO 430 IA=1,NSOUPO
  444. MSOUPO=IPCHP(IA)
  445. SEGACT MSOUPO
  446. MELEME=IGEOC
  447. SEGACT MELEME
  448. MPOVAL=IPOVAL
  449. SEGACT MPOVAL
  450. NBELEM=NUM(/2)
  451. N=VPOCHA(/1)
  452. NC=VPOCHA(/2)
  453. NC1=NOCOMP(/2)
  454. C
  455. DO 420 IB=1,NC1
  456. DO 390 IC=1,NNIN
  457. IF(NOCOMP(IB).NE.IPCOM(IC)) GO TO 390
  458. IF(NOHARM(IB).EQ.IPHAR(IC)) GO TO 400
  459. 390 CONTINUE
  460. 400 CONTINUE
  461. DO 411 ID=1,NBELEM
  462. DMOY(IB)=DMOY(IB)+ABS(VPOCHA(ID,IB)/NBELEM)
  463. 411 CONTINUE
  464. DO 410 ID=1,NBELEM
  465. KI=NOPOIN(NUM(1,ID))
  466. IGEO(KI)=NUM(1,ID)
  467. IF ( IBIN(IC,KI) .EQ. 1) THEN
  468. C la valeur au point est defini dans les deux champs
  469. XX1 = BB(IC,KI)
  470. XX2 = VPOCHA(ID,IB)
  471. DXX = ABS ( XX2 - XX1 )
  472. SXX = DMOY(IB)
  473. IF ( DXX .LE. (1.D-4*SXX).or.noer.eq.1) THEN
  474. BB(IC,KI) = ( XX1 + XX2 ) / 2.D0
  475. ELSE
  476. C les valeurs des champs au meme point sont différentes
  477. IF (IIMPI.EQ.123)
  478. & write (IOIMP,*) xx1,xx2,sxx
  479. CALL ERREUR (651)
  480. SEGSUP MTRAV,MTRA,MTR1,MTR4,MTR5
  481. GOTO 666
  482. ENDIF
  483. ELSE
  484. BB(IC,KI)=BB(IC,KI)+VPOCHA(ID,IB)
  485. * DMOY(IC) = DMOY(IC) +ABS(BB(IC,KI)/NNNOE)
  486. ENDIF
  487. IBIN(IC,KI)=1
  488. 410 CONTINUE
  489. 420 CONTINUE
  490. 430 CONTINUE
  491. 450 CONTINUE
  492. ITRAV=MTRAV
  493. C
  494. C reconstuction de la partition
  495. C
  496. CALL CRECHP(ITRAV,ICHPOI)
  497. C
  498. SEGSUP MTRAV,MTRA,MTR1,MTR4,MTR5
  499. IRET=ICHPOI
  500. MCHPOI=ICHPOI
  501. MCHPO1 = IP1
  502. MCHPO2 = IP2
  503. NAT1 = MCHPO1.JATTRI(/1)
  504. NAT2 = MCHPO2.JATTRI(/1)
  505. NAT=MAX(NAT1,NAT2,1)
  506. NSOUPO = IPCHP(/1)
  507. SEGADJ MCHPOI
  508. IRET=MCHPOI
  509. MTYPOI=MOT
  510. JATTRI(1) = 1
  511. C
  512. 666 CONTINUE
  513. RETURN
  514. END
  515.  
  516.  
  517.  
  518.  
  519.  

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