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

Last change on this file since 1009 was 353, checked in by acolaitis, 13 years ago

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