Index: trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90 =================================================================== --- trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90 (revision 1714) +++ trunk/LMDZ.GENERIC/libf/phystd/callcorrk.F90 (revision 1715) @@ -123,5 +123,5 @@ REAL*8 taucumi(L_LEVELS,L_NSPECTI,L_NGAUSS) - REAL*8 tauaero(L_LEVELS+1,naerkind) + REAL*8 tauaero(L_LEVELS,naerkind) REAL*8 nfluxtopv,nfluxtopi,nfluxtop,fluxtopvdn REAL*8 nfluxoutv_nu(L_NSPECTV) ! Outgoing band-resolved VI flux at TOA (W/m2). @@ -189,10 +189,10 @@ ! test on allocated necessary because of CLFvarying (two calls to callcorrk in physiq) - if(.not.allocated(QXVAER)) allocate(QXVAER(L_LEVELS+1,L_NSPECTV,naerkind)) - if(.not.allocated(QSVAER)) allocate(QSVAER(L_LEVELS+1,L_NSPECTV,naerkind)) - if(.not.allocated(GVAER)) allocate(GVAER(L_LEVELS+1,L_NSPECTV,naerkind)) - if(.not.allocated(QXIAER)) allocate(QXIAER(L_LEVELS+1,L_NSPECTI,naerkind)) - if(.not.allocated(QSIAER)) allocate(QSIAER(L_LEVELS+1,L_NSPECTI,naerkind)) - if(.not.allocated(GIAER)) allocate(GIAER(L_LEVELS+1,L_NSPECTI,naerkind)) + if(.not.allocated(QXVAER)) allocate(QXVAER(L_LEVELS,L_NSPECTV,naerkind)) + if(.not.allocated(QSVAER)) allocate(QSVAER(L_LEVELS,L_NSPECTV,naerkind)) + if(.not.allocated(GVAER)) allocate(GVAER(L_LEVELS,L_NSPECTV,naerkind)) + if(.not.allocated(QXIAER)) allocate(QXIAER(L_LEVELS,L_NSPECTI,naerkind)) + if(.not.allocated(QSIAER)) allocate(QSIAER(L_LEVELS,L_NSPECTI,naerkind)) + if(.not.allocated(GIAER)) allocate(GIAER(L_LEVELS,L_NSPECTI,naerkind)) !!! ALLOCATED instances are necessary because of CLFvarying (strategy to call callcorrk twice in physiq...) @@ -323,4 +323,5 @@ end do !iaer=1,naerkind. + ! How much light do we get ? do nw=1,L_NSPECTV @@ -439,5 +440,5 @@ ! Test / Correct for freaky s. s. albedo values. do iaer=1,naerkind - do k=1,L_LEVELS+1 + do k=1,L_LEVELS do nw=1,L_NSPECTV @@ -482,7 +483,5 @@ ! boundary conditions tauaero(1,iaer) = tauaero(2,iaer) - tauaero(L_LEVELS+1,iaer) = tauaero(L_LEVELS,iaer) !tauaero(1,iaer) = 0. - !tauaero(L_LEVELS+1,iaer) = 0. end do ! naerkind @@ -574,11 +573,5 @@ if(kastprof)then - - if(.not.global1d)then ! garde-fou/safeguard added by MT (to be removed in the future) - write(*,*) 'You have to fix mu0, ' - write(*,*) 'the cosinus of the solar angle' - stop - endif - + ! Initial values equivalent to mugaz. DO l=1,nlayer Index: trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90 =================================================================== --- trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90 (revision 1714) +++ trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90 (revision 1715) @@ -11,5 +11,4 @@ use init_print_control_mod, only: init_print_control use radinc_h, only: ini_radinc_h, naerkind - use radcommon_h, only: ini_radcommon_h use datafile_mod, only: datadir use comdiurn_h, only: sinlat, coslat, sinlon, coslon @@ -762,7 +761,4 @@ call ini_radinc_h(nlayer) - ! allocate "radcommon_h" arrays - call ini_radcommon_h() - ! allocate "comsoil_h" arrays call ini_comsoil_h(ngrid) Index: trunk/LMDZ.GENERIC/libf/phystd/optci.F90 =================================================================== --- trunk/LMDZ.GENERIC/libf/phystd/optci.F90 (revision 1714) +++ trunk/LMDZ.GENERIC/libf/phystd/optci.F90 (revision 1715) @@ -32,5 +32,5 @@ real*8 DTAUI(L_NLAYRAD,L_NSPECTI,L_NGAUSS) - real*8 DTAUKI(L_LEVELS+1,L_NSPECTI,L_NGAUSS) + real*8 DTAUKI(L_LEVELS,L_NSPECTI,L_NGAUSS) real*8 TAUI(L_NLEVRAD,L_NSPECTI,L_NGAUSS) real*8 TAUCUMI(L_LEVELS,L_NSPECTI,L_NGAUSS) @@ -42,9 +42,9 @@ ! for aerosols - real*8 QXIAER(L_LEVELS+1,L_NSPECTI,NAERKIND) - real*8 QSIAER(L_LEVELS+1,L_NSPECTI,NAERKIND) - real*8 GIAER(L_LEVELS+1,L_NSPECTI,NAERKIND) - real*8 TAUAERO(L_LEVELS+1,NAERKIND) - real*8 TAUAEROLK(L_LEVELS+1,L_NSPECTI,NAERKIND) + real*8 QXIAER(L_LEVELS,L_NSPECTI,NAERKIND) + real*8 QSIAER(L_LEVELS,L_NSPECTI,NAERKIND) + real*8 GIAER(L_LEVELS,L_NSPECTI,NAERKIND) + real*8 TAUAERO(L_LEVELS,NAERKIND) + real*8 TAUAEROLK(L_LEVELS,L_NSPECTI,NAERKIND) real*8 TAEROS(L_LEVELS,L_NSPECTI,NAERKIND) @@ -135,5 +135,5 @@ end do ! levels - + ! Spectral dependance of aerosol absorption do iaer=1,naerkind DO NW=1,L_NSPECTI @@ -147,7 +147,8 @@ do K=2,L_LEVELS - -! continuum absorption - DCONT = 0.0d0 + + DAERO=SUM(TAEROS(K,NW,1:naerkind)) ! aerosol absorption + + DCONT = 0.0d0 ! continuum absorption if(continuum.and.(.not.graybody))then @@ -225,7 +226,4 @@ endif -! aerosol absorption - DAERO=SUM(TAEROS(K,NW,1:naerkind)) - do ng=1,L_NGAUSS-1 @@ -286,5 +284,4 @@ end do - DTAUKI(L_LEVELS+1,1:L_NSPECTI,1:L_NGAUSS)=0.d0 !======================================================================= ! Now the full treatment for the layers, where besides the opacity @@ -293,6 +290,6 @@ do iaer=1,naerkind DO NW=1,L_NSPECTI - DO K=2,L_LEVELS+1 - TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER)*QSIAER(K,NW,IAER) + DO K=2,L_LEVELS + TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER)*QSIAER(K,NW,IAER) ! effect of scattering albedo ENDDO ENDDO @@ -300,8 +297,14 @@ DO NW=1,L_NSPECTI - DO L=1,L_NLAYRAD + DO L=1,L_NLAYRAD-1 K = 2*L+1 btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind)) END DO ! L vertical loop + + ! Last level + L = L_NLAYRAD + K = 2*L+1 + btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + END DO ! NW spectral loop @@ -309,5 +312,5 @@ DO NW=1,L_NSPECTI NG = L_NGAUSS - DO L=1,L_NLAYRAD + DO L=1,L_NLAYRAD-1 K = 2*L+1 @@ -334,4 +337,28 @@ END DO ! L vertical loop + + ! Last level + + L = L_NLAYRAD + K = 2*L+1 + DTAUI(L,nw,ng) = DTAUKI(K,NW,NG) ! + 1.e-50 + + atemp = 0. + if(DTAUI(L,NW,NG) .GT. 1.0D-9) then + do iaer=1,naerkind + atemp = atemp + GIAER(K,NW,IAER) * TAUAEROLK(K,NW,IAER) + end do + WBARI(L,nw,ng) = btemp(L,nw) / DTAUI(L,NW,NG) + else + WBARI(L,nw,ng) = 0.0D0 + DTAUI(L,NW,NG) = 1.0D-9 + endif + + if(btemp(L,nw) .GT. 0.0d0) then + cosbi(L,NW,NG) = atemp/btemp(L,nw) + else + cosbi(L,NW,NG) = 0.0D0 + end if + ! Now the other Gauss points, if needed. @@ -340,5 +367,5 @@ IF(TAUGSURF(NW,NG) .gt. TLIMIT) THEN - DO L=1,L_NLAYRAD + DO L=1,L_NLAYRAD-1 K = 2*L+1 DTAUI(L,nw,ng) = DTAUKI(K,NW,NG)+DTAUKI(K+1,NW,NG)! + 1.e-50 @@ -355,4 +382,21 @@ cosbi(L,NW,NG) = cosbi(L,NW,L_NGAUSS) END DO ! L vertical loop + + ! Last level + L = L_NLAYRAD + K = 2*L+1 + DTAUI(L,nw,ng) = DTAUKI(K,NW,NG)! + 1.e-50 + + if(DTAUI(L,NW,NG) .GT. 1.0D-9) then + + WBARI(L,nw,ng) = btemp(L,nw) / DTAUI(L,NW,NG) + + else + WBARI(L,nw,ng) = 0.0D0 + DTAUI(L,NW,NG) = 1.0D-9 + endif + + cosbi(L,NW,NG) = cosbi(L,NW,L_NGAUSS) + END IF Index: trunk/LMDZ.GENERIC/libf/phystd/optcv.F90 =================================================================== --- trunk/LMDZ.GENERIC/libf/phystd/optcv.F90 (revision 1714) +++ trunk/LMDZ.GENERIC/libf/phystd/optcv.F90 (revision 1715) @@ -24,5 +24,5 @@ ! ! THIS SUBROUTINE SETS THE OPTICAL CONSTANTS IN THE VISUAL - ! IT CALCUALTES FOR EACH LAYER, FOR EACH SPECRAL INTERVAL IN THE VISUAL + ! IT CALCULATES FOR EACH LAYER, FOR EACH SPECTRAL INTERVAL IN THE VISUAL ! LAYER: WBAR, DTAU, COSBAR ! LEVEL: TAU @@ -39,5 +39,5 @@ real*8 DTAUV(L_NLAYRAD,L_NSPECTV,L_NGAUSS) - real*8 DTAUKV(L_LEVELS+1,L_NSPECTV,L_NGAUSS) + real*8 DTAUKV(L_LEVELS,L_NSPECTV,L_NGAUSS) real*8 TAUV(L_NLEVRAD,L_NSPECTV,L_NGAUSS) real*8 TAUCUMV(L_LEVELS,L_NSPECTV,L_NGAUSS) @@ -48,9 +48,9 @@ ! for aerosols - real*8 QXVAER(L_LEVELS+1,L_NSPECTV,NAERKIND) - real*8 QSVAER(L_LEVELS+1,L_NSPECTV,NAERKIND) - real*8 GVAER(L_LEVELS+1,L_NSPECTV,NAERKIND) - real*8 TAUAERO(L_LEVELS+1,NAERKIND) - real*8 TAUAEROLK(L_LEVELS+1,L_NSPECTV,NAERKIND) + real*8 QXVAER(L_LEVELS,L_NSPECTV,NAERKIND) + real*8 QSVAER(L_LEVELS,L_NSPECTV,NAERKIND) + real*8 GVAER(L_LEVELS,L_NSPECTV,NAERKIND) + real*8 TAUAERO(L_LEVELS,NAERKIND) + real*8 TAUAEROLK(L_LEVELS,L_NSPECTV,NAERKIND) real*8 TAEROS(L_LEVELS,L_NSPECTV,NAERKIND) @@ -60,5 +60,4 @@ real*8 TAURAY(L_NSPECTV) real*8 TRAY(L_LEVELS,L_NSPECTV) - real*8 TRAYAER real*8 DPR(L_LEVELS), U(L_LEVELS) real*8 LCOEF(4), LKCOEF(L_LEVELS,4) @@ -66,4 +65,5 @@ real*8 taugsurf(L_NSPECTV,L_NGAUSS-1) real*8 DCONT,DAERO + real*8 DRAYAER double precision wn_cont, p_cont, p_air, T_cont, dtemp, dtempc double precision p_cross @@ -123,5 +123,5 @@ end do ! levels - + ! Spectral dependance of aerosol absorption do iaer=1,naerkind do NW=1,L_NSPECTV @@ -131,4 +131,6 @@ end do end do + + ! Rayleigh scattering do NW=1,L_NSPECTV do K=2,L_LEVELS @@ -142,8 +144,8 @@ do NW=1,L_NSPECTV - TRAYAER = TRAY(K,NW) - ! TRAYAER is Tau RAYleigh scattering, plus AERosol opacity + DRAYAER = TRAY(K,NW) + ! DRAYAER is Tau RAYleigh scattering, plus AERosol opacity do iaer=1,naerkind - TRAYAER = TRAYAER + TAEROS(K,NW,IAER) + DRAYAER = DRAYAER + TAEROS(K,NW,IAER) end do @@ -257,5 +259,5 @@ TAUGSURF(NW,NG) = TAUGSURF(NW,NG) + TAUGAS + DCONT DTAUKV(K,nw,ng) = TAUGAS & - + TRAYAER & ! TRAYAER includes all scattering contributions + + DRAYAER & ! DRAYAER includes all scattering contributions + DCONT ! For parameterized continuum aborption @@ -266,9 +268,9 @@ NG = L_NGAUSS - DTAUKV(K,nw,ng) = TRAY(K,NW) + DCONT ! For parameterized continuum absorption - - do iaer=1,naerkind - DTAUKV(K,nw,ng) = DTAUKV(K,nw,ng) + TAEROS(K,NW,IAER) - end do ! a bug was here! + DTAUKV(K,nw,ng) = DRAYAER + DCONT ! Scattering + parameterized continuum absorption + +! do iaer=1,naerkind +! DTAUKV(K,nw,ng) = DTAUKV(K,nw,ng) + TAEROS(K,NW,IAER) +! end do ! a bug was here! end do @@ -282,6 +284,6 @@ do iaer=1,naerkind DO NW=1,L_NSPECTV - DO K=2,L_LEVELS ! AS: shouldn't this be L_LEVELS+1 ? (see optci) - TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER) * QSVAER(K,NW,IAER) + DO K=2,L_LEVELS + TAUAEROLK(K,NW,IAER) = TAUAERO(K,IAER) * QSVAER(K,NW,IAER) ! effect of scattering albedo ENDDO ENDDO @@ -293,15 +295,15 @@ atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind))+SUM(GVAER(K+1,NW,1:naerkind) * TAUAEROLK(K+1,NW,1:naerkind)) btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) + SUM(TAUAEROLK(K+1,NW,1:naerkind)) - ctemp(L,NW) = btemp(L,NW) + 0.9999*(TRAY(K,NW) + TRAY(K+1,NW)) + ctemp(L,NW) = btemp(L,NW) + 0.9999*(TRAY(K,NW) + TRAY(K+1,NW)) ! JVO 2017 : does this 0.999 is really meaningful ? btemp(L,NW) = btemp(L,NW) + TRAY(K,NW) + TRAY(K+1,NW) COSBV(L,NW,1:L_NGAUSS) = atemp(L,NW)/btemp(L,NW) END DO ! L vertical loop - !last level - L = L_NLAYRAD - K = 2*L+1 - atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind)) + ! Last level + L = L_NLAYRAD + K = 2*L+1 + atemp(L,NW) = SUM(GVAER(K,NW,1:naerkind) * TAUAEROLK(K,NW,1:naerkind)) btemp(L,NW) = SUM(TAUAEROLK(K,NW,1:naerkind)) - ctemp(L,NW) = btemp(L,NW) + 0.9999*TRAY(K,NW) + ctemp(L,NW) = btemp(L,NW) + 0.9999*TRAY(K,NW) ! JVO 2017 : does this 0.999 is really meaningful ? btemp(L,NW) = btemp(L,NW) + TRAY(K,NW) COSBV(L,NW,1:L_NGAUSS) = atemp(L,NW)/btemp(L,NW) @@ -320,5 +322,5 @@ END DO ! L vertical loop - ! No vertical averaging on bottom layer + ! Last level L = L_NLAYRAD Index: trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90 =================================================================== --- trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90 (revision 1714) +++ trunk/LMDZ.GENERIC/libf/phystd/radcommon_h.F90 (revision 1715) @@ -126,5 +126,4 @@ real*8,save :: PTOP - real*8,save,allocatable :: TAUREF(:) real*8,parameter :: UBARI = 0.5D0 @@ -132,5 +131,5 @@ real*8,save :: gweight(L_NGAUSS) !$OMP THREADPRIVATE(QREFvis,QREFir,omegaREFvis,omegaREFir,& ! gweight read by master in sugas_corrk - !$OMP tstellar,planckir,PTOP,TAUREF) + !$OMP tstellar,planckir,PTOP) ! If the gas optical depth (top to the surface) is less than @@ -156,13 +155,3 @@ !$OMP THREADPRIVATE(glat,eclipse) -contains - - subroutine ini_radcommon_h - use radinc_h, only: L_LEVELS - implicit none - - allocate(TAUREF(L_LEVELS+1)) - - end subroutine ini_radcommon_h - end module radcommon_h