c**************************************************************** C SUBROUTINE k5777(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_tdetau c loi t(tau) de C van Veer issue du modèle Kurucz 5777, 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=143 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 /) REAL (kind=dp), PARAMETER, DIMENSION(nx) :: t1=(/ 1 3.7015D+03, 3.7213D+03, 3.7269D+03, 3.7401D+03, 3.7507D+03, 2 3.7633D+03, 3.7751D+03, 3.7877D+03, 3.8005D+03, 3.8138D+03, 3 3.8268D+03, 3.8397D+03, 3.8524D+03, 3.8649D+03, 3.8773D+03, 4 3.8894D+03, 3.9016D+03, 3.9135D+03, 3.9255D+03, 3.9371D+03, 5 3.9492D+03, 3.9608D+03, 3.9726D+03, 3.9842D+03, 3.9962D+03, 6 4.0078D+03, 4.0196D+03, 4.0311D+03, 4.0428D+03, 4.0541D+03, 7 4.0660D+03, 4.0772D+03, 4.0891D+03, 4.1006D+03, 4.1131D+03, 8 4.1246D+03, 4.1370D+03, 4.1490D+03, 4.1614D+03, 4.1739D+03, 9 4.1867D+03, 4.1995D+03, 4.2125D+03, 4.2256D+03, 4.2389D+03, 1 4.2523D+03, 4.2657D+03, 4.2793D+03, 4.2931D+03, 4.3071D+03, 2 4.3211D+03, 4.3351D+03, 4.3493D+03, 4.3632D+03, 4.3777D+03, 3 4.3918D+03, 4.4061D+03, 4.4204D+03, 4.4346D+03, 4.4488D+03, 4 4.4630D+03, 4.4773D+03, 4.4913D+03, 4.5059D+03, 4.5199D+03, 5 4.5344D+03, 4.5483D+03, 4.5626D+03, 4.5767D+03, 4.5911D+03, 6 4.6053D+03, 4.6198D+03, 4.6341D+03, 4.6486D+03, 4.6634D+03, 7 4.6781D+03, 4.6931D+03, 4.7082D+03, 4.7240D+03, 4.7400D+03, 8 4.7563D+03, 4.7732D+03, 4.7906D+03, 4.8086D+03, 4.8274D+03, 9 4.8469D+03, 4.8674D+03, 4.8889D+03, 4.9117D+03, 4.9358D+03, 1 4.9615D+03, 4.9891D+03, 5.0184D+03, 5.0515D+03, 5.0862D+03, 2 5.1236D+03, 5.1639D+03, 5.2075D+03, 5.2544D+03, 5.3055D+03, 3 5.3599D+03, 5.4204D+03, 5.4850D+03, 5.5560D+03, 5.6319D+03, 4 5.7136D+03, 5.8030D+03, 5.8996D+03, 6.0041D+03, 6.1207D+03, 5 6.2484D+03, 6.3873D+03, 6.5378D+03, 6.7024D+03, 6.8775D+03, 6 7.0526D+03, 7.2191D+03, 7.3727D+03, 7.5160D+03, 7.6498D+03, 7 7.7748D+03, 7.8927D+03, 8.0042D+03, 8.1106D+03, 8.2123D+03, 8 8.3105D+03, 8.4053D+03, 8.4974D+03, 8.5877D+03, 8.6754D+03, 9 8.7627D+03, 8.8476D+03, 8.9326D+03, 9.0158D+03, 9.0995D+03, 1 9.1814D+03, 9.2641D+03, 9.3462D+03, 9.4273D+03, 9.5110D+03, 2 9.5908D+03, 9.6756D+03, 9.7558D+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=5777.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 k5777' ; STOP ENDIF rad=.FALSE. ; tau_min=1.d-4 WRITE(*,1)tau_min,ro_ext0 ; WRITE(2,1)tau_min,ro_ext0 1 FORMAT(/,1x,'loi t(tau,teff,grav), k5777, non purement radiative',/, 1 'tau_min=',1pd10.3,' ro_ext=',1pd10.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 k5777