source: trunk/LMDZ.MARS/libf/phymars/lwmain.F @ 1655

Last change on this file since 1655 was 1266, checked in by aslmd, 11 years ago

LMDZ.MARS
IMPORTANT CHANGE

  • Remove all reference/use of nlayermx and dimphys.h
  • Made use of automatic arrays whenever arrays are needed with dimension nlayer
  • Remove lots of obsolete reference to dimensions.h
  • Converted iono.h and param_v4.h into corresponding modules

(with embedded subroutine to allocate arrays)
(no arrays allocated if thermosphere not used)

  • Deleted param.h and put contents into module param_v4_h
  • Adapted testphys1d, newstart, etc...
  • Made DATA arrays in param_read to be initialized by subroutine

fill_data_thermos in module param_v4_h

  • Optimized computations in paramfoto_compact (twice less dlog10 calculations)
  • Checked consistency before/after modification in debug mode
  • Checked performance is not impacted (same as before)
File size: 7.4 KB
Line 
1       subroutine lwmain (ig0,icount,kdlon,kflev
2     .                   ,dp,dt0,emis
3     .                   ,plev,tlev,tlay,aerosol,coolrate
4     .                   ,fluxground,fluxtop
5     .                   ,netrad
6     &                   ,QIRsQREF3d,omegaIR3d,gIR3d
7     &                   ,co2ice)
8
9c----------------------------------------------------------------------
10c     LWMAIN     organizes the LTE longwave calculations
11c     for layer 1 to layer "nlaylte" (stored in  "yomlw_h")
12c----------------------------------------------------------------------
13
14      use dimradmars_mod, only: ndlo2, nflev, nir, ndlon, nuco2
15      use dimradmars_mod, only: naerkind
16      use yomlw_h, only: nlaylte, xi
17      implicit none
18 
19#include "callkeys.h"
20#include "comg1d.h"
21
22c----------------------------------------------------------------------
23c         0.1   arguments
24c               ---------
25c                                                            inputs:
26c                                                            -------
27      integer ig0
28      integer icount
29      integer kdlon            ! part of ngrid
30      integer kflev            ! part of nlayer
31
32      real dp (ndlo2,kflev)         ! layer pressure thickness (Pa)
33      real dt0 (ndlo2)              ! surface temperature discontinuity (K)
34      real emis (ndlo2)             ! surface emissivity
35      real plev (ndlo2,kflev+1)     ! level pressure (Pa)
36      real tlev (ndlo2,kflev+1)     ! level temperature (K)
37      real tlay (ndlo2,kflev)       ! layer temperature (K)
38      real aerosol(ndlo2,kflev,naerkind)      !  aerosol extinction optical
39c                         depth at reference wavelength "longrefvis" set
40c                         in dimradmars_mod , in each layer, for one of
41c                         the "naerkind" kind of aerosol optical properties.
42
43
44c                                                            outputs:
45c                                                            --------
46      real coolrate(ndlo2,kflev)      ! cooling rate (K/s)
47      real fluxground(ndlo2)          ! downward ground flux (W/m2)
48      real fluxtop(ndlo2)             ! outgoing upward flux (W/m2) ("OLR")
49      real netrad (ndlo2,kflev)       ! radiative budget (W/m2)
50c     Aerosol optical properties
51      REAL :: QIRsQREF3d(ndlo2,kflev,nir,naerkind)
52      REAL :: omegaIR3d(ndlo2,kflev,nir,naerkind)
53      REAL :: gIR3d(ndlo2,kflev,nir,naerkind)
54
55c----------------------------------------------------------------------
56c         0.2   local arrays
57c               ------------
58
59      real aer_t (ndlon,nuco2,nflev+1)  ! transmission (aer)
60      real co2_u (ndlon,nuco2,nflev+1)  ! absorber amounts (co2)
61      real co2_up (ndlon,nuco2,nflev+1) ! idem scaled by the pressure (co2)
62
63      real bsurf (ndlon,nir)            ! surface spectral planck function
64      real btop (ndlon,nir)             ! top spectral planck function
65      real blev (ndlon,nir,nflev+1)     ! level   spectral planck function
66      real blay (ndlon,nir,nflev)       ! layer   spectral planck function
67      real dblay (ndlon,nir,nflev)      ! layer gradient spectral planck function
68      real dbsublay (ndlon,nir,2*nflev) ! layer gradient spectral planck function
69                                        ! in sub layers
70
71      real tautotal(ndlon,nflev,nir)  ! \   Total single scattering
72      real omegtotal(ndlon,nflev,nir) !  >  properties (Addition of the
73      real gtotal(ndlon,nflev,nir)    ! /   NAERKIND aerosols prop.)
74
75      real newcoolrate(ndlon,nflev) ! cooling rate (K/s) / with implicite scheme
76
77      REAL co2ice(ndlo2)           ! co2 ice surface layer (kg.m-2)
78      REAL emis_gaz(ndlo2)         ! emissivity for gaz computations
79
80      integer jk,jkk,ja,jl
81
82      logical firstcall
83      save firstcall
84      data firstcall/.true./
85
86
87c----------------------------------------------------------------------
88c         0.3   Initialisation
89c               --------------
90
91      if (firstcall) then
92
93        firstcall = .false.
94
95        xi (:,:,:,:)=0.
96
97      endif
98
99      DO jl=1 , kdlon
100         IF(co2ice(jl) .GT. 20.e-3) THEN
101             emis_gaz(jl)=1.
102         ELSE
103             emis_gaz(jl)=emis(jl)
104         ENDIF
105      ENDDO
106
107c----------------------------------------------------------------------
108c         1.0   planck function
109c               ---------------
110
111      call lwb ( kdlon, kflev, tlev, tlay, dt0
112     .         , bsurf, btop, blay, blev, dblay, dbsublay)
113
114c----------------------------------------------------------------------
115c         2.0   absorber amounts
116c               ----------------
117
118      call lwu ( kdlon, kflev
119     .         , dp, plev, tlay, aerosol
120     &         , QIRsQREF3d,omegaIR3d,gIR3d
121     .         , aer_t, co2_u, co2_up
122     .         , tautotal,omegtotal,gtotal)
123
124c----------------------------------------------------------------------
125c         3.0   transmission functions / exchange coefficiants
126c               ----------------------------------------------
127
128c                                                                distants
129c                                                                --------
130                    if( mod(icount-1,ilwd).eq.0) then
131
132c     print*, 'CALL of DISTANTS'
133      call lwxd ( ig0, kdlon, kflev, emis_gaz
134     .          , aer_t, co2_u, co2_up)
135
136                    endif
137c                                                              neighbours
138c                                                              ----------
139                    if( mod(icount-1,ilwn).eq.0) then
140
141c     print*, 'CALL of NEIGHBOURS'
142      call lwxn ( ig0, kdlon, kflev
143     .          , dp
144     .          , aer_t, co2_u, co2_up)
145
146                    endif
147c                                                              boundaries
148c                                                              ----------
149                    if( mod(icount-1,ilwb).eq.0) then
150
151c     print*, 'CALL of BOUNDARIES'
152      call lwxb ( ig0, kdlon, kflev, emis_gaz
153     .          , aer_t, co2_u, co2_up)
154
155                    endif
156
157c----------------------------------------------------------------------
158c         4.0   cooling rate
159c               ------------
160
161      call lwflux ( ig0, kdlon, kflev, dp
162     .            , bsurf, btop, blev, blay, dbsublay
163     .            , tlay, tlev, dt0      ! pour sortie dans g2d uniquement
164     .            , emis
165     .            , tautotal,omegtotal,gtotal
166     .            , coolrate, fluxground, fluxtop
167     .            , netrad)
168
169c     do jk = 1, nlaylte
170c       print*,coolrate(1,jk)
171c     enddo
172     
173c       do jkk = 0 , nlaylte+1
174c         do jk = 0 , nlaylte+1
175c           do ja = 1 , nuco2
176c             do jl = 1 , ngrid
177c      if (xi (jl,ja,jk,jkk) .LT. 0
178c    .       .OR. xi (jl,ja,jk,jkk) .GT. 1 ) then
179c                 print*,'xi bande',ja,jk,jkk,xi (jl,ja,jk,jkk)
180c      endif
181c             enddo
182c           enddo
183c         enddo
184c       enddo
185
186c----------------------------------------------------------------------
187c
188c          5.    shema semi-implicite  (lwi)
189c                ---------------------------
190c
191c
192      call lwi (ig0,kdlon,kflev,netrad,dblay,dp
193     .          , newcoolrate)
194c
195c  Verif que   (X sol,space) + somme(X i,sol) = 1
196c
197      do jkk = 1 , nlaylte
198        do jl = 1 , kdlon
199c     print*,'NEW et OLD coolrate :',jkk,newcoolrate(jl,jkk)
200c    .  ,coolrate(jl,jkk)
201      coolrate(jl,jkk) = newcoolrate(jl,jkk)
202        enddo
203      enddo
204c
205c----------------------------------------------------------------------
206
207      return
208      end
Note: See TracBrowser for help on using the repository browser.