source: trunk/LMDZ.MARS/libf/phymars/compute_dtau_mod.F90 @ 2156

Last change on this file since 2156 was 1985, checked in by mvals, 6 years ago

Mars GCM:
Correction of the dust storm scheme:

  • in rocketduststorm_mod.F90 the case of not daytime and no storm was not taken into account (it was actually in comment), which lead to NaNs? in the calculation of pdqrds (scheme=0, division by alpha=0.)

Useful:

  • in compute_dtau_mod.F90, the writediagfi of tauref_scenario has been added

MV

File size: 3.0 KB
RevLine 
[1974]1        MODULE compute_dtau_mod
2
3        IMPLICIT NONE
4
5        REAL,PARAMETER :: ti_injection=10/24. ! time of beginning injection
6        REAL,PARAMETER :: tf_injection=18/24. ! time of end injection
7        REAL,PARAMETER :: t_scenario=14/24.   ! time of tauref_scenario
8
9        CONTAINS
10
11        SUBROUTINE compute_dtau(ngrid,nlayer,                           &
12       &                        zday,pplev,tauref,                      &
13       &                        ptimestep,dustliftday,local_time)
14
15        USE geometry_mod, only: longitude_deg
16        USE time_phylmdz_mod, only: dtphys, daysec
17        USE comcstfi_h, only: g
18        USE tracer_mod, only: alpha_lift,igcm_dust_mass,igcm_dust_number
19       
20        IMPLICIT NONE
21       
22        include "callkeys.h"
23       
24        INTEGER, INTENT(in) :: ngrid
25        INTEGER, INTENT(in) :: nlayer
26        REAL, INTENT(in) :: zday ! date at lon=0, in fraction of sols
27        REAL, INTENT(in) :: pplev(ngrid,nlayer+1) ! pressure (Pa)
28        REAL, INTENT(in) :: tauref(ngrid) ! Computed dust opacity at 610Pa
29        REAL, INTENT(in) :: ptimestep
30        REAL, INTENT(in) :: local_time(ngrid)
31        REAL, INTENT(out) :: dustliftday(ngrid) ! Dust injection rate (s-1)
32       
33        INTEGER :: ig, l
34        INTEGER, SAVE :: nb_daystep ! nomber of step a day
35        REAL :: tauref_scenario(ngrid) ! from dust scenario
36        REAL :: zday_scenario
37        REAL,ALLOCATABLE,SAVE :: local_time_prev(:)
38        REAL,ALLOCATABLE,SAVE :: dtau(:) ! Dust opacity increment (at 610Pa)
39        REAL,PARAMETER :: odpref=610. !DOD reference pressure (Pa)
40       
41        LOGICAL, SAVE :: firstcall=.TRUE. ! signals first call to physics
42       
43       
44        IF(firstcall)THEN
45                ALLOCATE(local_time_prev(ngrid))
46                ALLOCATE(dtau(ngrid))
47                DO ig=1,ngrid
48                   local_time_prev(ig)=modulo(1.+(zday-ptimestep/daysec)&
49                                      -INT(zday-ptimestep/daysec)       &
50                                      +(longitude_deg(ig)/15)/24,1.)
51                   dtau(ig)=0.
52                ENDDO
53                nb_daystep=(daysec/dtphys)
54                firstcall=.FALSE.
55        ENDIF
56       
57        ! 1. Obtain tauref_scenario from dust scenario
58        zday_scenario=zday-modulo(zday,1.)
59        call read_dust_scenario(ngrid,nlayer,zday_scenario,pplev,       &
60       &                        tauref_scenario)
[1985]61       ! for diagnostics
62        call WRITEDIAGFI(ngrid,"tauref_scenario","tauref_scenario",     &
63       &                  "",2,tauref_scenario)
[1974]64
65        ! 2. Compute dtau
66        DO ig=1,ngrid
67         IF ((local_time(ig).ge.t_scenario).and.                        &
68       &          (local_time_prev(ig).lt.(t_scenario)))THEN
69                 dtau(ig)=tauref_scenario(ig)-tauref(ig)
70         ENDIF
71         IF (dtau(ig).LT.0) THEN
72                  dtau(ig)=0.
73         ENDIF
74         
75        ! 3. Use dtau to compute dustliftday
76         dustliftday(ig)=coeff_injection*                               &
77       &                (dtau(ig)*pplev(ig,1)/odpref)                   &
78       &                /(daysec*(tf_injection-ti_injection))
79        ENDDO
80         
81        ! 4. Save local time
82        local_time_prev(1:ngrid)=local_time(1:ngrid)
83       
84        end subroutine compute_dtau
85
86        end module compute_dtau_mod
Note: See TracBrowser for help on using the repository browser.