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

Last change on this file since 2265 was 2265, checked in by emillour, 5 years ago

Mars GCM:
Save "dtau", the opacity difference between model and target dust scenario
in the restartfi.nc file so that we have 1+1=2 when running with dust
injection schemes.
EM

File size: 4.6 KB
Line 
1       MODULE compute_dtau_mod
2
3        IMPLICIT NONE
4
5        include "callkeys.h"
6
7        REAL,SAVE :: ti_injection_sol ! time of beginning injection
8        REAL,SAVE :: tf_injection_sol ! time of end injection
9        REAL,PARAMETER :: t_scenario_sol=14/24.   ! time of tauref_scenario
10
11        REAL,SAVE,ALLOCATABLE :: dtau(:) ! Dust opacity difference (at 610Pa)
12                                         ! between GCM and dust scenario
13
14       CONTAINS
15
16        SUBROUTINE compute_dtau(ngrid,nlayer,                           &
17                                 zday,pplev,tauref,                     &
18                                 ptimestep,dustliftday,local_time)
19
20        USE geometry_mod, only: longitude_deg
21        USE time_phylmdz_mod, only: dtphys, daysec
22        USE comcstfi_h, only: g
23        USE tracer_mod, only: alpha_lift,igcm_dust_mass,igcm_dust_number
24        USE dimradmars_mod, only: tauvis
25       
26        IMPLICIT NONE
27       
28        include "callkeys.h"
29       
30        INTEGER, INTENT(in) :: ngrid
31        INTEGER, INTENT(in) :: nlayer
32        REAL, INTENT(in) :: zday ! date at lon=0, in fraction of sols
33        REAL, INTENT(in) :: pplev(ngrid,nlayer+1) ! pressure (Pa)
34        REAL, INTENT(in) :: tauref(ngrid) ! Computed dust opacity at 610Pa
35        REAL, INTENT(in) :: ptimestep
36        REAL, INTENT(in) :: local_time(ngrid)
37        REAL, INTENT(out) :: dustliftday(ngrid) ! Dust injection rate (s-1)
38       
39        INTEGER :: ig, l
40        INTEGER, SAVE :: nb_daystep ! nomber of step a day
41        REAL :: tauref_scenario(ngrid) ! from dust scenario
42        REAL :: zday_scenario
43        REAL,ALLOCATABLE,SAVE :: local_time_prev(:)
44        REAL,PARAMETER :: odpref=610. !DOD reference pressure (Pa)
45       
46        LOGICAL, SAVE :: firstcall=.TRUE. ! signals first call to physics
47       
48       
49        IF(firstcall)THEN
50                ALLOCATE(local_time_prev(ngrid))
51                DO ig=1,ngrid
52                   local_time_prev(ig)=modulo(1.+(zday-ptimestep/daysec)&
53                                      -INT(zday-ptimestep/daysec)       &
54                                      +(longitude_deg(ig)/15)/24,1.)
55                ENDDO
56                nb_daystep=(daysec/dtphys)
57                ! Local time in sol fraction
58                ti_injection_sol=ti_injection/24.
59                tf_injection_sol=tf_injection/24.
60                firstcall=.FALSE.
61        ENDIF
62       
63        ! 1. Obtain tauref_scenario from dust scenario at zday+1
64        if (iaervar.eq.1) then
65          tauref_scenario = tauvis
66        else
67          zday_scenario=zday-modulo(zday,1.) ! integer value of the day: the scenario opacity is measured at 14:00
68          zday_scenario=zday_scenario+1      ! opacity of the dust scenario is read the day after
69          call read_dust_scenario(ngrid,nlayer,zday_scenario,pplev,     &
70                                         tauref_scenario)
71        endif
72       ! for diagnostics
73        call WRITEDIAGFI(ngrid,"tauref_scenario","tauref_scenario",     &
74                          "",2,tauref_scenario)
75
76        ! 2. Compute dtau() and dustliftday()
77        DO ig=1,ngrid
78         IF ((local_time(ig).ge.t_scenario_sol).and.                    &
79                 (local_time_prev(ig).lt.(t_scenario_sol)))THEN
80                 dtau(ig)=tauref_scenario(ig)-tauref(ig)
81         ENDIF
82
83        ! Use dtau (when positiove) to compute dustliftday
84         IF (dtau(ig).LT.0) THEN
85             dustliftday(ig)=0.
86         ELSE
87             dustliftday(ig)=coeff_injection*                           &
88                        (dtau(ig)*pplev(ig,1)/odpref)                   &
89                        /(daysec*(tf_injection_sol-ti_injection_sol))
90         ENDIF
91        ENDDO ! of DO ig=1,ngrid
92
93       ! for diagnostics
94        call WRITEDIAGFI(ngrid,"dtau","opacity difference wrt scenario",&
95                          "",2,dtau)
96        call WRITEDIAGFI(ngrid,"dustliftday","dust injection rate",     &
97                          "s-1",2,dustliftday)
98         
99        ! 4. Save local time
100        local_time_prev(1:ngrid)=local_time(1:ngrid)
101       
102        end subroutine compute_dtau
103
104!=======================================================================
105! Initialization of the module variables
106
107        subroutine ini_compute_dtau_mod(ngrid)
108       
109          implicit none
110       
111          integer, intent(in) :: ngrid
112       
113          allocate(dtau(ngrid))
114       
115        end subroutine ini_compute_dtau_mod
116       
117        subroutine end_compute_dtau_mod
118       
119          implicit none
120       
121          if (allocated(dtau)) deallocate(dtau)
122
123          end subroutine end_compute_dtau_mod       
124
125       END MODULE compute_dtau_mod
Note: See TracBrowser for help on using the repository browser.