Télécharger devlfa.eso

Retour à la liste

Numérotation des lignes :

  1. C DEVLFA SOURCE BP208322 15/07/22 21:15:23 8586
  2. SUBROUTINE DEVLFA(Q1,Q2,FTOTA,NA1,IPALA,IPLIA,XPALA,XVALA,
  3. & NLIAA,XDT,NPAS,IND,FINERT,IVINIT,FTEST,FTOTA0)
  4. IMPLICIT INTEGER(I-N)
  5. IMPLICIT REAL*8(A-H,O-Z)
  6. *--------------------------------------------------------------------*
  7. * *
  8. * Operateur DYNE : algorithme de Fu - de Vogelaere *
  9. * ________________________________________________ *
  10. * *
  11. * Calcul des forces de choc base A. *
  12. * *
  13. * Param}tres: *
  14. * *
  15. * e Q1(.,.) Vecteur des deplacements generalises. *
  16. * e Q2(.,.) Vecteur des vitesses generalisees. *
  17. * es FTOTA Forces ext{rieures totalis{es sur la base A. *
  18. * e NA1 Nombre total d'inconnues en base A. *
  19. * e IPALA Renseigne sur le type de la liaison. *
  20. * e IPLIA Tableau contenant les num{ros "DYNE" de la liaison. *
  21. * e XPALA Tableau contenant les param}tres de la liaison. *
  22. * es XVALA Tableau contenant les variables internes des liaisons *
  23. * e NLIAA Nombre de liaisons sur la base A. *
  24. * e XDT Valeurs des pas de temps *
  25. * e NPAS Numero du pas de temps *
  26. * e IND Indice du demi-pas de temps *
  27. * = 2 si 1er demi-pas = 1 si 2eme demi-pas *
  28. * es FINERT Forces d'inertie *
  29. * e IVINIT =1 si vitesses initiales, =0 sinon *
  30. * *
  31. * Auteur, date de cr{ation: *
  32. * *
  33. * Lionel VIVAN, le 20 aout 1989. *
  34. * E. de LANGRE, 15 fevrier 91 correction jeu negatif *
  35. * Lionel VIVAN, 12 mars 1991 ajout de la liaison FLUIDE *
  36. * Denis ROBERT, 30 avril 1992 ajout liaison POLYNOMIALE *
  37. * *
  38. * remarque: *
  39. * ========= *
  40. * Si jeu n{gatif (cas particulier de la base A ou il n'y a pas de *
  41. * normale), on renverse les variables avec XNORM. *
  42. * *
  43. *--------------------------------------------------------------------*
  44. *
  45. INTEGER IPALA(NLIAA,*),IPLIA(NLIAA,*)
  46. REAL*8 XPALA(NLIAA,*),Q1(NA1,*),Q2(NA1,*),FTOTA(NA1,*)
  47. REAL*8 XVALA(NLIAA,4,*),XDT(*),FINERT(NA1,*)
  48. PARAMETER (XZERO = 0.D0)
  49. REAL*8 FTest(nA1,4)
  50. REAL*8 ftOTa0(NA1,4)
  51. *
  52. * Boucle sur le nombre de liaisons
  53. *
  54. XFIN = 0.D0
  55. PDT = XDT(NPAS)
  56. PDTS2 = 0.5D0 * PDT
  57. DO 10 I = 1,NLIAA
  58. ITYP = IPALA(I,1)
  59. icond= IPALa(I,2)
  60. iannul= 0
  61. * --- cas des liaisons conditionnelles
  62. if (icond .eq. 1 ) then
  63. DO 101 j = 4,20
  64. jliai = ipala(i,j)
  65. jpliai = abs ( jliai)
  66. if ( jliai . EQ. 0 ) then
  67. goto 101
  68. else
  69. jtyp = ipala(jpliai,1)
  70. do 102 ik = 1,nA1
  71. do 103 jk = 1,4
  72. ftest(ik,jk) = 0d0
  73. ftota0 (ik,jk) = ftota(ik,jk)
  74. 103 continue
  75. 102 continue
  76. * CALCUL DES LIAISONS POUR TEST
  77. *
  78. * ****** choc {l{mentaire POINT_PLAN sans amortissement
  79. *
  80. IF (JTYP.EQ.1) THEN
  81. XRAID = XPALA(jpliai,1)
  82. XJEU = XPALA(jpliai,2)
  83. XNORM = 1.
  84. IF (XJEU.LT.0D0) THEN
  85. XNORM = -1.
  86. XJEU = -1*XJEU
  87. ENDIF
  88. INA1 = IPLIA(jpliai,1)
  89. XDEP = XNORM*Q1(INA1,IND)
  90. CALL DYCHEL(XDEP,XRAID,XJEU, XFL,iannul)
  91. XVALA(jpliai,IND,1) = XNORM*XFL
  92. XVALA(jpliai,IND,4) = XNORM*XDEP
  93. FTest(INA1,IND) =Ftest(INA1,IND) + XNORM*XFL
  94. *
  95. * ------ choc {l{mentaire POINT_PLAN avec amortissement
  96. *
  97. ELSE IF (JTYP.EQ.2) THEN
  98. XRAID = XPALA(jpliai,1)
  99. XJEU = XPALA(jpliai,2)
  100. XAMO = XPALA(jpliai,3)
  101. XNORM = 1.
  102. IPERM = 0
  103. IF (XJEU.LT.XZERO) THEN
  104. XNORM = -1.
  105. XJEU = -1*XJEU
  106. ENDIF
  107. INA1 = IPLIA(jpliai,1)
  108. XDEP = XNORM * Q1(INA1,IND)
  109. IF ((NPAS.EQ.1).AND.(IND.EQ.3)) THEN
  110. XVIT = XNORM * Q2(INA1,IND)
  111. ELSE
  112. IND2 = IND + 1
  113. XDEPM1 = XNORM * Q1(INA1,IND2)
  114. XVIT = (XDEP - XDEPM1) / PDTS2
  115. ENDIF
  116. XVALA(jpliai,IND,3) = XNORM*XVIT
  117. CALL DYCHAM(XDEP,XVIT,XRAID,XJEU,XAMO,XFL,IPERM,iannul)
  118. *+*
  119. IF (XDEP.GT.0.D0 .AND. XFL.GT.0.D0) XFL = 0.D0
  120. IF (XDEP.LT.0.D0 .AND. XFL.LT.0.D0) XFL = 0.D0
  121. *+*
  122. XVALA(jpliai,IND,1) = XNORM*XFL
  123. XVALA(jpliai,IND,4) = XNORM*XDEP
  124. FTest(INA1,IND) =Ftest(INA1,IND) + XNORM*XFL
  125. *
  126. * ------ choc {l{mentaire POINT_PLAN_FLUIDE
  127. *
  128. ELSE IF (JTYP.EQ.3) THEN
  129. XINER = XPALA(jpliai,1)
  130. XCONV = XPALA(jpliai,2)
  131. XVISC = XPALA(jpliai,3)
  132. XPCEL = XPALA(jpliai,4)
  133. XPCRA = XPALA(jpliai,5)
  134. XJEU = XPALA(jpliai,6)
  135. INA1 = IPLIA(jpliai,1)
  136. XDEP = Q1(INA1,IND)
  137. IF ((NPAS.EQ.1).AND.(IND.EQ.3)) THEN
  138. XVIT = Q2(INA1,IND)
  139. ELSE
  140. IND2 = IND + 1
  141. XDEPM1 = Q1(INA1,IND2)
  142. XVIT = (XDEP - XDEPM1) / PDTS2
  143. ENDIF
  144. IF (XJEU.GT.0D0) THEN
  145. XDH= XJEU - XDEP
  146. XNORM = 1.
  147. ELSE
  148. XDH= XDEP - XJEU
  149. XNORM = -1.
  150. ENDIF
  151. * Calcul de la masse ajout{e
  152. XXIN = -XINER / XDH
  153. FINERT(INA1,IND) = FINERT(INA1,IND) + XXIN
  154. * Calcul de la force de convection
  155. CALL DYFCON(XDH,XDEP,XVIT,XJEU,XCONV,XFCO,iannul)
  156. * Calcul de la force de viscosit{
  157. CALL DYFVIS(XDH,XDEP,XVIT,XJEU,XVISC,XFVI,iannul)
  158. * Calcul de la force de perte de charge
  159. CALL DYFPDC(XDH,XDEP,XVIT,XJEU,XPCEL,XPCRA,XFPE,iannul)
  160. XFL = (XFCO* XNORM ) + XFVI + XFPE
  161. XVALA(jpliai,IND,1) = XDEP
  162. XVALA(jpliai,IND,2) = XVIT
  163. XVALA(jpliai,IND,3) = XXIN
  164. XVALA(jpliai,IND,4) = XFCO*XNORM
  165. XVALA(jpliai,IND,5) = XFVI
  166. XVALA(jpliai,IND,6) = XFPE
  167. FTest(INA1,IND) =Ftest(iNA1,IND) + XFL
  168. *
  169. * ------ force {l{mentaire de COUPLAGE EN VITESSE
  170. *
  171. ELSE IF (JTYP.EQ.4) THEN
  172. INA1 = IPLIA(jpliai,1)
  173. INA2 = IPLIA(jpliai,2)
  174. XDEP = Q1(INA2,IND)
  175. IF ((NPAS.EQ.1).AND.(IND.EQ.3)) THEN
  176. XVIT = Q2(INA2,IND)
  177. ELSE
  178. IND2 = IND + 1
  179. XDEPM1 = Q1(INA2,IND2)
  180. XVIT = (XDEP - XDEPM1) / PDTS2
  181. ENDIF
  182. XCPLGE = XPALA(jpliai,1)
  183. CALL DYCPLV(XVIT,XCPLGE,XFL,iannul)
  184. XVALA(jpliai,IND,3) = XVIT
  185. XVALA(jpliai,IND,1) = XFL
  186. XVALA(jpliai,IND,4) = XDEP
  187. FTest(INA1,IND) =Ftest(INA1,IND) + XFL
  188. *
  189. * ------ force {l{mentaire de COUPLAGE EN DEPLACEMENT
  190. *
  191. ELSE IF (JTYP.EQ.5) THEN
  192. INA1 = IPLIA(jpliai,1)
  193. INA2 = IPLIA(jpliai,2)
  194. XDEP = Q1(INA2,IND)
  195. XCPLGE = XPALA(jpliai,1)
  196. CALL DYCPLD(XDEP,XCPLGE,XFL,iannul)
  197. cbp on mutliplie eventuellement par une fonction temporelle
  198. jfonct = ipala(jpliai,3)
  199. if(jfonct.ge.1) then
  200. XFREQ = XPALA(jpliai,2)
  201. XTIME = DBLE(NPAS-1)*PDT
  202. if(IND.EQ.2) XTIME=XTIME-PDTS2
  203. if(jfonct.eq.1) XFONCT = COS(XFREQ*XTIME)
  204. if(jfonct.eq.2) XFONCT = SIN(XFREQ*XTIME)
  205. XFL = XFL * XFONCT
  206. endif
  207. XVALA(jpliai,IND,1) = XFL
  208. XVALA(jpliai,IND,4) = XDEP
  209. FTest(INA1,IND) =Ftest(INA1,IND) + XFL
  210. *
  211. * ------ force elementaire de type POLYNOMIALE
  212. *
  213. ELSE IF (JTYP.EQ.6) THEN
  214. * WRITE(IOIMP,*)'DEVLFA : IPLIA(',jpliai,'1)=',IPLIA(jpliai,1)
  215. * nombre de modes "origine"
  216. NMOD = IPALA(I,2)
  217. CALL DYPOL1(Q1,Q2,NA1,IPLIA,XPALA,XVALA,NLIAA,IND,XDT,
  218. & NPAS,jpliai,NMOD,FTest,IVINIT,iannul)
  219. *
  220. * ------ choc {l{mentaire ...
  221. *
  222. * ELSE IF (JTYP.EQ. ) THEN
  223. * .......
  224. * .......
  225. *
  226. ENDIF
  227. ***** tests pour voir si on annulera
  228. xff = 0d0
  229. do 104 ik = 1,na1
  230. do 105 jk = 1,4
  231. xff = xff + ( ftest(ik,jk) ** 2)
  232. 105 continue
  233. 104 continue
  234. xff = xff ** .5
  235. if ( ((xff .le. 1e-20 ) .and. ( jliai .gt. 0) )
  236. & .OR. ((xff .gt. 1e-20 ) .and. ( jliai .lt. 0) ) )
  237. & then
  238. iannul = 1
  239. endif
  240. endif
  241. 101 continue
  242. * --- fin du cas des liaisons conditionnelles
  243. endif
  244. *
  245. * ------ choc {l{mentaire POINT_PLAN sans amortissement
  246. *
  247. IF (ITYP.EQ.1) THEN
  248. XRAID = XPALA(I,1)
  249. XJEU = XPALA(I,2)
  250. XNORM = 1.
  251. IF (XJEU.LT.0D0) THEN
  252. XNORM = -1.
  253. XJEU = -1*XJEU
  254. ENDIF
  255. INA1 = IPLIA(I,1)
  256. XDEP = XNORM*Q1(INA1,IND)
  257. CALL DYCHEL(XDEP,XRAID,XJEU, XFL,iannul)
  258. XVALA(I,IND,1) = XNORM*XFL
  259. XVALA(I,IND,4) = XNORM*XDEP
  260. FTOTA(INA1,IND) = FTOTA(INA1,IND) + XNORM*XFL
  261. *
  262. * ------ choc {l{mentaire POINT_PLAN avec amortissement
  263. *
  264. ELSE IF (ITYP.EQ.2) THEN
  265. XRAID = XPALA(I,1)
  266. XJEU = XPALA(I,2)
  267. XAMO = XPALA(I,3)
  268. XNORM = 1.
  269. IPERM = 0
  270. IF (XJEU.LT.XZERO) THEN
  271. XNORM = -1.
  272. XJEU = -1*XJEU
  273. ENDIF
  274. INA1 = IPLIA(I,1)
  275. XDEP = XNORM * Q1(INA1,IND)
  276. IF ((NPAS.EQ.1).AND.(IND.EQ.3)) THEN
  277. XVIT = XNORM * Q2(INA1,IND)
  278. ELSE
  279. IND2 = IND + 1
  280. XDEPM1 = XNORM * Q1(INA1,IND2)
  281. XVIT = (XDEP - XDEPM1) / PDTS2
  282. ENDIF
  283. XVALA(I,IND,3) = XNORM*XVIT
  284. CALL DYCHAM(XDEP,XVIT,XRAID,XJEU,XAMO,XFL,IPERM,iannul)
  285. *+*
  286. IF (XDEP.GT.0.D0 .AND. XFL.GT.0.D0) XFL = 0.D0
  287. IF (XDEP.LT.0.D0 .AND. XFL.LT.0.D0) XFL = 0.D0
  288. *+*
  289. XVALA(I,IND,1) = XNORM*XFL
  290. XVALA(I,IND,4) = XNORM*XDEP
  291. FTOTA(INA1,IND) = FTOTA(INA1,IND) + XNORM*XFL
  292. *
  293. * ------ choc {l{mentaire POINT_PLAN_FLUIDE
  294. *
  295. ELSE IF (ITYP.EQ.3) THEN
  296. XINER = XPALA(I,1)
  297. XCONV = XPALA(I,2)
  298. XVISC = XPALA(I,3)
  299. XPCEL = XPALA(I,4)
  300. XPCRA = XPALA(I,5)
  301. XJEU = XPALA(I,6)
  302. INA1 = IPLIA(I,1)
  303. XDEP = Q1(INA1,IND)
  304. IF ((NPAS.EQ.1).AND.(IND.EQ.3)) THEN
  305. XVIT = Q2(INA1,IND)
  306. ELSE
  307. IND2 = IND + 1
  308. XDEPM1 = Q1(INA1,IND2)
  309. XVIT = (XDEP - XDEPM1) / PDTS2
  310. ENDIF
  311. IF (XJEU.GT.0D0) THEN
  312. XDH= XJEU - XDEP
  313. XNORM = 1.
  314. ELSE
  315. XDH= XDEP - XJEU
  316. XNORM = -1.
  317. ENDIF
  318. * Calcul de la masse ajout{e
  319. XXIN = -XINER / XDH
  320. FINERT(INA1,IND) = FINERT(INA1,IND) + XXIN
  321. * Calcul de la force de convection
  322. CALL DYFCON(XDH,XDEP,XVIT,XJEU,XCONV,XFCO,iannul)
  323. * Calcul de la force de viscosit{
  324. CALL DYFVIS(XDH,XDEP,XVIT,XJEU,XVISC,XFVI,iannul)
  325. * Calcul de la force de perte de charge
  326. CALL DYFPDC(XDH,XDEP,XVIT,XJEU,XPCEL,XPCRA,XFPE,
  327. & iannul)
  328. XFL = (XFCO* XNORM ) + XFVI + XFPE
  329. XVALA(I,IND,1) = XDEP
  330. XVALA(I,IND,2) = XVIT
  331. XVALA(I,IND,3) = XXIN
  332. XVALA(I,IND,4) = XFCO*XNORM
  333. XVALA(I,IND,5) = XFVI
  334. XVALA(I,IND,6) = XFPE
  335. FTOTA(INA1,IND) = FTOTA(INA1,IND) + XFL
  336. *
  337. * ------ force {l{mentaire de COUPLAGE EN VITESSE
  338. *
  339. ELSE IF (ITYP.EQ.4) THEN
  340. INA1 = IPLIA(I,1)
  341. INA2 = IPLIA(I,2)
  342. XDEP = Q1(INA2,IND)
  343. IF ((NPAS.EQ.1).AND.(IND.EQ.3)) THEN
  344. XVIT = Q2(INA2,IND)
  345. ELSE
  346. IND2 = IND + 1
  347. XDEPM1 = Q1(INA2,IND2)
  348. XVIT = (XDEP - XDEPM1) / PDTS2
  349. ENDIF
  350. XCPLGE = XPALA(I,1)
  351. CALL DYCPLV(XVIT,XCPLGE,XFL,iannul)
  352. XVALA(I,IND,3) = XVIT
  353. XVALA(I,IND,1) = XFL
  354. XVALA(I,IND,4) = XDEP
  355. FTOTA(INA1,IND) = FTOTA(INA1,IND) + XFL
  356. *
  357. * ------ force {l{mentaire de COUPLAGE EN DEPLACEMENT
  358. *
  359. ELSE IF (ITYP.EQ.5) THEN
  360. INA1 = IPLIA(I,1)
  361. INA2 = IPLIA(I,2)
  362. XDEP = Q1(INA2,IND)
  363. XCPLGE = XPALA(I,1)
  364. jfonct = ipala(I,3)
  365. XFREQ = XPALA(I,2)
  366. CALL DYCPLD(XDEP,XCPLGE,XFL,iannul)
  367. cbp on mutliplie eventuellement par une fonction temporelle
  368. jfonct = ipala(I,3)
  369. if(jfonct.ge.1) then
  370. XFREQ = XPALA(I,2)
  371. XTIME = DBLE(NPAS-1)*PDT
  372. if(IND.EQ.2) XTIME=XTIME-PDTS2
  373. if(jfonct.eq.1) XFONCT = COS(XFREQ*XTIME)
  374. if(jfonct.eq.2) XFONCT = SIN(XFREQ*XTIME)
  375. XFL = XFL * XFONCT
  376. endif
  377. XVALA(I,IND,1) = XFL
  378. XVALA(I,IND,4) = XDEP
  379. FTOTA(INA1,IND) = FTOTA(INA1,IND) + XFL
  380. *
  381. * ------ force elementaire de type POLYNOMIALE
  382. *
  383. ELSE IF (ITYP.EQ.6) THEN
  384. * WRITE(IOIMP,*)'DEVLFA : IPLIA(',I,'1)=',IPLIA(I,1)
  385. * nombre de modes "origine"
  386. NMOD = IPALA(I,2)
  387. CALL DYPOL1(Q1,Q2,NA1,IPLIA,XPALA,XVALA,NLIAA,IND,XDT,
  388. & NPAS,I,NMOD,FTOTA,IVINIT,iannul)
  389. *
  390. * ------ choc {l{mentaire ...
  391. *
  392. * ELSE IF (ITYP.EQ. ) THEN
  393. * .......
  394. * .......
  395. *
  396. ENDIF
  397.  
  398.  
  399. *
  400. * la suite n'est plus utile car on passe iannul aux s_p de calcul des
  401. * forces de liaisons.
  402.  
  403. * si la liaisn était annulée on l'annule
  404. * if ( ( icond.eq. 1 ) .and. ( iannul.eq.1)) then
  405. * on annulle l'increment de ftotb
  406. * do 112 ik = 1,na1
  407. * do 113 jk = 1,4
  408. * ftota (ik,jk) = ftota0(ik,jk)
  409. * 113 continue
  410. * 112 continue
  411.  
  412. * end if
  413. 10 CONTINUE
  414. * END DO
  415. *
  416. END
  417.  
  418.  
  419.  
  420.  
  421.  
  422.  
  423.  

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