Télécharger biosav.eso

Retour à la liste

Numérotation des lignes :

  1. C BIOSAV SOURCE CHAT 05/11/29 21:15:02 5153
  2. SUBROUTINE BIOSAV
  3. *----------------------------------------------------------
  4. *
  5. * CALCUL DES CHAMPS DE BIOT ET SAVART
  6. *
  7. *----------------------------------------------------------
  8. IMPLICIT INTEGER(I-N)
  9. IMPLICIT REAL*8(A-H,O-Z)
  10. -INC CCOPTIO
  11. -INC CCREEL
  12. -INC SMELEME
  13. -INC SMCOORD
  14. -INC TMTRAV
  15. -INC SMCHPOI
  16. INTEGER ENT
  17. DIMENSION C(3),PN1(3),RESP(3)
  18. DIMENSION VX(3),VZ(3),VY(3),VI(3),PP3(3)
  19. DIMENSION ZF1(1),ZALF(1),ZDENS(1),ZBET(1)
  20. DIMENSION P1(3),P2(3),P3(3),ROTA(3,3)
  21. SEGMENT ICPR(XCOOR(/1)/(IDIM+1))
  22. SEGMENT SANGLE
  23. REAL*8 TETM(NT4)
  24. REAL*8 TETI(NT4)
  25. REAL*8 DTE(NT4)
  26. ENDSEGMENT
  27. *
  28. CHARACTER*4 MCLE(4),CHSECT(1),CHOIX(2)
  29. LOGICAL INDU,POTE
  30. *
  31. DATA EPSIL3/1.D-12/,EPSIL2/1.D-8/
  32. *
  33. DATA MCLE/'CERC','ARC ','BARR','FIL '/
  34. DATA CHSECT/'TRAP'/
  35. DATA CHOIX/'INDU','POTE'/
  36. *
  37. CALL LIROBJ('CHPOINT ',MCHPO1,0,IRETOU)
  38. IF (IERR.NE.0) RETURN
  39. IF (MCHPO1.GT.0) THEN
  40. CALL LIROBJ('MAILLAGE',IPT1,1,IRETOU)
  41. IF (IERR.NE.0) RETURN
  42.  
  43. * OPTION CALCUL PAR INTEGRALE ELLIPTIQUE (PHR 2002)
  44.  
  45. CALL BIOKAM (MCHPO1,IPT1,ENT)
  46. IF (IERR.NE.0) RETURN
  47. MCHPOI = ENT
  48. GOTO 9990
  49. ENDIF
  50. * CHOIX DE L'OPTION DE CALCUL
  51. * (INDUCTION PAR DEFAUT)
  52. *
  53. INDU=.TRUE.
  54. POTE=.FALSE.
  55. CALL LIRMOT(CHOIX,2,IMC,0)
  56. IF (IERR.NE.0) RETURN
  57. IF(IMC.EQ.2) THEN
  58. POTE=.TRUE.
  59. INDU=.FALSE.
  60. ENDIF
  61. *
  62. * LECTURE D'UN MAILLAGE
  63. *
  64. CALL LIROBJ('MAILLAGE',MELEME,1,IRETOU)
  65. IF (IERR.NE.0) RETURN
  66. *
  67. * LECTURE D'UN MOT CLE
  68. *
  69. CALL LIRMOT(MCLE,4,IMC,1)
  70. IF (IERR.NE.0) RETURN
  71. *
  72. * IMC = 1 ON A LU CERCLE
  73. * IMC = 2 ON A LU ARC
  74. * IMC = 3 ON A LU BARRE
  75. * IMC = 4 ON A LU FIL
  76. *
  77. *
  78. * LECTURE DE 2 ou 3 POINTS
  79. *
  80. CALL LIROBJ('POINT ',IP1,1,IRETOU)
  81. CALL LIROBJ('POINT ',IP2,1,IRETOU)
  82. IP3=0
  83. IF(IMC.LT.4) THEN
  84. CALL LIROBJ('POINT ',IP3,1,IRETOU)
  85. ELSE
  86. C CAS DU FIL RECTILIGNE
  87. GO TO 10
  88. ENDIF
  89. IF (IERR.NE.0) RETURN
  90. *
  91. * LECTURE DE 4 OU 5 FLOTTANTS SELON LES CAS
  92. *
  93. CALL LIRREE(F1,1,IRETOU)
  94. CALL LIRREE(F2,1,IRETOU)
  95. F3=0.D0
  96. IF(IMC.LT.3) CALL LIRREE(F3,1,IRETOU)
  97. *
  98. * LECTURE DU TYPE DE SECTION
  99. *
  100. ISECT=0
  101. CALL LIRMOT(CHSECT,1,ISECT,0)
  102. IF(ISECT.EQ.1) THEN
  103. CALL LIRREE(PENT1,1,IRETOU)
  104. CALL LIRREE(PENT2,1,IRETOU)
  105. IF(IERR.NE.0) RETURN
  106. ELSE
  107. ISECT=3
  108. ENDIF
  109. 10 CONTINUE
  110. CALL LIRREE(DENS,1,IRETOU)
  111. CALL LIRREE(XIU,1,IRETOU)
  112. C ON VIENT DE LIRE MU0 DANS L UNITE UTILISATEUR
  113. RUMKSA= 4.D0 * XPI * 1.D-7
  114. RAP=(XIU / RUMKSA)
  115. IU= 1
  116. IF (IERR.NE.0) RETURN
  117. ZF1(1)=F1
  118. ZDENS(1)=DENS
  119. *
  120. * EXEMPLE DE RECUPERATION DES COORDONNEES DU POINT P1
  121. *
  122. IREF=(IDIM+1)*IP1
  123. P1(1)=XCOOR(IREF-IDIM)
  124. P1(2)=XCOOR(IREF-IDIM+1)
  125. P1(3)=0.D0
  126. IF (IDIM.EQ.3) P1(3)=XCOOR(IREF-IDIM+2)
  127. IREF=(IDIM+1)*IP2
  128. P2(1)=XCOOR(IREF-IDIM)
  129. P2(2)=XCOOR(IREF-IDIM+1)
  130. P2(3)=0.D0
  131. IF (IDIM.EQ.3) P2(3)=XCOOR(IREF-IDIM+2)
  132. IF(IP3.NE.0) THEN
  133. IREF=(IDIM+1)*IP3
  134. P3(1)=XCOOR(IREF-IDIM)
  135. P3(2)=XCOOR(IREF-IDIM+1)
  136. P3(3)=0.D0
  137. IF (IDIM.EQ.3) P3(3)=XCOOR(IREF-IDIM+2)
  138. ENDIF
  139. *
  140. * -----------------------
  141. * | LES CALCULS .... |
  142. * -----------------------
  143. *
  144. * ON COMMENCE PAR EXPLORER LE MAILLAGE
  145. *
  146. SEGACT MELEME
  147. IF (ITYPEL.NE.1.OR.LISOUS(/1).NE.0) CALL CHANGE(MELEME,1)
  148. *
  149. * ON VERIFIE QUE CHAQUE POINT N'APPARAIT QU'UNE SEULE FOIS
  150. *
  151. NNNOE=NUM(/2)
  152. SEGINI ICPR
  153. ICON=0
  154. DO 12 I=1,NNNOE
  155. IKI=NUM(1,I)
  156. IF(ICPR(IKI).NE.0) GO TO 12
  157. ICON=ICON+1
  158. ICPR(IKI)=ICON
  159. 12 CONTINUE
  160. SEGSUP ICPR
  161. *
  162. IF(ICON.NE.NNNOE) THEN
  163. CALL ERREUR(75)
  164. SEGDES MELEME
  165. RETURN
  166. ENDIF
  167. *
  168. * ON INITIALISE LES NOMS DES COMPOSANTES ET LES HARMONIQUES
  169. * POUR L'INSTANT LIMITE AU 3D ET 2D PLAN
  170. *
  171. *
  172. NNIN=IDIM
  173. SEGINI MTRAV
  174. NHRM=NIFOUR
  175. IF(NIFOUR.EQ.1) NHRM=IFOUR
  176. INCO(1)='BX '
  177. IF(POTE) INCO(1)='AX '
  178. INCO(2)='BY '
  179. IF(POTE) INCO(2)='AY '
  180. NHAR(1)=NHRM
  181. NHAR(2)=NHRM
  182. IF(IDIM.EQ.3) THEN
  183. INCO(3)='BZ '
  184. IF(POTE) INCO(3)='AZ '
  185. NHAR(3)=NHRM
  186. ENDIF
  187. C
  188. GOTO (100,100,200,300),IMC
  189. 100 CONTINUE
  190. C BOBINE CIRCULAIRE
  191.  
  192. C ETABLISSEMENT DU REPERE LOCAL
  193. C AXE DES Z (P2 - P1 ) * (P3 - P1)
  194. DO 801 I=1,3
  195. VX(I)=P2(I)- P1(I)
  196. VI(I)=P3(I)- P1(I)
  197. 801 CONTINUE
  198. CALL PROVEC(VX,VI,VZ)
  199. CALL PROVEC(VZ,VX,VY)
  200. CALL ROTREP (ROTA,VX,VY,VZ)
  201. IF(IERR.NE.0) THEN
  202. SEGDES MELEME
  203. SEGSUP MTRAV
  204. RETURN
  205. ENDIF
  206. NB=1
  207. IF(F1.EQ.0.D0) RETURN
  208. ALF=F2/F1
  209. BET=F3/(2.D0*F1)
  210. ZALF(1)=ALF
  211. ZBET(1)=BET
  212. IC= 3
  213. DO 802 J=1,3
  214. PP3(J)= P3(J) - P1(J)
  215. 802 CONTINUE
  216. RAYO2=SDOT(3,VX,1,VX,1)
  217. PROD=SDOT(3,VX,1,VI,1)
  218. RAYO2=RAYO2+EPSIL3
  219. TM=ACOS(PROD/RAYO2)
  220. GOTO 350
  221. 200 CONTINUE
  222. C CAS DES BARRES ON VA TOURNER / OX POUR AMENER OY // P1P2
  223. GY=0.D0
  224. DO 803 I=1,3
  225. VX(I)=P2(I)- P1(I)
  226. VI(I)=P3(I)- P1(I)
  227. GY=GY +( VX(I)*VX(I))
  228. C(I)= (P2(I)+ P1(I))/2.D0
  229. 803 CONTINUE
  230. GY= SQRT( GY)
  231. CALL PROVEC(VX,VI,VZ)
  232. CALL PROVEC(VZ,VX,VY)
  233. CALL ROTREP (ROTA,VX,VY,VZ)
  234. IF(IERR.NE.0) THEN
  235. SEGDES MELEME
  236. SEGSUP MTRAV
  237. RETURN
  238. ENDIF
  239. DO 804 I= 1,3
  240. VZ(I)=0.D0
  241. DO 8041 J= 1,3
  242. VZ(I)= VZ(I) + ROTA(I,J) * C(J)
  243. 8041 CONTINUE
  244. 804 CONTINUE
  245. * CALL ROTATI (C,1,1,1)
  246. GO TO 350
  247. C
  248. C CAS DU FIL
  249. 300 CONTINUE
  250. C
  251. 350 CONTINUE
  252. UPSI =1.D-9
  253. *
  254. * ON BOUCLE SUR LES NOEUDS
  255. *
  256. NT4=8000
  257. SEGINI SANGLE
  258. DO 30 IPT=1,NNNOE
  259. L=NUM(1,IPT)
  260. IREF= (IDIM+1)*L
  261. * L EST LE NUMERO DU NOEUD DANS LE TABLEAU DES COORDONNEES
  262. DO 805 J=1,3
  263. PN1(J)=XCOOR(IREF-IDIM-1+J)
  264. 805 CONTINUE
  265. *
  266. GOTO (400,400,500,600),IMC
  267. 400 CONTINUE
  268. * CAS DES CERCLES ET DES ARC
  269. C MESURE DU POINT DANS LE REPERE LOCAL
  270. DO 806 J=1,3
  271. PN1(J)=PN1(J) - P1(J)
  272. 806 CONTINUE
  273. DO 807 I= 1,3
  274. VZ(I)=0.D0
  275. VY(I)=0.D0
  276. DO 8071 J= 1,3
  277. VY(I)= VY(I) + ROTA(I,J) * PN1(J)
  278. VZ(I)= VZ(I) + ROTA(I,J) * PP3(J)
  279. 8071 CONTINUE
  280. 807 CONTINUE
  281. * CALL ROTATI (PP3,1,1,1)
  282. * CALL ROTATI (PN1,1,1,1)
  283. C ON CHANGE DE RERERENTIEL
  284. C MESURE DANS LE REPERE LOCAL DES POINTS DE CALCUL DE HS
  285.  
  286. ZP= VY(3)
  287. RP= SQRT ( VY(1)**2 + VY(2)**2 )
  288. IF(RP.EQ.0.D0) THEN
  289. TMIN=0.D0
  290. COST=1.D0
  291. SINT=0.D0
  292. RP=EPSIL2
  293. ELSE
  294. TMIN=ACOS(VY(1)/(RP+EPSIL3))
  295. COST=COS(TMIN)
  296. SINT=VY(2)/(RP+EPSIL3)
  297. ENDIF
  298. DO 808 J= 1,3
  299. C(J)=0.D0
  300. 808 CONTINUE
  301. C
  302. IF(BET.EQ.0.D0.AND.ALF.EQ.1.D0) THEN
  303. C CAS DU FIL CIRCULAIRE
  304. RAYON=F1
  305. COUR=DENS
  306. C ON CONSTRUIT LE REPERE LOCAL
  307. sigt=1.d0
  308. if( abs(sint).gt.1.e-6)SIGT=SIGN(1.D0,-SINT)
  309. TMIN=SIGT*(TMIN-XPI)+XPI
  310. IF(IMC.EQ.1) TM=2.D0*XPI
  311. TMAX=TMIN+TM
  312. CALL TEGRAL(SANGLE)
  313. IF(INDU) THEN
  314. CALL BFILCI(SANGLE,COUR,RAYON,TMIN,TMAX,RP,ZP,BX,BY,BZ)
  315. ENDIF
  316. IF(POTE) THEN
  317. CALL AFILCI(SANGLE,COUR,RAYON,TMIN,TMAX,RP,ZP,BX,BY,BZ)
  318. ENDIF
  319. RESP(1)=BX*COST-BY*SINT
  320. RESP(2)=BX*SINT+BY*COST
  321. RESP(3)=BZ
  322. ELSE IF(BET.EQ.0.D0.AND.ISECT.EQ.3) THEN
  323. C CAS DE LA PLAQUE CIRCULAIRE
  324. RINT=F1
  325. REXT=F2
  326. DR=REXT-RINT
  327. COUR=DENS*DR
  328. C ON CONSTRUIT LE REPERE LOCAL
  329. sigt=1.d0
  330. if( abs(sint).gt.1.e-6) SIGT=SIGN(1.D0,-SINT)
  331. TMIN=SIGT*(TMIN-XPI)+XPI
  332. IF(IMC.EQ.1) TM=2.D0*XPI
  333. TMAX=TMIN+TM
  334. CALL TEGRAL(SANGLE)
  335. IF(INDU) THEN
  336. CALL BPLQCI(SANGLE,COUR,RINT,REXT,TMIN,TMAX,RP,ZP,
  337. & BX,BY,BZ)
  338. ENDIF
  339. IF(POTE) THEN
  340. CALL APLQCI(SANGLE,COUR,RINT,REXT,TMIN,TMAX,RP,ZP,
  341. & BX,BY,BZ)
  342. ENDIF
  343. RESP(1)=BX*COST-BY*SINT
  344. RESP(2)=BX*SINT+BY*COST
  345. RESP(3)=BZ
  346. ELSE IF(ALF.EQ.1D0) THEN
  347. C CAS DU CYLINDRE
  348. RAYON=F1
  349. H=F3
  350. COUR=DENS*H
  351. C ON CONSTRUIT LE REPERE LOCAL
  352. sigt=1.d0
  353. if( abs(sint).gt.1.e-6)SIGT=SIGN(1.D0,-SINT)
  354. TMIN=SIGT*(TMIN-XPI)+XPI
  355. IF(IMC.EQ.1) TM=2.D0*XPI
  356. TMAX=TMIN+TM
  357. CALL TEGRAL(SANGLE)
  358. IF(INDU) THEN
  359. CALL BCYLCI(SANGLE,COUR,RAYON,H,TMIN,TMAX,RP,ZP,
  360. & BX,BY,BZ)
  361. ENDIF
  362. IF(POTE) THEN
  363. CALL ACYLCI(SANGLE,COUR,RAYON,H,TMIN,TMAX,RP,ZP,
  364. & BX,BY,BZ)
  365. ENDIF
  366. RESP(1)=BX*COST-BY*SINT
  367. RESP(2)=BX*SINT+BY*COST
  368. RESP(3)=BZ
  369. ELSE
  370. IF (IMC.EQ.1.AND.ISECT.EQ.3) THEN
  371. IF(INDU) THEN
  372. CALL CHAMBO(IU,IC,NB,C(3),ZF1,ZALF,ZBET,ZDENS,ZP,RP,BZ
  373. $ ,BR,BP)
  374. IF (RP.LE.UPSI) THEN
  375. RESP(1)=0.D0
  376. RESP(2)=0.D0
  377. ELSE
  378. RESP(1)=BR * VY(1) / RP
  379. RESP(2)=BR * VY(2) / RP
  380. END IF
  381. RESP(3)=BZ
  382. GO TO 700
  383. ENDIF
  384. IF(POTE) THEN
  385. RINT=F1
  386. REXT=F2
  387. HCEN=F3
  388. WIDTH=F2-F1
  389. AIRE=F3*WIDTH
  390. COUR=DENS*AIRE
  391. DELTA=0.D0
  392. C ON CONSTRUIT LE REPERE LOCAL
  393. SIGT=1.D0
  394. if(ABS(SINT).gt.1.e-8)SIGT=SIGN(1.D0,-SINT)
  395. TMIN=SIGT*(TMIN-XPI)+XPI
  396. TM=2.D0*XPI
  397. TMAX=TMIN+TM
  398. PENT1=0.D0
  399. PENT2=0.D0
  400. CALL TEGRAL(SANGLE)
  401. CALL AARCTR(SANGLE,COUR,HCEN,HCEN,RINT,REXT,
  402. * PENT1,PENT2,TMIN,TMAX,
  403. * RP,ZP,BX,BY,BZ)
  404. RESP(1)=BX*COST-BY*SINT
  405. RESP(2)=BX*SINT+BY*COST
  406. RESP(3)=BZ
  407. GO TO 700
  408. ENDIF
  409. ENDIF
  410. IF(ISECT.EQ.3) THEN
  411. C CAS DES SECTIONS RECTANGULAIRES
  412. C ON CALCULE POUR UN ARC ---> ANGLES
  413. A1=0.D0
  414. IF (ABS(VZ(1)).LT.UPSI) THEN
  415. A2=XPI/2.D0
  416. ELSE
  417. AA= VZ(2)/VZ(1)
  418. A2=ATAN ( AA )
  419. IF ( A2.LT.0.D0) A2 = XPI + A2
  420. END IF
  421. IF(INDU) THEN
  422. CALL CHAMAR(IU,IC,C,A1,A2,F1,ALF,BET,DENS,
  423. * VY(1),VY(2),VY(3),BX,BY,BZ)
  424. RESP(1)=BX
  425. RESP(2)=BY
  426. RESP(3)=BZ
  427. ENDIF
  428. IF(POTE) THEN
  429. RINT=F1
  430. REXT=F2
  431. HCEN=F3
  432. WIDTH=F2-F1
  433. AIRE=F3*WIDTH
  434. COUR=DENS*AIRE
  435. PENT1=0.D0
  436. PENT2=0.D0
  437. C ON CONSTRUIT LE REPERE LOCAL
  438. sigt=1.d0
  439. if( abs(sint).gt.1.e-6)SIGT=SIGN(1.d0,-SINT)
  440. TMIN=SIGT*(TMIN-XPI)+XPI
  441. TMAX=TMIN+TM
  442. CALL TEGRAL(SANGLE)
  443. CALL AARCTR(SANGLE,COUR,HCEN,HCEN,RINT,REXT,
  444. * PENT1,PENT2,TMIN,TMAX,
  445. * RP,ZP,BX,BY,BZ)
  446. RESP(1)=BX*COST-BY*SINT
  447. RESP(2)=BX*SINT+BY*COST
  448. RESP(3)=BZ
  449. ENDIF
  450. ELSEIF(ISECT.EQ.1) THEN
  451. C CAS DES SECTIONS TRAPEZOIDALES
  452. RINT=F1
  453. REXT=F2
  454. HCEN=F3
  455. WIDTH=F2-F1
  456. AIRE=F3*WIDTH
  457. COUR=DENS*AIRE
  458. DELTA=0.5D0*WIDTH*(PENT2-PENT1)
  459. HINT=HCEN-DELTA
  460. HEXT=HCEN+DELTA
  461. C ON CONSTRUIT LE REPERE LOCAL
  462. SIGT=1.d0
  463. if(abs(sint).gt.1.e-06)SIGT=SIGN(1.D0,-SINT)
  464. TMIN=SIGT*(TMIN-XPI)+XPI
  465. IF(IMC.EQ.1) TM=2.D0*XPI
  466. TMAX=TMIN+TM
  467. CALL TEGRAL(SANGLE)
  468. IF(HCEN.EQ.0.D0) THEN
  469. C CAS DU TRONC DE CONE
  470. TETA=ATAN(PENT1)
  471. COUR=DENS*WIDTH/COS(TETA)
  472. IF(INDU) THEN
  473. CALL BTRCON(SANGLE,COUR,RINT,REXT,PENT1,TMIN,TMAX,
  474. & RP,ZP,BX,BY,BZ)
  475. ENDIF
  476. IF(POTE) THEN
  477. CALL ATRCON(SANGLE,COUR,RINT,REXT,PENT1,TMIN,TMAX,
  478. & RP,ZP,BX,BY,BZ)
  479. ENDIF
  480. ELSE
  481. IF(INDU) THEN
  482. CALL ARCTRA(SANGLE,COUR,HINT,HEXT,RINT,REXT,
  483. * PENT1,PENT2,TMIN,TMAX,
  484. * RP,ZP,BX,BY,BZ)
  485. ENDIF
  486. IF(POTE) THEN
  487. CALL AARCTR(SANGLE,COUR,HINT,HEXT,RINT,REXT,
  488. * PENT1,PENT2,TMIN,TMAX,
  489. * RP,ZP,BX,BY,BZ)
  490. ENDIF
  491. ENDIF
  492. RESP(1)=BX*COST-BY*SINT
  493. RESP(2)=BX*SINT+BY*COST
  494. RESP(3)=BZ
  495. ENDIF
  496. END IF
  497. IF(IERR.NE.0) THEN
  498. SEGDES MELEME
  499. SEGSUP MTRAV
  500. RETURN
  501. ENDIF
  502. C
  503. GOTO 700
  504. *
  505. 500 CONTINUE
  506. C CAS DES BARRES ON CHANGE DE RERERENTIEL
  507. DO 809 I= 1,3
  508. VY(I)=0.D0
  509. DO 8091 J= 1,3
  510. VY(I)= VY(I) + ROTA(I,J) * PN1(J)
  511. 8091 CONTINUE
  512. 809 CONTINUE
  513. * CALL ROTATI (PN1,1,1,1)
  514. C MESURE DANS LE REPERE LOCAL DES POINTS DE CALCUL DE HS
  515. IF(ISECT.EQ.3) THEN
  516. C CAS DES SECTIONS RETANGULAIRES
  517. IF(F2.EQ.0.D0) THEN
  518. C CAS DE LA PLAQUE RECTANGULAIRE
  519. DO 813 I=1,3
  520. VY(I)=0.D0
  521. DO 8131 J=1,3
  522. VY(I)=VY(I)+ROTA(I,J)*(PN1(J)-P1(J))
  523. 8131 CONTINUE
  524. 813 CONTINUE
  525. DY=F1
  526. REXT=0.5D0*DY
  527. RINT=-REXT
  528. COUR=DENS*DY
  529. YMIN=-VY(1)
  530. YMAX=YMIN+GY
  531. XP=-VY(2)
  532. ZP=VY(3)
  533. IF(INDU) THEN
  534. CALL BPLQDR(COUR,RINT,REXT,YMIN,YMAX,XP,ZP,BXP,BYP,BZP
  535. $ )
  536. ENDIF
  537. IF(POTE) THEN
  538. CALL APLQDR(COUR,RINT,REXT,YMIN,YMAX,XP,ZP,BXP,BYP,BZP
  539. $ )
  540. ENDIF
  541. RESP(1)=BYP
  542. RESP(2)=-BXP
  543. RESP(3)= BZP
  544. ELSE
  545. IF(INDU) THEN
  546. CALL CHABAR(IU,VZ,GY,F1,F2,DENS,VY(1),VY(2),VY(3),
  547. * BXP,BYP,BZP)
  548. RESP(1)=BXP
  549. RESP(2)=BYP
  550. RESP(3)=BZP
  551. ENDIF
  552. IF(POTE) THEN
  553. DO 814 I=1,3
  554. VY(I)=0.D0
  555. DO 8141 J=1,3
  556. VY(I)=VY(I)+ROTA(I,J)*(PN1(J)-P1(J))
  557. 8141 CONTINUE
  558. 814 CONTINUE
  559. DY=F1
  560. DZ=F2
  561. REXT=0.5D0*DY
  562. RINT=-REXT
  563. AIRE=DY*DZ
  564. COUR=DENS*AIRE
  565. PENT1=0.D0
  566. PENT2=0.D0
  567. YMIN=-VY(1)
  568. YMAX=YMIN+GY
  569. XP=-VY(2)
  570. ZP=VY(3)
  571. CALL ASEGTR(COUR,DZ,DZ,RINT,REXT,PENT1,PENT2,
  572. * YMIN,YMAX,XP,ZP,BXP,BYP,BZP)
  573. RESP(1)=BYP
  574. RESP(2)=-BXP
  575. RESP(3)= BZP
  576. ENDIF
  577. ENDIF
  578. ELSE IF(ISECT.EQ.1) THEN
  579. C CAS DES SECTIONS TRAPEZOIDALES
  580. DO 812 I=1,3
  581. VY(I)=0.D0
  582. DO 8121 J=1,3
  583. VY(I)=VY(I)+ROTA(I,J)*(PN1(J)-P1(J))
  584. 8121 CONTINUE
  585. 812 CONTINUE
  586. DY=F1
  587. DZ=F2
  588. REXT=0.5D0*DY
  589. RINT=-REXT
  590. AIRE=DY*DZ
  591. COUR=DENS*AIRE
  592. DELTA=0.5D0*DY*(PENT2-PENT1)
  593. HINT=DZ-DELTA
  594. HEXT=DZ+DELTA
  595. YMIN=-VY(1)
  596. YMAX=YMIN+GY
  597. XP=-VY(2)
  598. ZP=VY(3)
  599. IF(INDU) THEN
  600. CALL SEGTRA(COUR,HINT,HEXT,RINT,REXT,PENT1,PENT2,
  601. * YMIN,YMAX,XP,ZP,BXP,BYP,BZP)
  602. ENDIF
  603. IF(POTE) THEN
  604. CALL ASEGTR(COUR,HINT,HEXT,RINT,REXT,PENT1,PENT2,
  605. * YMIN,YMAX,XP,ZP,BXP,BYP,BZP)
  606. ENDIF
  607. RESP(1)=BYP
  608. RESP(2)=-BXP
  609. RESP(3)= BZP
  610. ENDIF
  611. GO TO 700
  612. C
  613. C CAS DU FIL
  614. 600 CONTINUE
  615. GY=0.D0
  616. DO 601 I=1,3
  617. VY(I)=P2(I)-P1(I)
  618. VI(I)=PN1(I)-P1(I)
  619. GY=GY+(VY(I)*VY(I))
  620. 601 CONTINUE
  621. GY=SQRT(GY)
  622. CALL PROVEC(VI,VY,VZ)
  623. CALL PROVEC(VY,VZ,VX)
  624. CALL ROTREP (ROTA,VX,VY,VZ)
  625. IF(IERR.NE.0) THEN
  626. SEGDES MELEME
  627. SEGSUP MTRAV
  628. RETURN
  629. ENDIF
  630. XP=0.D0
  631. YMIN=0.D0
  632. DO 602 J= 1,3
  633. XP=XP+ROTA(1,J)*VI(J)
  634. YMIN=YMIN-ROTA(2,J)*VI(J)
  635. 602 CONTINUE
  636. YMAX=YMIN+GY
  637. COUR=DENS
  638. IF(INDU) THEN
  639. CALL BFILDR(COUR,YMIN,YMAX,XP,BZP)
  640. RESP(1)=0.D0
  641. RESP(2)=0.D0
  642. RESP(3)=BZP
  643. ENDIF
  644. IF(POTE) THEN
  645. CALL AFILDR(COUR,YMIN,YMAX,XP,AYP)
  646. RESP(1)=0.D0
  647. RESP(2)=AYP
  648. RESP(3)=0.D0
  649. ENDIF
  650. CC
  651. 700 CONTINUE
  652. *
  653. C RETOUR DU CHAMP DANS LE REPERE GENERAL
  654.  
  655. DO 810 I= 1,3
  656. VY(I)=0.D0
  657. DO 8101 J= 1,3
  658. VY(I)= VY(I) + ROTA(J,I) * RESP(J) * RAP
  659. 8101 CONTINUE
  660. 810 CONTINUE
  661. * CALL RET_REF (RESP)
  662. * REMPLISSAGE DU SEGMENT DE TRAVAIL MTRAV
  663. *
  664. IGEO(IPT)=L
  665. DO 811 J= 1,3
  666. BB(J,IPT) = VY(J)
  667. IBIN(J,IPT) = 1
  668. 811 CONTINUE
  669. C FIN BOUCLE SUR LES POINTS
  670. *
  671. 30 CONTINUE
  672. SEGSUP SANGLE
  673. SEGDES,MELEME
  674. *
  675. CALL CRECHP(MTRAV,MCHPOI)
  676. SEGSUP MTRAV
  677. *
  678. * ATTRIBUTION D'UNE NATURE DISCRETE AU CHAMP QUI SORT
  679. *
  680.  
  681. 9990 SEGACT MCHPOI*MOD
  682. JATTRI(1) = 2
  683. SEGDES MCHPOI
  684. *
  685. *
  686. CALL ECROBJ('CHPOINT ',MCHPOI)
  687. RETURN
  688. END
  689.  
  690.  
  691.  
  692.  
  693.  
  694.  
  695.  
  696.  
  697.  
  698.  
  699.  
  700.  
  701.  
  702.  
  703.  
  704.  
  705.  
  706.  

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