Télécharger devfb4.eso

Retour à la liste

Numérotation des lignes :

  1. C DEVFB4 SOURCE CHAT 05/01/12 22:45:30 5004
  2. SUBROUTINE DEVFB4(ITYP,FTOTB,XPTB,IPALB,IPLIB,XPALB,XVALB,NLIAB,
  3. & NPLB,IND,IND1,INDM1,NPA,NPAM1,IND2,PDT,PDTS2,
  4. & FEXPSM,NPC1,I,IERRD,iannul)
  5. IMPLICIT INTEGER(I-N)
  6. IMPLICIT REAL*8(A-H,O-Z)
  7. *--------------------------------------------------------------------*
  8. * *
  9. * Op{rateur DYNE : algorithme de Fu - de Vogelaere *
  10. * ________________________________________________ *
  11. * *
  12. * Calcul des forces de choc sur base B pour les liaisons *
  13. * de type PROFIL_PROFIL_... *
  14. * *
  15. * Param}tres: *
  16. * *
  17. * e ITYP type de la liaison. *
  18. * es FTOTB Forces ext{rieures totalis{es sur la base B. *
  19. * e XPTB Tableau des d{placements des points *
  20. * e IPALB Renseigne sur la liaison. *
  21. * e IPLIB Tableau contenant les num{ros "DYNE" de la liaison. *
  22. * e XPALB Tableau contenant les param}tres de la liaison. *
  23. * es XVALB Tableau contenant les variables internes de liaisons. *
  24. * e NLIAB Nombre de liaisons sur la base B. *
  25. * e NPLB Nombre total de points intervenant dans les liaisons. *
  26. * e IND Indice du pas. *
  27. * e I num{ro de la liaison. *
  28. * es IERRD num{ro d'erreur. *
  29. * *
  30. * *
  31. * Auteur, date de cr{ation: *
  32. * *
  33. * Lionel VIVAN : le 1 f{vrier 1990 *
  34. * *
  35. *--------------------------------------------------------------------*
  36. *
  37. INTEGER IPALB(NLIAB,*),IPLIB(NLIAB,*)
  38. REAL*8 XPALB(NLIAB,*),XPTB(NPLB,4,*),FTOTB(NPLB,*)
  39. REAL*8 XVALB(NLIAB,4,*),XFOR(3),FEXPSM(NPLB,NPC1,2,*)
  40. LOGICAL LPETPP
  41. PARAMETER ( ZERO = 0.D0 , PRECIS = 1.D-15 )
  42. *
  43. * --- choc {l{mentaire PROFIL_PROFIL_INTERIEUR
  44. * --- choc {l{mentaire PROFIL_PROFIL_EXTERIEUR
  45. *
  46. IF (ITYP.EQ.31 .OR. ITYP.EQ.32) THEN
  47. NPOI = IPLIB(I,1)
  48. IDIM = IPALB(I,3)
  49. NOMBN1 = IPALB(I,4)
  50. NOMBN2 = IPALB(I,5)
  51. XRAID = XPALB(I,1)
  52. XSECT = XPALB(I,2)
  53. XEXPO = XPALB(I,3)
  54. ID1 = 3
  55. ID2 = ID1 + IDIM
  56. ID3 = ID1 + 2*IDIM
  57. ID4 = ID1 + 3*IDIM
  58. ID5 = ID1 + 4*IDIM
  59. ID6 = ID1 + 5*IDIM
  60. ID7 = ID6 + IDIM*NOMBN1
  61. ID8 = ID7 + IDIM*NOMBN2
  62. IP1 = 5
  63. IP2 = IP1 + NOMBN1
  64. IP3 = IP2
  65. IP4 = IP1 + NOMBN1 + NOMBN1*NOMBN2
  66. SOMFN = ZERO
  67. *
  68. * Recherche du nombre de point qui ont travers{ le maillage fixe,
  69. * donc s'il y a eu choc.
  70. IPOIN = 0
  71. ITIN = ID7
  72. DO 10 IN = 1,NOMBN2
  73. IP5 = IP4 + NOMBN1*(IN - 1)
  74. * D{placement du point IN du profil mobile dans le plan des profils
  75. XN1 = ZERO
  76. YN1 = ZERO
  77. DO 12 ID = 1,IDIM
  78. XXX = XPALB(I,ITIN+ID) + XPTB(NPOI,IND,ID)
  79. & - XPALB(I,ID2+ID)
  80. XN1 = XN1 + ( XXX * XPALB(I,ID3+ID) )
  81. YN1 = YN1 + ( XXX * XPALB(I,ID4+ID) )
  82. 12 CONTINUE
  83. * end do
  84. * On boucle sur le nombre d'{l{ment du profil fixe pour savoir si
  85. * le point IN @ travers{ le profil fixe
  86. NUME1 = ID6
  87. NUME2 = ID6 + IDIM
  88. ICOEF = ID8
  89. DO 14 IE = 1,NOMBN1
  90. IPALB(I,IP5+IE) = 0
  91. IPOSI = IPALB(I,IP3+IE)
  92. ICAS = IPALB(I,IP1+IE)
  93. XAIJ = XPALB(I,ICOEF+1)
  94. XBIJ = XPALB(I,ICOEF+2)
  95. IF (ICAS.EQ.1) THEN
  96. IF (XN1.GT.XAIJ) THEN
  97. IPOS2 = 1
  98. ELSE IF (XN1.LT.XAIJ) THEN
  99. IPOS2 = -1
  100. ELSE
  101. IPOS2 = 0
  102. ENDIF
  103. ELSE IF (ICAS.EQ.2) THEN
  104. IF (YN1.GT.XBIJ) THEN
  105. IPOS2 = 1
  106. ELSE IF (YN1.LT.XBIJ) THEN
  107. IPOS2 = -1
  108. ELSE
  109. IPOS2 = 0
  110. ENDIF
  111. ELSE
  112. RES = ( XAIJ * XN1 ) + XBIJ
  113. IF (YN1.GT.RES) THEN
  114. IPOS2 = 1
  115. ELSE IF (YN1.LT.RES) THEN
  116. IPOS2 = -1
  117. ELSE
  118. IPOS2 = 0
  119. ENDIF
  120. ENDIF
  121. IF (IPOSI.NE.IPOS2) THEN
  122. * Le point IN a travers{ le maillage fixe,
  123. * quel {l{ment a {t{ travers{ par ce point ?
  124. IF (IE.EQ.NOMBN1) NUME2 = ID6
  125. XE1 = ZERO
  126. YE1 = ZERO
  127. XE2 = ZERO
  128. YE2 = ZERO
  129. DO 16 ID = 1,IDIM
  130. XXX = XPALB(I,NUME1+ID) - XPALB(I,ID2+ID)
  131. YYY = XPALB(I,NUME2+ID) - XPALB(I,ID2+ID)
  132. XE1 = XE1 + ( XXX * XPALB(I,ID3+ID) )
  133. YE1 = YE1 + ( XXX * XPALB(I,ID4+ID) )
  134. XE2 = XE2 + ( YYY * XPALB(I,ID3+ID) )
  135. YE2 = YE2 + ( YYY * XPALB(I,ID4+ID) )
  136. 16 CONTINUE
  137. * end do
  138. IF (XE2.GT.XE1) THEN
  139. IF (XE1.LE.XN1 .AND. XN1.LT.XE2) THEN
  140. IPOIN = IPOIN + 1
  141. IPALB(I,IP5+IE) = 1
  142. GOTO 18
  143. ENDIF
  144. ELSE
  145. IF (XE2.LT.XN1 .AND. XN1.LE.XE1) THEN
  146. IPOIN = IPOIN + 1
  147. IPALB(I,IP5+IE) = 1
  148. GOTO 18
  149. ENDIF
  150. ENDIF
  151. IF (YE2.GT.YE1) THEN
  152. IF (YE1.LE.YN1 .AND. YN1.LT.YE2) THEN
  153. IPOIN = IPOIN + 1
  154. IPALB(I,IP5+IE) = 1
  155. GOTO 18
  156. ENDIF
  157. ELSE
  158. IF (YE2.LT.YN1 .AND. YN1.LE.YE1) THEN
  159. IPOIN = IPOIN + 1
  160. IPALB(I,IP5+IE) = 1
  161. GOTO 18
  162. ENDIF
  163. ENDIF
  164. ENDIF
  165. NUME1 = NUME2
  166. NUME2 = NUME2 + IDIM
  167. ICOEF = ICOEF + 2
  168. 14 CONTINUE
  169. * end do
  170. 18 CONTINUE
  171. IP3 = IP3 + NOMBN1
  172. ITIN = ITIN + IDIM
  173. 10 CONTINUE
  174. * end do
  175. *
  176. * Si le nombre de points ayant travers{s le profil fixe est {gale
  177. * au nombre de points du maillage mobile, alors on imprime une
  178. * erreur et on arrete le calcul.
  179. IF (IPOIN.EQ.NOMBN2) THEN
  180. IERRD = 1
  181. RETURN
  182. ENDIF
  183. *
  184. * s'il y a eu choc,
  185. * il y a IPOIN qui ont travers{ le maillage fixe
  186. *
  187. IF (IPOIN.NE.0) THEN
  188. *
  189. * Il y a eu choc :
  190. * - recherche du premier point qui a travers{ le maillage fixe, d'o:
  191. * d{termination du point P, intersection maillage fixe-mobile
  192. * - recherche du dernier point qui a travers{ le maillage fixe, d'o:
  193. * d{termination du point PP, intersection maillage fixe-mobile
  194. * - @ partir des points P et PP, d{termination de la normale de choc
  195. * - calcul de l'aire de choc, @ partir du point P, des points du
  196. * maillage fixe, du point PP et des points du maillage mobile qui
  197. * n'ont pas travers{ le maillage fixe
  198. * - calcul de la force de choc.
  199. *
  200. ICOMP = 0
  201. IN2 = 0
  202. 20 CONTINUE
  203. * On boucle sur les points du profil mobile
  204. IN2 = IN2 + 1
  205. IP5 = IP4 + NOMBN1*(IN2 - 1)
  206. IF (IN2.GT.NOMBN2) GOTO 22
  207. LPETPP = .FALSE.
  208. DO 24 IE1 = 1,NOMBN1
  209. IETAT = IPALB(I,IP5+IE1)
  210. * On regarde si le point IN2 a travers{ un {l{ment du profil fixe
  211. IF (IETAT.EQ.1) THEN
  212. ICOMP = ICOMP + 1
  213. IF (IN2.EQ.1) THEN
  214. * Si on traite le premier point, on cherche les points en partant du
  215. * dernier
  216. DO 26 INN2 = NOMBN2,1,-1
  217. IP6 = IP4 + NOMBN1*(INN2 - 1)
  218. DO 28 IEE1 = 1,NOMBN1
  219. IETA2 = IPALB(I,IP6+IEE1)
  220. IF (IETA2.EQ.1) THEN
  221. ICOMP = ICOMP + 1
  222. IEEE1 = IEE1
  223. GOTO 26
  224. ENDIF
  225. 28 CONTINUE
  226. * end do
  227. * Le point INN2 n'a pas travers{ le maillage fixe
  228. * Le point (INN2-1) et la point INN2 forment le segment de droite
  229. * pour d{finir le point d'intersection P avec l'{l{ment IEEE1
  230. IF (INN2.EQ.NOMBN2) THEN
  231. IPOIN1 = NOMBN2
  232. IPOIN2 = 1
  233. IELEM1 = IE1
  234. ELSE
  235. IPOIN1 = INN2
  236. IPOIN2 = INN2 + 1
  237. IELEM1 = IEEE1
  238. ENDIF
  239. GOTO 30
  240. 26 CONTINUE
  241. * end do
  242. 30 CONTINUE
  243. ELSE
  244. IPOIN1 = IN2 - 1
  245. IPOIN2 = IN2
  246. IELEM1 = IE1
  247. ENDIF
  248. * D{termination du point P
  249. * Le point IPOIN1 n'a pas travers{ le maillage fixe
  250. * Le point IPOIN2 a travers{ le maillage fixe par l'{l{ment IELEM1
  251. CALL DYNE31(IPOIN1,IPOIN2,IELEM1,XPALB,IPALB,XPTB,
  252. & NLIAB,NPLB,I,NPOI,IND,ID1,IP1,XP,YP)
  253. *
  254. * Recherche du point PP, on regarde les points suivant
  255. IF (IN2.EQ.NOMBN2) THEN
  256. IPOIN3 = IN2
  257. IPOIN4 = 1
  258. IELEM2 = IE1
  259. GOTO 32
  260. ENDIF
  261. IN22 = IN2 + 1
  262. IEEE1 = IE1
  263. DO 34 INN2 = IN22,NOMBN2
  264. IP7 = IP4 + NOMBN1*(INN2 - 1)
  265. DO 36 IEE1 = 1,NOMBN1
  266. IETA2 = IPALB(I,IP7+IEE1)
  267. IF (IETA2.EQ.1) THEN
  268. ICOMP = ICOMP + 1
  269. IEEE1 = IEE1
  270. IF (INN2.EQ.NOMBN2) THEN
  271. IPOIN3 = INN2
  272. IPOIN4 = 1
  273. IELEM2 = IEE1
  274. ENDIF
  275. GOTO 34
  276. ENDIF
  277. 36 CONTINUE
  278. * end do
  279. * Le point INN2 n'a pas travers{ le maillage fixe
  280. * Le point (INN2-1) et la point INN2 forment le segment de droite
  281. * pour d{finir le point d'intersection P avec l'{l{ment IEEE1
  282. IPOIN3 = INN2 - 1
  283. IPOIN4 = INN2
  284. IELEM2 = IEEE1
  285. GOTO 32
  286. 34 CONTINUE
  287. * end do
  288. 32 CONTINUE
  289. * D{termination du point PP
  290. * Le point IPOIN3 a travers{ le maillage fixe par l'{l{ment IELEM2
  291. * Le point IPOIN4 n'a pas travers{ le maillage fixe
  292. CALL DYNE31(IPOIN3,IPOIN4,IELEM2,XPALB,IPALB,XPTB,
  293. & NLIAB,NPLB,I,NPOI,IND,ID1,IP1,XPP,YPP)
  294. LPETPP = .TRUE.
  295. GOTO 40
  296. ENDIF
  297. 24 CONTINUE
  298. * end do
  299. 40 CONTINUE
  300. *
  301. IF ( LPETPP ) THEN
  302. * On a d{termin{ les points P et PP, calcul de la normale
  303. XNORM1 = ( XP + XPP ) * 0.5D0
  304. YNORM1 = ( YP + YPP ) * 0.5D0
  305. XXX = XP - XPP
  306. YYY = YP - YPP
  307. IPT3 = ID7 + IDIM*(IPOIN3-1)
  308. XN4 = ZERO
  309. YN4 = ZERO
  310. DO 42 ID = 1,IDIM
  311. XXX = XPALB(I,IPT3+ID) + XPTB(NPOI,IND,ID)
  312. & - XPALB(I,ID2+ID)
  313. XN4 = XN4 + ( XXX * XPALB(I,ID3+ID) )
  314. YN4 = YN4 + ( XXX * XPALB(I,ID4+ID) )
  315. 42 CONTINUE
  316. * END DO
  317. IF (ABS(YYY).LT.PRECIS) THEN
  318. YNN = YN4 - YNORM1
  319. XNORM = ZERO
  320. YNORM = YNN / ABS(YNN)
  321. ELSE IF (ABS(XXX).LT.PRECIS) THEN
  322. XNN = XN4 - XNORM1
  323. XNORM = XNN / ABS(XNN)
  324. YNORM = ZERO
  325. ELSE
  326. * le premier vecteur est P-PP
  327. * la normale de choc sera le produit vectoriel P_PP, N
  328. XPPP = XPP - XP
  329. YPPP = YPP - YP
  330. PS = SQRT( XPPP*XPPP + YPPP*YPPP )
  331. XNORM = YPPP * XPALB(I,ID5+3) / PS
  332. YNORM = -XPPP * XPALB(I,ID5+3) / PS
  333. ENDIF
  334. * Calcul de l'aire
  335. CALL DYNE32(IPOIN1,IPOIN4,IELEM1,IELEM2,XP,YP,XPP,YPP,
  336. & XPALB,IPALB,XPTB,NLIAB,NPLB,I,NPOI,ID1,IP1,IND,SOMAIR)
  337. AIRE = XSECT - SOMAIR
  338. * Calcul de la force de choc
  339. * Liaison conditionnelle
  340. IF (iannul.EQ.0) THEN
  341. FN = -XRAID * ( AIRE ** XEXPO )
  342. ELSE
  343. FN = 0
  344. ENDIF
  345. SOMFN = SOMFN + FN
  346. * On revient dans le repere global
  347. XFOR1 = FN * XNORM
  348. XFOR2 = FN * YNORM
  349. DO 38 ID = 1,IDIM
  350. XXX = XFOR1*XPALB(I,ID3+ID) + XFOR2*XPALB(I,ID4+ID)
  351. FTOTB(NPOI,ID) = FTOTB(NPOI,ID) + XXX
  352. 38 CONTINUE
  353. * end do
  354. IN2 = INN2
  355. ENDIF
  356. IF (ICOMP.EQ.IPOIN) GOTO 22
  357. GOTO 20
  358. 22 CONTINUE
  359. ENDIF
  360. DO 100 ID = 1,IDIM
  361. XVALB(I,IND,3+ID) = XPTB(NPOI,IND,ID)
  362. 100 CONTINUE
  363. * end do
  364. XVALB(I,IND,1) = SOMFN
  365. *
  366. * --- choc ...........
  367. *
  368. * ELSE IF (ITYP.EQ. ) THEN
  369. * .......
  370. * .......
  371. *
  372. ENDIF
  373. *
  374. END
  375.  
  376.  
  377.  

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