Télécharger queshi.eso

Retour à la liste

Numérotation des lignes :

  1. C QUESHI SOURCE BP208322 13/11/13 21:15:21 7865
  2. subroutine queshi(shifti,ienc,ifin,frequ,nmodma,shipro,iexis,
  3. $icle,nn,Northo)
  4. ccc $icle,nn,fpet,fgra)
  5. *
  6. implicit real*8(a-h,o-z)
  7. implicit integer(i-n)
  8. C-INC CCREEL
  9. -INC CCOPTIO
  10. *
  11. * routine pour calculer le shift a appliquer
  12. * appelée par strate
  13. *
  14. * creation : chat 11.2009
  15. * modifs : bp 01.2010 (traitements de qq cas particuliers)
  16. *
  17. * icle=1 on cherche à partir de ifin au dessus
  18. * icle=2 on cherche au dessous
  19. *
  20. *
  21. dimension frequ(*),iexis(*)
  22. *
  23. ideca= 4
  24. imul=1
  25. ifin2=nn
  26. if (icle.eq.2) then
  27. imul = -1
  28. ifin2=1
  29. endif
  30. *
  31. *---- on commence par regarder si on est pas en presence d'un trou ----*
  32. *
  33. itrou=0
  34. do iou=ifin+imul,ifin2,imul
  35. if (iexis(iou).ne.0) then
  36. cbp itrou=(iou-ifin)*imul
  37. itrou=((iou-ifin)*imul) - 1
  38. go to 2
  39. endif
  40. enddo
  41. 2 continue
  42. c write(6,*) 'iou,ifin,ifin2,itrou,ienc=',iou,ifin,ifin2,itrou,ienc
  43.  
  44.  
  45. *---- on est en presence d'un trou -----------------------------------*
  46. if (itrou.ne.0) then
  47.  
  48. cbp reporté + loin if (itrou.lt.ienc) then
  49.  
  50. 22 continue
  51. if (icle.eq.1) then
  52. fpet=frequ(ifin)
  53. fgra=frequ(iou)
  54. else
  55. fpet=frequ(iou)
  56. fgra=frequ(ifin)
  57. endif
  58. cbp : ajout pour le cas ou delta =0 car on a raté un mode et mal numeroté
  59. c et on cherche entre 2 frequences identiques trouvées avec 2 cycles
  60. df1 = abs(fpet - fgra) / (max(1.D-6,(abs(fpet)+abs(fgra))))
  61. if (df1.lt.1.D-4) then
  62. if (iou.ne.ifin2.and.iexis(iou+imul).ne.0) then
  63. iou=iou+imul
  64. goto 22
  65. elseif(iexis(ifin-imul).ne.0) then
  66. ifin=ifin-imul
  67. goto 22
  68. else
  69. write(6,*) 'ERREUR QUESHI : Contactez support '
  70. stop
  71. endif
  72. endif
  73.  
  74. if (itrou.lt.ienc) then
  75. fgra2=fgra*fgra
  76. fpet2=fpet*fpet
  77. c b2= 4.d0*fpet2
  78. c quatac=-2.d0*(fgra2-fpet2)
  79. c delta=b2 - quatac
  80. c if ( delta.lt.0.D0) then
  81. c call erreur(5)
  82. c return
  83. c endif
  84. c aa = sqrt(delta)
  85. c x2= (-2.d0 * fpet + aa )/ 2.d0
  86. c shipro=x2+fpet
  87. cbp : simplification
  88. fmoy2 = (fgra2 + fpet2) / 2.D0
  89. shipro = sqrt(fmoy2)
  90. cbp : cas d un petit trou, hyp. modes de multiplicite=2 => ienc*2
  91. c ienc=itrou + 1
  92. c ienc = 2 * (itrou + 1)
  93. c puisqu il y a dans strate: nmod=min ( ienc+2, nmodma), on peut reduire:
  94. ienc = 2 * itrou
  95. Northo = ienc - itrou
  96. c
  97. if (iimpi.ge.6) then
  98. WRITE(6,*) '----------------- Queshi -----------------'
  99. WRITE(IOIMP,1991) itrou,ienc,shipro
  100. 1991 FORMAT(2X,'PETIT TROU DE ',I5,' MODES'
  101. C ,'=> ON EN REDEMANDE ',I5,' AUTOUR DE ',F10.3)
  102. endif
  103. c
  104. return
  105.  
  106. else
  107. * on est en presence d'un trou et itrou>ienc donc on avait visé beaucoup
  108. * trop loin
  109. fpet=frequ(ifin)
  110. friou=frequ(iou)
  111. fgra= fpet + ( friou -fpet)/itrou*ienc
  112. Northo = 0
  113. go to 3
  114. endif
  115.  
  116. *---- pas de trou ----------------------------------------------------*
  117. else
  118. cbp : on ne passe ici que pour un cas improbable
  119. if ( icle.eq.2) then
  120. if (ifin-ienc.le.0) then
  121. WRITE(6,*) '!!! on ne passe ici que pour un cas improbable !!!'
  122. ienc=ifin - 1
  123. shipro=0.d0
  124. Northo = ienc
  125. return
  126. endif
  127. endif
  128.  
  129. endif
  130.  
  131. *---- calcul du shift a partir de la pente de y=freq(numero_mode) ----*
  132. * (cas sans trou)
  133. fpet= frequ(ifin)
  134. idec5 = ifin- (ideca*imul)
  135. idec10 = ifin- 2*(ideca*imul)
  136. idec20 = ifin- 4*(ideca*imul)
  137. pentm=0.D0
  138. pent5=0.D0
  139. pent10=0.D0
  140. pent20=0.D0
  141. cbp13/10/2010 : pentm = pente moyenne
  142. ndec = 0
  143. cbp plutot que 3,2,1 on prend 1,2,3 pour une moyenne + equilibrée
  144. icoe1 = 1
  145. icoe2 = 2
  146. icoe3 = 3
  147. if (idec5.gt.0) then
  148. if (iexis(idec5).ne.0) then
  149. pent5 = fpet - frequ(idec5)
  150. pentm = pentm + (icoe1*(pent5/ideca))
  151. ndec = ndec + icoe1
  152. endif
  153. endif
  154. if (idec10.gt.0) then
  155. if (iexis(idec10).ne.0) then
  156. pent10 = fpet - frequ(idec10)
  157. pentm = pentm + (icoe2*(pent10 /(2.*ideca)))
  158. ndec = ndec + icoe2
  159. endif
  160. endif
  161. if (idec20.gt.0) then
  162. if (iexis(idec20).ne.0) then
  163. pent20 = fpet - frequ(idec20)
  164. pentm = pentm + (icoe3*pent20 / (4.*ideca))
  165. ndec = ndec + icoe3
  166. endif
  167. endif
  168. cbp if (ndec.eq.0) then
  169. cbp write(6,*) 'ERREUR dans queshi.eso : contactez le support'
  170. cbp return
  171. cbp endif
  172. cbp fgra=fpet + (pentm * ienc)
  173. cbp 13/10/2010 :
  174. PENMIN = 3.E-4*fpet
  175. if ((ndec.ne.0.).and.(abs(pentm)).gt.PENMIN) then
  176. pentm = pentm/ndec
  177. else
  178. c pentm = imul * XPI * abs(fpet-shifti)
  179. pentm = imul * abs(fpet-shifti)
  180. write(6,*) 'DIFFICULTE dans QUESHI :',
  181. & ' calcul de la pente majorée ',pentm
  182. endif
  183. c on change completement la strategie en se placant juste apres fpet
  184. c fgra=fpet + (pentm * (min(ienc,2.)))
  185. c on revient a une strategie moins aggressive
  186. fgra=fpet + (pentm*ienc)
  187. if (iimpi.ge.6) then
  188. WRITE(6,*) 'pent5,pent10,pent20,pentm,fpet,fgra='
  189. WRITE(6,*) pent5,pent10,pent20,pentm,fpet,fgra,imul
  190. endif
  191. Northo = max(4,(ienc/10))
  192.  
  193. * calcul du shift (cas sans trou ou avec un grand trou)
  194. 3 continue
  195. if (icle.eq.2) then
  196. az=fpet
  197. fpet=fgra
  198. fgra=az
  199. endif
  200. c fgra2= fgra*fgra
  201. c fpet2= fpet*fpet
  202. c b2= 4.d0*fpet2
  203. c quatac= -2.d0*( fgra2-fpet2)
  204. c delta=b2-quatac
  205. c if (delta.lt.0.D0) then
  206. c call erreur(5)
  207. c return
  208. c endif
  209. c aa=sqrt(delta)
  210. c x2= (-2.d0 * fpet + aa)/2.D0
  211. c shipro=x2 + fpet
  212. cbp : simplification
  213. * fmoy2 = (fgra2 + fpet2) / 2.D0
  214. * shipro = sqrt(fmoy2)
  215. shipro = (fgra + fpet) / 2.D0
  216. c WRITE(6,*) 'fpet,fgra,shipro=',fpet,fgra,shipro
  217. ccc fgra = (shipro + fgra)/2.D0
  218. c
  219. if (iimpi.ge.6) then
  220. WRITE(6,*) '----------------- Queshi -----------------'
  221. WRITE(IOIMP,1993) ienc,shipro
  222. 1993 FORMAT(2X,'=> ON DEMANDE ',I5,' MODES AUTOUR DE ',F12.3)
  223. endif
  224. c
  225. return
  226. c
  227. end
  228.  
  229.  
  230.  
  231.  
  232.  
  233.  
  234.  
  235.  

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