Télécharger queshi.eso

Retour à la liste

Numérotation des lignes :

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

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