source: trunk/LMDZ.VENUS/libf/phyvenus/lwi.F @ 1121

Last change on this file since 1121 was 1017, checked in by emillour, 11 years ago

Common dynamics: (and collateral adaptations in Venus physics)
Improved cpdet routines in and additional sponge mode:

  • Additionnal sponge mode (trigered with "callsponge" flag), in line with the one used in the Generic and Martian GCM. This sponge is called whenever there is a dissipation step.
  • Improvement of the cpdet routines : created routines tpot2t_glo_p and t2tpot_glo_p which handle fields on the whole dynamics (scaler) grid, which are more efficient than calling tpot2t_p or t2tpot_p with slabs of data (generated use of intermediate copies of these chunks of data at every call)
  • Turned cpdet.F into a module cpdet_mod.F90 (and correspondingly adapted all routines in the Venus physics).

EM
File size: 4.4 KB
Line 
1      subroutine lwi(nl,netrad,deltapsi,deltap,temp,coolrate)
2
3      use dimphy
4      use cpdet_mod, only: cpdet
5      implicit none
6
7
8#include "dimensions.h"
9#include "YOMCST.h"
10#include "timerad.h"
11     
12CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
13C
14C                             -   lwi    -   
15C
16C     PURPOSE:       Schema semi - implicite
17C
18C
19CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
20
21
22c************************************************************************
23c
24c        0.    Declarations
25c              ------------
26c
27c-------------------------------------------------------------------------
28c        0.1   Arguments
29c              ---------
30c Inputs
31c.......
32
33      integer nl
34      real    deltapsi(0:nl+1,0:nl+1) ! Dpsi/DT = sum(nu)[ksi DB/DT]
35      real    netrad(0:nl)     ! radiative budget (W/m**2)
36      real    deltap(nl) !epaisseur de la couche en pression (Pa)
37! ADAPTATION GCM POUR CP(T)
38      real    temp(nl)   !temperature de la couche
39     
40c Outputs
41c........
42
43      real    coolrate(nl)
44     
45c-------------------------------------------------------------------------
46c        0.2   local arrays
47c              ------------
48c
49      real    di(nl)
50     .      , spec(nl)
51     .      , hi(nl)
52     .      , bi(nl)
53
54      real    ci(nl)
55     .      , ai(nl)
56
57      real    deltat
58
59      real    semi,semit,denom
60
61      integer i,jl,band
62
63c************************************************************************
64c
65c        1.    Initialisations
66c              ---------------
67c
68c-----------------------------------------------------------------------
69        deltat = dtimerad
70c       print*,'SEMI = ',semi, '(expl:0  semi-implicite:0.5  impl:1)'
71c        semi  = 0.5
72c       semi = 0.5
73        semi = 1.
74        semit = semi * deltat
75
76c       print*,'dtimerad,deltat,semit:',dtimerad,deltat,semit
77
78c************************************************************************
79c
80c        2.   
81c              ---------------
82c
83c-------------------------------------------------------------------------
84c        2.1   Calcul des di
85c              -------------
86c
87
88      do i = 1 , nl-1
89        spec(i) =
90     .      ( deltapsi(i,nl+1)
91     .      + deltapsi(i,i+1)
92     .      + deltapsi(i,i-1) )
93     
94        di(i) =  1. + semit * (RG/deltap(i)/cpdet(temp(i))) * spec(i)
95c     print*,i,' di(i)=',di(i)
96      enddo
97         
98
99c couche nl
100c ------------
101c      , on enleve i,i+1 sinon on a 2 fois le cooling2space
102 
103        spec(nl) =
104     .      ( deltapsi(nl,nl+1)
105     .      + deltapsi(nl,nl-1) )
106     
107        di(nl) =  1. + semit * (RG/deltap(nl)/cpdet(temp(nl)))*spec(nl)
108
109c-------------------------------------------------------------------------
110c        2.2   Calcul des hi
111c              -------------
112c
113
114      do i = 1 , nl-1
115     
116        spec(i) = deltapsi(i+1,i)
117
118        hi(i) = - semit * (RG/deltap(i)/cpdet(temp(i))) * spec(i)
119
120      enddo
121                       
122c     print*,'hi(i)',hi(i)
123
124c-------------------------------------------------------------------------
125c        2.3   Calcul des bi
126c              -------------
127c
128
129
130      do i = 2 , nl
131     
132        spec(i) = deltapsi(i-1,i)
133       
134        bi(i) = - semit * (RG/deltap(i)/cpdet(temp(i))) * spec(i)
135       
136      enddo
137
138c     print*,'bi(i)',bi(i)
139
140c couche 1
141c --------
142c  tant qu'on a pas un calcul propre de deltab(0) qui tienne compte de
143c    la discontinuite de temperature au sol , on met  b1 = 0
144
145     
146        bi(1) = 0
147
148c-------------------------------------------------------------------------
149c        2.4   
150c              -------------
151c
152
153c couche nl
154c ------------
155
156c     ci(nl) = ((RG/CP) * netrad(nl) / deltap(nl))
157c    .                   / di(nl)
158      ci(nl) = netrad(nl) *RG/cpdet(temp(nl)) / deltap(nl)
159
160      ai(nl) = - bi(nl) / di(nl)
161
162      do i = nl-1 , 1 , -1
163        denom = di(i) + hi(i) * ai(i+1)
164
165c       ci(i) = (  (RG/CP) * netrad(i) / deltap(i)
166c    .               - hi(i) * ci(i+1)  )  / denom
167        ci(i) = netrad(i) *RG/cpdet(temp(i)) / deltap(i)
168 
169        ai(i) = -bi(i) / denom
170      enddo
171
172c-------------------------------------------------------------------------
173c        2.5   
174c              -------------
175c
176
177      coolrate(1) = ci(1)
178
179      do i = 2 , nl
180           coolrate(i) = ci(i) + ai(i) * coolrate(i-1)
181c       print*,i,' coolrate(i)=',coolrate(i)
182      enddo
183
184c-------------------------------------------------------------------------
185
186      RETURN
187      END
Note: See TracBrowser for help on using the repository browser.