c**************************************************************** SUBROUTINE tdetau_piau(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 lois t(tau) de L.Piau 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 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 c nombre de points INTEGER, PARAMETER :: mx=4, nx=163 c épaisseur optique REAL (kind=dp), PARAMETER, DIMENSION(nx) :: taus=(/ 1 1.07504E-03, 2.51873E-03, 4.36329E-03, 6.61127E-03, 9.39312E-03, 2 1.10752E-02, 1.31208E-02, 1.56305E-02, 1.88616E-02, 2.32120E-02, 3 2.92046E-02, 3.74503E-02, 4.87520E-02, 6.53412E-02, 7.47222E-02, 4 8.57656E-02, 9.88327E-02, 1.13312E-01, 1.29631E-01, 1.48272E-01, 5 1.69670E-01, 1.93360E-01, 2.20212E-01, 2.51395E-01, 2.87442E-01, 6 3.29084E-01, 3.79390E-01, 4.45051E-01, 5.44726E-01, 7.43902E-01, 7 1.17371E+00, 2.03815E+00, 3.58251E+00, 6.09214E+00, 9.79602E+00, 8 1.49105E+01, 2.12075E+01, 2.86514E+01, 3.72406E+01, 4.71399E+01, 9 5.79066E+01, 6.97118E+01, 8.26743E+01, 9.70919E+01, 1.12355E+02, 1 1.28784E+02, 1.46533E+02, 1.66024E+02, 1.86445E+02, 2.08241E+02, 2 2.31588E+02, 2.57004E+02, 2.83450E+02, 3.11559E+02, 3.41519E+02, 3 3.73977E+02, 4.07617E+02, 4.43226E+02, 4.81031E+02, 5.21807E+02, 4 5.63938E+02, 6.08452E+02, 6.55547E+02, 7.06154E+02, 7.58282E+02, 5 8.13143E+02, 8.70981E+02, 9.32903E+02, 9.96563E+02, 1.06352E+03, 6 1.13397E+03, 1.20927E+03, 1.28661E+03, 1.36788E+03, 1.45325E+03, 7 1.54447E+03, 1.63820E+03, 1.73693E+03, 1.84062E+03, 1.95145E+03, 8 2.06544E+03, 2.18558E+03, 2.31174E+03, 2.44685E+03, 2.58613E+03, 9 2.73336E+03, 2.88804E+03, 3.05362E+03, 3.22458E+03, 3.40512E+03, 1 3.59485E+03, 3.79812E+03, 4.00853E+03, 4.23100E+03, 4.46485E+03, 2 4.71476E+03, 4.97380E+03, 5.24736E+03, 5.53515E+03, 5.84336E+03, 3 6.16376E+03, 6.50272E+03, 6.85939E+03, 7.24025E+03, 7.63656E+03, 4 8.05551E+03, 8.49676E+03, 8.96958E+03, 9.46305E+03, 9.98572E+03, 5 1.05364E+04, 1.11244E+04, 1.17383E+04, 1.23878E+04, 1.30727E+04, 6 1.38067E+04, 1.45750E+04, 1.53891E+04, 1.62477E+04, 1.71641E+04, 7 1.81241E+04, 1.91407E+04, 2.02143E+04, 2.13652E+04, 2.25741E+04, 8 2.38564E+04, 2.52113E+04, 2.66586E+04, 2.81801E+04, 2.97952E+04, 9 3.15050E+04, 3.33411E+04, 3.52772E+04, 3.73359E+04, 3.95165E+04, 1 4.18502E+04, 4.43128E+04, 4.69338E+04, 4.97158E+04, 5.27087E+04, 2 5.58752E+04, 5.92493E+04, 6.28312E+04, 6.66699E+04, 7.07329E+04, 3 7.50651E+04, 7.96722E+04, 8.46329E+04, 8.98946E+04, 9.55067E+04, 4 1.01475E+05, 1.07875E+05, 1.14665E+05, 1.21912E+05, 1.29633E+05, 5 1.37950E+05, 1.46793E+05, 1.56233E+05, 1.66290E+05, 1.77081E+05, 6 1.88554E+05, 2.00814E+05, 2.13893E+05 /) c températures REAL (kind=dp), PARAMETER, DIMENSION(nx) :: t1=(/ 1 4.39711E+03, 4.39877E+03, 4.40033E+03, 4.40175E+03, 4.40355E+03, 2 4.40388E+03, 4.40624E+03, 4.40936E+03, 4.41364E+03, 4.41926E+03, 3 4.42636E+03, 4.43569E+03, 4.44873E+03, 4.46962E+03, 4.49135E+03, 4 4.51937E+03, 4.55473E+03, 4.60249E+03, 4.66190E+03, 4.73254E+03, 5 4.81294E+03, 4.89912E+03, 4.99060E+03, 5.08770E+03, 5.18990E+03, 6 5.29742E+03, 5.41398E+03, 5.54679E+03, 5.70760E+03, 5.91993E+03, 7 6.19145E+03, 6.53522E+03, 6.94528E+03, 7.36536E+03, 7.76671E+03, 8 8.12402E+03, 8.40563E+03, 8.63935E+03, 8.83753E+03, 9.01048E+03, 9 9.15945E+03, 9.29763E+03, 9.42330E+03, 9.53897E+03, 9.64306E+03, 1 9.74152E+03, 9.83506E+03, 9.92409E+03, 1.00056E+04, 1.00830E+04, 2 1.01563E+04, 1.02275E+04, 1.02945E+04, 1.03596E+04, 1.04224E+04, 3 1.04843E+04, 1.05434E+04, 1.06013E+04, 1.06572E+04, 1.07123E+04, 4 1.07650E+04, 1.08167E+04, 1.08668E+04, 1.09163E+04, 1.09640E+04, 5 1.10109E+04, 1.10569E+04, 1.11027E+04, 1.11470E+04, 1.11911E+04, 6 1.12346E+04, 1.12783E+04, 1.13208E+04, 1.13631E+04, 1.14048E+04, 7 1.14465E+04, 1.14872E+04, 1.15276E+04, 1.15674E+04, 1.16073E+04, 8 1.16463E+04, 1.16855E+04, 1.17244E+04, 1.17637E+04, 1.18024E+04, 9 1.18413E+04, 1.18800E+04, 1.19191E+04, 1.19576E+04, 1.19964E+04, 1 1.20350E+04, 1.20741E+04, 1.21129E+04, 1.21521E+04, 1.21913E+04, 2 1.22312E+04, 1.22711E+04, 1.23115E+04, 1.23520E+04, 1.23930E+04, 3 1.24336E+04, 1.24744E+04, 1.25151E+04, 1.25561E+04, 1.25968E+04, 4 1.26380E+04, 1.26791E+04, 1.27209E+04, 1.27625E+04, 1.28042E+04, 5 1.28456E+04, 1.28872E+04, 1.29282E+04, 1.29695E+04, 1.30108E+04, 6 1.30528E+04, 1.30947E+04, 1.31370E+04, 1.31793E+04, 1.32220E+04, 7 1.32645E+04, 1.33074E+04, 1.33504E+04, 1.33942E+04, 1.34381E+04, 8 1.34826E+04, 1.35273E+04, 1.35729E+04, 1.36185E+04, 1.36647E+04, 9 1.37112E+04, 1.37583E+04, 1.38056E+04, 1.38535E+04, 1.39017E+04, 1 1.39508E+04, 1.40001E+04, 1.40501E+04, 1.41006E+04, 1.41518E+04, 2 1.42035E+04, 1.42557E+04, 1.43084E+04, 1.43616E+04, 1.44152E+04, 3 1.44693E+04, 1.45238E+04, 1.45792E+04, 1.46350E+04, 1.46915E+04, 4 1.47485E+04, 1.48061E+04, 1.48643E+04, 1.49232E+04, 1.49830E+04, 5 1.50439E+04, 1.51059E+04, 1.51685E+04, 1.52319E+04, 1.52958E+04, 6 1.53601E+04, 1.54245E+04, 1.54897E+04 /) REAL (kind=dp), SAVE, DIMENSION(nx+mx) :: taut REAL (kind=dp), SAVE, DIMENSION(1,nx) :: ts REAL (kind=dp), DIMENSION(mx,0:mx-1) :: fx REAL (kind=dp), PARAMETER :: teffc=5800.d0, ro_ext0=1.30917e-08 REAL (kind=dp), SAVE :: cte1 INTEGER, SAVE :: knot, l=1 LOGICAL, SAVE :: init=.TRUE. c--------------------------------------------------------------------- 2000 FORMAT(8es10.3) c----------initialisations début------------------ 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 tdetau_piau' ; STOP ENDIF rad=.FALSE. ; tau_min=3.d-3 WRITE(*,1)tau_min,ro_ext0 ; WRITE(2,1)tau_min,ro_ext0 1 FORMAT('loi t(tau,teff) 3D (tdetau_piau) non purement radiative',/, 1 'tau_min=',es10.3,', ro_ext=',es10.3,/, 2 'ATTENTION: seulement valable pour les modèles SOLAIRES') ENDIF c----------------initialisations fin ----------------------------- 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 tdetau_piau