Télécharger fuchpo.eso

Retour à la liste

Numérotation des lignes :

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

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