source: palm/trunk/SOURCE/ls_forcing.f90 @ 1382

Last change on this file since 1382 was 1382, checked in by boeske, 7 years ago

minor changes in profile data output of lsf tendencies, variables renamed

  • Property svn:keywords set to Id
File size: 20.8 KB
Line 
1 MODULE ls_forcing_mod
2
3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2014 Leibniz Universitaet Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22! Renamed variables which store large scale forcing tendencies
23! pt_lsa -> td_lsa_lpt, pt_subs -> td_sub_lpt,
24! q_lsa  -> td_lsa_q,   q_subs  -> td_sub_q,
25! high|lowpt_lsa -> high|low_td_lsa_lpt, ...
26!
27! Former revisions:
28! -----------------
29! $Id: ls_forcing.f90 1382 2014-04-30 12:15:41Z boeske $
30!
31! 1365 2014-04-22 15:03:56Z boeske
32! Usage of large scale forcing for pt and q enabled:
33! Added new subroutine ls_advec for horizontal large scale advection and large
34! scale subsidence,
35! error message in init_ls_forcing specified,
36! variable t renamed nt
37!
38! 1353 2014-04-08 15:21:23Z heinze
39! REAL constants provided with KIND-attribute
40!
41! 1320 2014-03-20 08:40:49Z raasch
42! ONLY-attribute added to USE-statements,
43! kind-parameters added to all INTEGER and REAL declaration statements,
44! kinds are defined in new module kinds,
45! comment fields (!:) to be used for variable explanations added to
46! all variable declaration statements
47!
48! 1318 2014-03-17 13:35:16Z raasch
49! module interfaces removed
50!
51! 1299 2014-03-06 13:15:21Z heinze
52! Ensure a zero large scale vertical velocity at the surface
53! Bugfix: typo in case of boundary condition in if-clause
54!
55! 1276 2014-01-15 13:40:41Z heinze
56! Use LSF_DATA also in case of Dirichlet bottom boundary condition for scalars
57!
58! 1249 2013-11-06 10:45:47Z heinze
59! remove call of user module
60! reformatting
61!
62! 1241 2013-10-30 11:36:58Z heinze
63! Initial revision
64!
65! Description:
66! ------------
67! Calculates large scale forcings (geostrophic wind and subsidence velocity) as
68! well as surfaces fluxes dependent on time given in an external file (LSF_DATA).
69! Code is based in parts on DALES and UCLA-LES.
70!--------------------------------------------------------------------------------!
71
72    PRIVATE
73    PUBLIC init_ls_forcing, ls_forcing_surf, ls_forcing_vert, ls_advec
74    SAVE
75
76    INTERFACE ls_advec
77       MODULE PROCEDURE ls_advec
78       MODULE PROCEDURE ls_advec_ij
79    END INTERFACE ls_advec
80
81 CONTAINS
82
83    SUBROUTINE init_ls_forcing
84
85       USE arrays_3d,                                                          &
86           ONLY:  p_surf, pt_surf, q_surf, qsws_surf, shf_surf, td_lsa_lpt,    &
87                  td_lsa_q, td_sub_lpt, td_sub_q, time_surf, time_vert,        &
88                  ug_vert, vg_vert, wsubs_vert, zu
89
90       USE control_parameters,                                                 &
91           ONLY:  end_time, lsf_surf, lsf_vert, message_string, nlsf
92
93       USE indices,                                                            &
94           ONLY:  ngp_sums_ls, nzb, nz, nzt
95
96       USE kinds
97
98       USE statistics,                                                         &
99           ONLY:  sums_ls_l
100
101
102       IMPLICIT NONE
103
104       CHARACTER(100) ::  chmess      !:
105       CHARACTER(1)   ::  hash        !:
106
107       INTEGER(iwp) ::  ierrn         !:
108       INTEGER(iwp) ::  finput = 90   !:
109       INTEGER(iwp) ::  k             !:
110       INTEGER(iwp) ::  nt             !:
111
112       REAL(wp) ::  fac               !:
113       REAL(wp) ::  highheight        !:
114       REAL(wp) ::  highug_vert       !:
115       REAL(wp) ::  highvg_vert       !:
116       REAL(wp) ::  highwsubs_vert    !:
117       REAL(wp) ::  lowheight         !:
118       REAL(wp) ::  lowug_vert        !:
119       REAL(wp) ::  lowvg_vert        !:
120       REAL(wp) ::  lowwsubs_vert     !:
121       REAL(wp) ::  high_td_lsa_lpt   !:
122       REAL(wp) ::  low_td_lsa_lpt    !:
123       REAL(wp) ::  high_td_lsa_q     !:
124       REAL(wp) ::  low_td_lsa_q      !:
125       REAL(wp) ::  high_td_sub_lpt   !:
126       REAL(wp) ::  low_td_sub_lpt    !:
127       REAL(wp) ::  high_td_sub_q     !:
128       REAL(wp) ::  low_td_sub_q      !:
129       REAL(wp) ::  r_dummy           !:
130
131       ALLOCATE( p_surf(0:nlsf), pt_surf(0:nlsf), q_surf(0:nlsf),              &
132                 qsws_surf(0:nlsf), shf_surf(0:nlsf),                          &
133                 td_lsa_lpt(nzb:nzt+1,0:nlsf), td_lsa_q(nzb:nzt+1,0:nlsf),     &
134                 td_sub_lpt(nzb:nzt+1,0:nlsf), td_sub_q(nzb:nzt+1,0:nlsf),     &
135                 time_vert(0:nlsf), time_surf(0:nlsf),                         &
136                 ug_vert(nzb:nzt+1,0:nlsf), vg_vert(nzb:nzt+1,0:nlsf),         &
137                 wsubs_vert(nzb:nzt+1,0:nlsf) )
138
139       p_surf = 0.0_wp; pt_surf = 0.0_wp; q_surf = 0.0_wp; qsws_surf = 0.0_wp
140       shf_surf = 0.0_wp; time_vert = 0.0_wp; td_lsa_lpt = 0.0_wp
141       td_lsa_q = 0.0_wp; td_sub_lpt = 0.0_wp; td_sub_q = 0.0_wp
142       time_surf = 0.0_wp; ug_vert = 0.0_wp; vg_vert = 0.0_wp
143       wsubs_vert = 0.0_wp
144
145!
146!--    Array for storing large scale forcing and nudging tendencies at each
147!--    timestep for data output
148       ALLOCATE( sums_ls_l(nzb:nzt+1,0:7) )
149       sums_ls_l = 0.0_wp
150
151       ngp_sums_ls = (nz+2)*6
152
153       OPEN ( finput, FILE='LSF_DATA', STATUS='OLD', &
154              FORM='FORMATTED', IOSTAT=ierrn )
155
156       IF ( ierrn /= 0 )  THEN
157          message_string = 'file LSF_DATA does not exist'
158          CALL message( 'ls_forcing', 'PA0368', 1, 2, 0, 6, 0 )
159       ENDIF
160
161       ierrn = 0
162!
163!--    First three lines of LSF_DATA contain header
164       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
165       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
166       READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
167
168       IF ( ierrn /= 0 )  THEN
169          message_string = 'errors in file LSF_DATA'
170          CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
171       ENDIF
172
173!
174!--    Surface values are read in
175       nt     = 0
176       ierrn = 0
177
178       DO WHILE ( time_surf(nt) < end_time )
179          nt = nt + 1
180          READ ( finput, *, IOSTAT = ierrn ) time_surf(nt), shf_surf(nt),      &
181                                             qsws_surf(nt), pt_surf(nt),       &
182                                             q_surf(nt), p_surf(nt)
183
184          IF ( ierrn < 0 )  THEN
185            WRITE ( message_string, * ) 'No time dependent surface variables ',&
186                              'in&LSF_DATA for end of run found'
187
188             CALL message( 'ls_forcing', 'PA0370', 1, 2, 0, 6, 0 )
189          ENDIF
190       ENDDO
191
192
193       IF ( time_surf(1) > end_time )  THEN
194          WRITE ( message_string, * ) 'No time dependent surface variables in ',&
195                                     '&LSF_DATA for end of run found - ',      &
196                                     'lsf_surf is set to FALSE'
197          CALL message( 'ls_forcing', 'PA0371', 0, 0, 0, 6, 0 )
198          lsf_surf = .FALSE.
199       ENDIF
200
201!
202!--    Go to the end of the list with surface variables
203       DO WHILE ( ierrn == 0 )
204          READ ( finput, *, IOSTAT = ierrn ) r_dummy
205       ENDDO
206
207!
208!--    Profiles of ug, vg and w_subs are read in (large scale forcing)
209
210       nt = 0
211       DO WHILE ( time_vert(nt) < end_time )
212          nt = nt + 1
213          hash = "#"
214          ierrn = 1 ! not zero
215!
216!--       Search for the next line consisting of "# time",
217!--       from there onwards the profiles will be read
218          DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 
219             READ ( finput, *, IOSTAT=ierrn ) hash, time_vert(nt)
220             IF ( ierrn < 0 )  THEN
221                WRITE( message_string, * ) 'No time dependent vertical profiles',&
222                                 ' in&LSF_DATA for end of run found'
223                CALL message( 'ls_forcing', 'PA0372', 1, 2, 0, 6, 0 )
224             ENDIF
225          ENDDO
226
227          IF ( nt == 1 .AND. time_vert(nt) > end_time ) EXIT
228
229          READ ( finput, *, IOSTAT=ierrn ) lowheight, lowug_vert, lowvg_vert,  &
230                                           lowwsubs_vert, low_td_lsa_lpt,      &
231                                           low_td_lsa_q, low_td_sub_lpt,       &
232                                           low_td_sub_q
233          IF ( ierrn /= 0 )  THEN
234             message_string = 'errors in file LSF_DATA'
235             CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
236          ENDIF
237
238          READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert,            &
239                                           highvg_vert, highwsubs_vert,        &
240                                           high_td_lsa_lpt, high_td_lsa_q,     &
241                                           high_td_sub_lpt, high_td_sub_q
242     
243          IF ( ierrn /= 0 )  THEN
244             message_string = 'errors in file LSF_DATA'
245             CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
246          ENDIF
247
248
249          DO  k = nzb, nzt+1
250             IF ( highheight < zu(k) )  THEN
251                lowheight      = highheight
252                lowug_vert     = highug_vert
253                lowvg_vert     = highvg_vert
254                lowwsubs_vert  = highwsubs_vert
255                low_td_lsa_lpt = high_td_lsa_lpt
256                low_td_lsa_q   = high_td_lsa_q
257                low_td_sub_lpt = high_td_sub_lpt
258                low_td_sub_q   = high_td_sub_q
259
260                ierrn = 0
261                READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert,      &
262                                                 highvg_vert, highwsubs_vert,  &
263                                                 high_td_lsa_lpt,              &
264                                                 high_td_lsa_q,                &
265                                                 high_td_sub_lpt, high_td_sub_q
266
267                IF ( ierrn /= 0 )  THEN
268                   WRITE( message_string, * ) 'zu(nzt+1) = ', zu(nzt+1), 'm ', &
269                        'is higher than the maximum height in LSF_DATA which ',&
270                        'is ', lowheight, 'm. Interpolation on PALM ',         &
271                        'grid is not possible.'
272                   CALL message( 'ls_forcing', 'PA0395', 1, 2, 0, 6, 0 )
273                ENDIF
274
275             ENDIF
276
277!
278!--          Interpolation of prescribed profiles in space
279             fac = (highheight-zu(k))/(highheight - lowheight)
280
281             ug_vert(k,nt)    = fac * lowug_vert                               &
282                                + ( 1.0_wp - fac ) * highug_vert
283             vg_vert(k,nt)    = fac * lowvg_vert                               &
284                                + ( 1.0_wp - fac ) * highvg_vert
285             wsubs_vert(k,nt) = fac * lowwsubs_vert                            &
286                                + ( 1.0_wp - fac ) * highwsubs_vert
287
288             td_lsa_lpt(k,nt) = fac * low_td_lsa_lpt                           &
289                                + ( 1.0_wp - fac ) * high_td_lsa_lpt
290             td_lsa_q(k,nt)   = fac * low_td_lsa_q                             &
291                                + ( 1.0_wp - fac ) * high_td_lsa_q
292             td_sub_lpt(k,nt) = fac * low_td_sub_lpt                           &
293                                + ( 1.0_wp - fac ) * high_td_sub_lpt
294             td_sub_q(k,nt)   = fac * low_td_sub_q                             &
295                                + ( 1.0_wp - fac ) * high_td_sub_q
296
297          ENDDO
298
299       ENDDO
300
301!
302!--    Large scale vertical velocity has to be zero at the surface
303       wsubs_vert(nzb,:) = 0.0_wp
304 
305       IF ( time_vert(1) > end_time )  THEN
306          WRITE ( message_string, * ) 'Time dependent large scale profile ',   &
307                             'forcing from&LSF_DATA sets in after end of ' ,   &
308                             'simulation - lsf_vert is set to FALSE'
309          CALL message( 'ls_forcing', 'PA0373', 0, 0, 0, 6, 0 )
310          lsf_vert = .FALSE.
311       ENDIF
312
313       CLOSE( finput )
314
315
316    END SUBROUTINE init_ls_forcing
317
318
319    SUBROUTINE ls_forcing_surf ( time )
320
321       USE arrays_3d,                                                          &
322           ONLY:  p_surf, pt_surf, q_surf, qsws, qsws_surf, shf, shf_surf,     &
323                  time_surf, time_vert, ug, ug_vert, vg, vg_vert
324
325       USE control_parameters,                                                 &
326           ONLY:  bc_q_b, ibc_pt_b, ibc_q_b, pt_surface, q_surface,            &
327                  surface_pressure
328
329       USE kinds
330
331       IMPLICIT NONE
332
333       INTEGER(iwp) ::  nt                     !:
334
335       REAL(wp)             :: fac            !:
336       REAL(wp), INTENT(in) :: time           !:
337
338!
339!--    Interpolation in time of LSF_DATA at the surface
340       nt = 1
341       DO WHILE ( time > time_surf(nt) )
342          nt = nt + 1
343       ENDDO
344       IF ( time /= time_surf(nt) )  THEN
345         nt = nt - 1
346       ENDIF
347
348       fac = ( time -time_surf(nt) ) / ( time_surf(nt+1) - time_surf(nt) )
349
350       IF ( ibc_pt_b == 0 )  THEN
351!
352!--       In case of Dirichlet boundary condition shf must not
353!--       be set - it is calculated via MOST in prandtl_fluxes
354          pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
355
356       ELSEIF ( ibc_pt_b == 1 )  THEN
357!
358!--       In case of Neumann boundary condition pt_surface is needed for
359!--       calculation of reference density
360          shf        = shf_surf(nt) + fac * ( shf_surf(nt+1) - shf_surf(nt) )
361          pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
362
363       ENDIF
364
365       IF ( ibc_q_b == 0 )  THEN
366!
367!--       In case of Dirichlet boundary condition qsws must not
368!--       be set - it is calculated via MOST in prandtl_fluxes
369          q_surface = q_surf(nt) + fac * ( q_surf(nt+1) - q_surf(nt) )
370
371       ELSEIF ( ibc_q_b == 1 )  THEN
372
373          qsws = qsws_surf(nt) + fac * ( qsws_surf(nt+1) - qsws_surf(nt) )
374
375       ENDIF
376
377       surface_pressure = p_surf(nt) + fac * ( p_surf(nt+1) - p_surf(nt) )
378
379    END SUBROUTINE ls_forcing_surf
380
381
382    SUBROUTINE ls_forcing_vert ( time )
383
384       USE arrays_3d,                                                          &
385           ONLY:  time_vert, ug, ug_vert, vg, vg_vert, w_subs, wsubs_vert
386
387       USE control_parameters,                                                 &
388           ONLY:  large_scale_subsidence
389
390       USE kinds
391
392       IMPLICIT NONE
393
394       INTEGER(iwp) ::  nt                     !:
395
396       REAL(wp)             ::  fac           !:
397       REAL(wp), INTENT(in) ::  time          !:
398
399!
400!--    Interpolation in time of LSF_DATA for ug, vg and w_subs
401       nt = 1
402       DO WHILE ( time > time_vert(nt) )
403          nt = nt + 1
404       ENDDO
405       IF ( time /= time_vert(nt) )  THEN
406         nt = nt - 1
407       ENDIF
408
409       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
410
411       ug     = ug_vert(:,nt) + fac * ( ug_vert(:,nt+1) - ug_vert(:,nt) )
412       vg     = vg_vert(:,nt) + fac * ( vg_vert(:,nt+1) - vg_vert(:,nt) )
413
414       IF ( large_scale_subsidence )  THEN
415          w_subs = wsubs_vert(:,nt)                                            &
416                   + fac * ( wsubs_vert(:,nt+1) - wsubs_vert(:,nt) )
417       ENDIF
418
419    END SUBROUTINE ls_forcing_vert
420
421
422!------------------------------------------------------------------------------!
423! Call for all grid points
424!------------------------------------------------------------------------------!
425    SUBROUTINE ls_advec ( time, prog_var )
426
427       USE arrays_3d,                                                          &
428           ONLY:  td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, tend, time_vert       
429       
430       USE control_parameters,                                                 &
431           ONLY:  large_scale_subsidence, use_subsidence_tendencies
432       
433       USE indices
434       
435       USE kinds
436
437       IMPLICIT NONE
438
439       CHARACTER (LEN=*) ::  prog_var   !:
440
441       REAL(wp), INTENT(in)  :: time    !:
442       REAL(wp) :: fac                  !: 
443
444       INTEGER(iwp) ::  i               !:
445       INTEGER(iwp) ::  j               !:
446       INTEGER(iwp) ::  k               !:
447       INTEGER(iwp) ::  nt               !:
448
449!
450!--    Interpolation in time of LSF_DATA
451       nt = 1
452       DO WHILE ( time > time_vert(nt) )
453          nt = nt + 1
454       ENDDO
455       IF ( time /= time_vert(nt) )  THEN
456         nt = nt - 1
457       ENDIF
458
459       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
460
461!
462!--    Add horizontal large scale advection tendencies of pt and q
463       SELECT CASE ( prog_var )
464
465          CASE ( 'pt' )
466
467             DO  i = nxl, nxr
468                DO  j = nys, nyn
469                   DO  k = nzb_u_inner(j,i)+1, nzt
470                      tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) + fac *     &
471                                    ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )
472                   ENDDO
473                ENDDO
474             ENDDO
475
476          CASE ( 'q' )
477
478             DO  i = nxl, nxr
479                DO  j = nys, nyn
480                   DO  k = nzb_u_inner(j,i)+1, nzt
481                      tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) + fac *       &
482                                    ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) )
483                   ENDDO
484                ENDDO
485             ENDDO
486
487       END SELECT
488
489!
490!--    Subsidence of pt and q with prescribed subsidence tendencies
491       IF ( large_scale_subsidence .AND. use_subsidence_tendencies )  THEN
492
493          SELECT CASE ( prog_var )
494
495             CASE ( 'pt' )
496
497                DO  i = nxl, nxr
498                   DO  j = nys, nyn
499                      DO  k = nzb_u_inner(j,i)+1, nzt
500                         tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *  &
501                                       ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) )
502                      ENDDO
503                   ENDDO
504                ENDDO
505 
506             CASE ( 'q' )
507
508                DO  i = nxl, nxr
509                   DO  j = nys, nyn
510                      DO  k = nzb_u_inner(j,i)+1, nzt
511                         tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *    &
512                                       ( td_sub_q(k,nt+1) - td_sub_q(k,nt) )
513                      ENDDO
514                   ENDDO
515                ENDDO
516
517          END SELECT
518
519       ENDIF
520
521    END SUBROUTINE ls_advec
522
523
524!------------------------------------------------------------------------------!
525! Call for grid point i,j
526!------------------------------------------------------------------------------!
527    SUBROUTINE ls_advec_ij ( i, j, time, prog_var )
528
529       USE arrays_3d,                                                          &
530           ONLY:  td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, tend, time_vert       
531       
532       USE control_parameters,                                                 &
533           ONLY:  large_scale_subsidence, use_subsidence_tendencies
534       
535       USE indices
536       
537       USE kinds
538
539       IMPLICIT NONE
540
541       CHARACTER (LEN=*) ::  prog_var   !:
542
543       REAL(wp), INTENT(in)  :: time    !:
544       REAL(wp) :: fac                  !:
545
546       INTEGER(iwp) ::  i               !:
547       INTEGER(iwp) ::  j               !:
548       INTEGER(iwp) ::  k               !:
549       INTEGER(iwp) ::  nt               !:
550
551!
552!--    Interpolation in time of LSF_DATA
553       nt = 1
554       DO WHILE ( time > time_vert(nt) )
555          nt = nt + 1
556       ENDDO
557       IF ( time /= time_vert(nt) )  THEN
558         nt = nt - 1
559       ENDIF
560
561       fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
562
563!
564!--    Add horizontal large scale advection tendencies of pt and q
565       SELECT CASE ( prog_var )
566
567          CASE ( 'pt' )
568
569             DO  k = nzb_u_inner(j,i)+1, nzt
570                tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt)                   &
571                              + fac * ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )
572             ENDDO
573
574          CASE ( 'q' )
575
576             DO  k = nzb_u_inner(j,i)+1, nzt
577                tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt)                     &
578                              + fac * ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) )
579             ENDDO
580
581       END SELECT
582
583!
584!--    Subsidence of pt and q with prescribed profiles
585       IF ( large_scale_subsidence .AND. use_subsidence_tendencies )  THEN
586
587          SELECT CASE ( prog_var )
588
589             CASE ( 'pt' )
590
591                DO  k = nzb_u_inner(j,i)+1, nzt
592                   tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *        &
593                                 ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) )
594                ENDDO
595 
596             CASE ( 'q' )
597
598                DO  k = nzb_u_inner(j,i)+1, nzt
599                   tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *          &
600                                 ( td_sub_q(k,nt+1) - td_sub_q(k,nt) )
601                ENDDO
602
603          END SELECT
604
605       ENDIF
606
607    END SUBROUTINE ls_advec_ij
608
609
610 END MODULE ls_forcing_mod
Note: See TracBrowser for help on using the repository browser.