source: trunk/LMDZ.MARS/libf/phymars/lwxd.F @ 1619

Last change on this file since 1619 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: 9.3 KB
Line 
1      subroutine lwxd (ig0,kdlon,kflev,emis
2     .                ,aer_t,co2_u,co2_up)
3
4c----------------------------------------------------------------------
5c     LWXD   computes transmission function and exchange coefficiants
6c                        for distant layers
7c                          (co2 / aerosols)
8c                       (bands 1 and 2 of co2)
9c----------------------------------------------------------------------
10
11c              |---|---|---|---|---|---|---|---|
12c   kflev+1    |   |   |   |   |   |   |   | 0 |  (space)
13c              |---|---|---|---|---|---|---|---|
14c    kflev     |   |***|***|***|***|   | 0 |   |
15c              |---|---|---|---|---|---|---|---|
16c     ...      |   |***|***|***|   | 0 |   |   |
17c              |---|---|---|---|---|---|---|---|
18c      4       |   |***|***|   | 0 |   |***|   |
19c              |---|---|---|---|---|---|---|---|
20c      3       |   |***|   | 0 |   |***|***|   |
21c              |---|---|---|---|---|---|---|---|
22c      2       |   |   | 0 |   |   |***|***|   |
23c              |---|---|---|---|---|---|---|---|
24c      1       |   | 0 |   |   |***|***|***|   |
25c              |---|---|---|---|---|---|---|---|
26c      0       | 0 |   |   |***|***|***|***|   |  (ground)
27c              |---|---|---|---|---|---|---|---|
28c                0   1   2   3   4  ...  k |k+1
29c             (ground)                    (space)
30
31c  (*)  xi computed in this subroutine
32c----------------------------------------------------------------------
33
34      use dimradmars_mod, only: ndlon, nuco2, nflev, ndlo2
35      use yomlw_h, only: nlaylte, xi, xi_emis
36      implicit none
37
38#include "callkeys.h"
39
40c----------------------------------------------------------------------
41c         0.1   arguments
42c               ---------
43c                                                            inputs:
44c                                                            -------
45      integer ig0
46      integer kdlon      ! part of ngrid
47      integer kflev      ! part of nalyer
48 
49      real emis (ndlo2)                  ! surface emissivity
50      real aer_t (ndlo2,nuco2,kflev+1)   ! transmission (aer)
51      real co2_u (ndlo2,nuco2,kflev+1)   ! absorber amounts (co2)
52      real co2_up (ndlo2,nuco2,kflev+1)  ! idem scaled by the pressure (co2)
53
54c----------------------------------------------------------------------
55c         0.2   local arrays
56c               ------------
57
58      integer ja,jl,jk,jkk,ndim
59!      parameter(ndim = ndlon*nuco2*(nflev+2)*(nflev+2))
60
61
62      real zu (ndlon,nuco2)
63      real zup (ndlon,nuco2)
64      real zt_co2 (ndlon,nuco2)
65      real zt_aer (ndlon,nuco2)
66
67      real ksi (ndlon,nuco2,0:nflev+1,0:nflev+1)
68      real ksi_emis (ndlon,nuco2,0:nflev+1,0:nflev+1)
69      real trans (ndlon,nuco2,0:nflev+1,0:nflev+1)
70      real trans_emis (ndlon,nuco2,0:nflev+1,0:nflev+1)
71
72c----------------------------------------------------------------------
73      ndim = ndlon*nuco2*(nflev+2)*(nflev+2)
74      call zerophys(ndim,ksi_emis)
75c----------------------------------------------------------------------
76c         1.0   Transmission functions
77c               ----------------------
78
79c----------------------------------------------------------------------
80c        1.1   Direct transmission
81c              -------------------
82
83      do jk = 1 , nlaylte+1
84        do jkk = jk , nlaylte+1
85
86          do ja = 1 , nuco2
87            do jl = 1 , kdlon
88c                                                                   co2
89c                                                                   ---
90              zu(jl,ja) =  co2_u(jl,ja,jk)  - co2_u(jl,ja,jkk)
91              zup(jl,ja) = co2_up(jl,ja,jk) - co2_up(jl,ja,jkk)
92c                                                                   aer
93c                                                                   ---
94              zt_aer(jl,ja)= aer_t(jl,ja,jk)
95     .                      /aer_t(jl,ja,jkk)
96
97            enddo
98          enddo
99
100          call lwtt(kdlon,zu,zup,nuco2,zt_co2)
101c                                                            co2 and aer
102c                                                            -----------
103          do ja = 1 , nuco2
104            do jl = 1 , kdlon
105              trans(jl,ja,jk,jkk) = zt_co2(jl,ja) * zt_aer(jl,ja)
106            enddo
107          enddo
108c                                                       trans reciprocity
109c                                                       -----------------
110          do ja = 1 , nuco2
111            do jl = 1 , kdlon
112              trans(jl,ja,jkk,jk) = trans(jl,ja,jk,jkk)
113c           if (trans(jl,ja,jk,jkk) .LT. 0 ) then
114c           print*,'trans bande',ja,jk,jkk,trans(jl,ja,jk,jkk)
115c           endif
116c           if (trans(jl,ja,jk,jkk) .GT. 1) then
117c           print*,'trans bande',ja,jk,jkk,trans(jl,ja,jk,jkk)
118c           trans(jl,ja,jk,jkk)=1
119c           print*,'trans bande',ja,jk,jkk,trans(jl,ja,jk,jkk)
120c           endif
121
122            enddo
123          enddo
124     
125        enddo
126      enddo
127
128c----------------------------------------------------------------------
129c         1.2   Transmission with reflexion
130c               ---------------------------
131
132      do jk = 1 , nlaylte+1
133        do jkk = jk , nlaylte+1
134
135      if (callemis) then
136          do ja = 1 , nuco2
137            do jl = 1 , kdlon
138c                                                                   co2
139c                                                                   ---
140              zu(jl,ja)  = 2 * co2_u(jl,ja,1)  - co2_u(jl,ja,jk)
141     .                                         - co2_u(jl,ja,jkk)
142              zup(jl,ja) = 2 * co2_up(jl,ja,1) - co2_up(jl,ja,jk)
143     .                                         - co2_up(jl,ja,jkk)
144c                                                                   aer
145c                                                                   ---
146                zt_aer(jl,ja)  = aer_t(jl,ja,1)
147     .                         * aer_t(jl,ja,1)
148     .                         / aer_t(jl,ja,jk)
149     .                         / aer_t(jl,ja,jkk)
150            enddo
151          enddo
152
153          call lwtt(kdlon,zu,zup,nuco2,zt_co2)
154c                                                            co2 and aer
155c                                                            -----------
156          do ja = 1 , nuco2
157            do jl = 1 , kdlon
158              trans_emis(jl,ja,jk,jkk) = zt_co2(jl,ja)
159     .                                 * zt_aer(jl,ja)
160            enddo
161          enddo
162
163      else
164
165          do ja = 1 , nuco2
166            do jl = 1 , kdlon
167              trans_emis(jl,ja,jk,jkk) = 1.
168            enddo
169          enddo
170
171      endif
172c                                                       trans reciprocity
173c                                                       -----------------
174          do ja = 1 , nuco2
175            do jl = 1 , kdlon
176              trans_emis(jl,ja,jkk,jk) = trans_emis(jl,ja,jk,jkk)
177c         if (trans_emis(jl,ja,jk,jkk) .LT. 0
178c    .                .OR.  trans_emis(jl,ja,jk,jkk) .GT. 1) then
179c     print*,'trans_emis bande',ja,jk,jkk,trans_emis(jl,ja,jk,jkk)
180c         endif
181            enddo
182          enddo
183
184        enddo
185      enddo
186
187c----------------------------------------------------------------------
188c         2.0   Exchange Coefficiants
189c               ---------------------
190
191      do jk = 1 , nlaylte-2
192        do jkk = jk+2 , nlaylte
193          do ja = 1 , nuco2
194            do jl = 1 , kdlon
195
196      ksi(jl,ja,jk,jkk) =
197     .            trans(jl,ja,jk+1,jkk)   - trans(jl,ja,jk,jkk)
198     .          - trans(jl,ja,jk+1,jkk+1) + trans(jl,ja,jk,jkk+1)
199
200      ksi_emis(jl,ja,jk,jkk) =
201     .   trans_emis(jl,ja,jk,jkk)   - trans_emis(jl,ja,jk+1,jkk)
202     . - trans_emis(jl,ja,jk,jkk+1) + trans_emis(jl,ja,jk+1,jkk+1)
203
204c       if (ksi(jl,ja,jk,jkk) .LT. 0 ) then
205c           print*,'ksi bande',ja,jk,jkk,ksi(jl,ja,jk,jkk)
206c           ksi(jl,ja,jk,jkk)=0
207c           print*,'ksi bande',ja,jk,jkk,ksi(jl,ja,jk,jkk)
208c       endif
209c       if (ksi(jl,ja,jk,jkk) .GT. 1) then
210c           print*,'ksi bande',ja,jk,jkk,ksi(jl,ja,jk,jkk)
211c           ksi(jl,ja,jk,jkk)=1
212c           print*,'ksi bande',ja,jk,jkk,ksi(jl,ja,jk,jkk)
213c       endif
214
215c       if (ksi_emis(jl,ja,jk,jkk) .LT. 0
216c    .                .OR.  ksi_emis(jl,ja,jk,jkk) .GT. 1) then
217c     print*,'ksi_emis bande',ja,jk,jkk,ksi_emis(jl,ja,jk,jkk)
218c       endif
219
220      xi(ig0+jl,ja,jk,jkk) = ksi(jl,ja,jk,jkk)
221     .      + ksi_emis(jl,ja,jk,jkk) * (1 - emis(jl))
222
223c                                                        ksi reciprocity
224c                                                        ---------------
225      ksi(jl,ja,jkk,jk)      = ksi(jl,ja,jk,jkk)
226      ksi_emis(jl,ja,jkk,jk) = ksi_emis(jl,ja,jk,jkk)
227      xi(ig0+jl,ja,jkk,jk)   = xi(ig0+jl,ja,jk,jkk)
228
229            enddo
230          enddo
231        enddo
232      enddo
233
234c----------------------------------------------------------------------
235c         2.1   Save xi_emis for neighbours (lwxn.F)
236c               -----------------------------------
237
238      do jk = 1 , nlaylte-1
239        do ja = 1 , nuco2
240          do jl = 1 , kdlon
241
242c     ksi_emis(jl,ja,jk,jk+1) =
243c    .   trans_emis(jl,ja,jk,jk+1)   - trans_emis(jl,ja,jk+1,jk+1)
244c    . - trans_emis(jl,ja,jk,jk+2) + trans_emis(jl,ja,jk+1,jk+2)
245
246      xi_emis(ig0+jl,ja,jk) =
247     .                 ksi_emis(jl,ja,jk,jk+1) * (1-emis(jl))
248
249          enddo
250        enddo
251      enddo
252
253c----------------------------------------------------------------------
254      return
255      end
Note: See TracBrowser for help on using the repository browser.