c**************************************************************** SUBROUTINE k5750(tau,teff,t,dtsdtau,dtsdteff,dtsdg, 1 ro_ext,dro_grav,dro_teff,f_tau,df_tau,d2f_tau) c routine private du module mod_atm c loi t(tau) de C van Veer issue du modèle Kurucz 5750, 4.47, X=.7, Z=.02 c mis sous la forme t=Teff(3/4(tau+q(tau)))**1/4, l/Hp=1.78 c Auteur: P. Morel, Département J.D. Cassini, O.C.A. c entrées : c tau : profondeur optique Rosseland c teff : température effective c grav : gravité c sorties : c t : température c dtsd* : dérivées t/ tau, teff, grav c dtsd* : dérivées t/ tau, teff, grav c dro_** : dérivées ro_ext/ teff, grav c f_tau, df_tau, df_tau2 : f, d f / d tau, d2 f / d2 tau c--------------------------------------------------------------------- USE mod_kind USE mod_numerique, ONLY : bsp1ddn, no_croiss IMPLICIT NONE REAL (kind=dp), INTENT(in) :: tau, teff REAL (kind=dp), INTENT(out) :: df_tau, dro_grav, dro_teff, dtsdg, 1 dtsdtau, dtsdteff, f_tau, d2f_tau, ro_ext, t INTEGER, PARAMETER :: mx=4, nx=144 REAL (kind=dp), PARAMETER, DIMENSION(nx) :: taus=(/ 1 1.3335D-07, 1.5399D-07, 1.7783D-07, 2.0535D-07, 2.3714D-07, 2 2.7384D-07, 3.1623D-07, 3.6517D-07, 4.2170D-07, 4.8697D-07, 3 5.6234D-07, 6.4938D-07, 7.4989D-07, 8.6596D-07, 1.0000D-06, 4 1.1548D-06, 1.3335D-06, 1.5399D-06, 1.7783D-06, 2.0535D-06, 5 2.3714D-06, 2.7384D-06, 3.1623D-06, 3.6517D-06, 4.2170D-06, 6 4.8697D-06, 5.6234D-06, 6.4938D-06, 7.4989D-06, 8.6596D-06, 7 1.0000D-05, 1.1548D-05, 1.3335D-05, 1.5399D-05, 1.7783D-05, 8 2.0535D-05, 2.3714D-05, 2.7384D-05, 3.1623D-05, 3.6517D-05, 9 4.2170D-05, 4.8697D-05, 5.6234D-05, 6.4938D-05, 7.4989D-05, 1 8.6596D-05, 1.0000D-04, 1.1548D-04, 1.3335D-04, 1.5399D-04, 2 1.7783D-04, 2.0535D-04, 2.3714D-04, 2.7384D-04, 3.1623D-04, 3 3.6517D-04, 4.2170D-04, 4.8697D-04, 5.6234D-04, 6.4938D-04, 4 7.4989D-04, 8.6596D-04, 1.0000D-03, 1.1548D-03, 1.3335D-03, 5 1.5399D-03, 1.7783D-03, 2.0535D-03, 2.3714D-03, 2.7384D-03, 6 3.1623D-03, 3.6517D-03, 4.2170D-03, 4.8697D-03, 5.6234D-03, 7 6.4938D-03, 7.4989D-03, 8.6596D-03, 1.0000D-02, 1.1548D-02, 8 1.3335D-02, 1.5399D-02, 1.7783D-02, 2.0535D-02, 2.3714D-02, 9 2.7384D-02, 3.1623D-02, 3.6517D-02, 4.2170D-02, 4.8697D-02, 1 5.6234D-02, 6.4938D-02, 7.4989D-02, 8.6596D-02, 1.0000D-01, 2 1.1548D-01, 1.3335D-01, 1.5399D-01, 1.7783D-01, 2.0535D-01, 3 2.3714D-01, 2.7384D-01, 3.1623D-01, 3.6517D-01, 4.2170D-01, 4 4.8697D-01, 5.6234D-01, 6.4938D-01, 7.4989D-01, 8.6596D-01, 5 1.0000D+00, 1.1548D+00, 1.3335D+00, 1.5399D+00, 1.7783D+00, 6 2.0535D+00, 2.3714D+00, 2.7384D+00, 3.1623D+00, 3.6517D+00, 7 4.2170D+00, 4.8697D+00, 5.6234D+00, 6.4938D+00, 7.4989D+00, 8 8.6596D+00, 1.0000D+01, 1.1548D+01, 1.3335D+01, 1.5399D+01, 9 1.7783D+01, 2.0535D+01, 2.3714D+01, 2.7384D+01, 3.1623D+01, 1 3.6517D+01, 4.2170D+01, 4.8697D+01, 5.6234D+01, 6.4938D+01, 2 7.4989D+01, 8.6596D+01, 1.0000D+02, 1.1548D+02 /) REAL (kind=dp), PARAMETER, DIMENSION(nx) :: t1=(/ 1 3.6798D+03, 3.6968D+03, 3.7035D+03, 3.7154D+03, 3.7262D+03, 2 3.7379D+03, 3.7494D+03, 3.7613D+03, 3.7731D+03, 3.7850D+03, 3 3.7966D+03, 3.8084D+03, 3.8202D+03, 3.8319D+03, 3.8435D+03, 4 3.8549D+03, 3.8662D+03, 3.8775D+03, 3.8886D+03, 3.8999D+03, 5 3.9110D+03, 3.9220D+03, 3.9330D+03, 3.9440D+03, 3.9551D+03, 6 3.9660D+03, 3.9772D+03, 3.9885D+03, 4.0001D+03, 4.0117D+03, 7 4.0235D+03, 4.0355D+03, 4.0475D+03, 4.0600D+03, 4.0724D+03, 8 4.0851D+03, 4.0978D+03, 4.1106D+03, 4.1236D+03, 4.1366D+03, 9 4.1501D+03, 4.1630D+03, 4.1769D+03, 4.1906D+03, 4.2045D+03, 1 4.2185D+03, 4.2325D+03, 4.2467D+03, 4.2611D+03, 4.2753D+03, 2 4.2896D+03, 4.3038D+03, 4.3182D+03, 4.3324D+03, 4.3469D+03, 3 4.3611D+03, 4.3759D+03, 4.3903D+03, 4.4050D+03, 4.4196D+03, 4 4.4343D+03, 4.4489D+03, 4.4639D+03, 4.4789D+03, 4.4938D+03, 5 4.5086D+03, 4.5236D+03, 4.5384D+03, 4.5533D+03, 4.5682D+03, 6 4.5834D+03, 4.5984D+03, 4.6135D+03, 4.6287D+03, 4.6441D+03, 7 4.6596D+03, 4.6751D+03, 4.6912D+03, 4.7075D+03, 4.7241D+03, 8 4.7410D+03, 4.7584D+03, 4.7764D+03, 4.7949D+03, 4.8142D+03, 9 4.8342D+03, 4.8552D+03, 4.8771D+03, 4.9003D+03, 4.9248D+03, 1 4.9508D+03, 4.9785D+03, 5.0085D+03, 5.0405D+03, 5.0749D+03, 2 5.1121D+03, 5.1521D+03, 5.1954D+03, 5.2421D+03, 5.2928D+03, 3 5.3472D+03, 5.4069D+03, 5.4715D+03, 5.5416D+03, 5.6175D+03, 4 5.6988D+03, 5.7869D+03, 5.8832D+03, 5.9883D+03, 6.1043D+03, 5 6.2344D+03, 6.3783D+03, 6.5423D+03, 6.7128D+03, 6.8728D+03, 6 7.0225D+03, 7.1645D+03, 7.2982D+03, 7.4256D+03, 7.5472D+03, 7 7.6638D+03, 7.7759D+03, 7.8840D+03, 7.9887D+03, 8.0901D+03, 8 8.1891D+03, 8.2856D+03, 8.3801D+03, 8.4728D+03, 8.5640D+03, 9 8.6538D+03, 8.7426D+03, 8.8305D+03, 8.9175D+03, 9.0038D+03, 1 9.0897D+03, 9.1749D+03, 9.2602D+03, 9.3446D+03, 9.4297D+03, 2 9.5139D+03, 9.5987D+03, 9.6841D+03, 9.7674D+03 /) REAL (kind=dp), SAVE, DIMENSION(nx+2*mx) :: taut REAL (kind=dp), SAVE, DIMENSION(1,nx) :: ts REAL (kind=dp), DIMENSION(mx,0:mx-1) :: fx REAL (kind=dp), PARAMETER :: teffc=5750.d0, ro_ext0=3.55d-9 REAL (kind=dp), SAVE :: cte1 INTEGER, SAVE :: knot, l=1 LOGICAL, SAVE :: init=.TRUE. c--------------------------------------------------------------------- IF(init)THEN init=.FALSE. ; cte1=(3.d0/4.d0)**0.25d0 ; ts=RESHAPE(t1,(/1,nx/)) ts=(ts/teffc)**4*4.d0/3.d0 CALL bsp1ddn(1,ts,taus,taut,nx,mx,knot,.FALSE.,taus(1),l,fx) IF(no_croiss)THEN PRINT*,'Arrêt 1 dans k5750' ; STOP ENDIF rad=.FALSE. ; tau_min=1.d-4 WRITE(*,1)tau_min,ro_ext0 ; WRITE(2,1)tau_min,ro_ext0 1 FORMAT(/,'loi t(tau,teff,grav), k5750, non purement radiative',/, 1 'tau_min=',es10.3,' ro_ext=',es10.3,/) ENDIF CALL bsp1ddn(1,ts,taus,taut,nx,mx,knot,.TRUE.,tau,l,fx) f_tau=fx(1,0) ; df_tau=fx(1,1) ;d2f_tau=fx(1,2) t=teff*cte1*f_tau**0.25d0 ; dtsdtau=t/4.d0/f_tau*df_tau dtsdteff=t/teff dtsdg=0.d0 ; ro_ext=ro_ext0 ; dro_grav=0.d0 ; dro_teff=0.d0 RETURN END SUBROUTINE k5750