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

Last change on this file since 3026 was 2932, checked in by romain.vande, 21 months ago

Mars PCM:
Add a new routine to write output called write_output.
It needs to be imported (for example : use write_output_mod, only: write_output)
Then, it requires only 4 arguments : the name of the variable, its title, its units and the variable itself.
It detects the dimension of the variable and decide to output either in diagfi or diagsoil.
It is also compatible with XIOS (xml file needs to be adapted)
Writediagfi and writediagsoil routines are still available in case.
RV

File size: 5.3 KB
Line 
1       MODULE compute_dtau_mod
2
3        IMPLICIT NONE
4
5        REAL,SAVE :: ti_injection_sol ! time of beginning injection
6        REAL,SAVE :: tf_injection_sol ! time of end injection
7
8        REAL,SAVE,ALLOCATABLE :: dtau(:) ! Dust opacity difference (at 610Pa)
9                                         ! between GCM and dust scenario
10                                         
11!$OMP THREADPRIVATE(ti_injection_sol,tf_injection_sol,dtau)
12
13       CONTAINS
14
15        SUBROUTINE compute_dtau(ngrid,nlayer,                           &
16                                 zday,pplev,tau_pref_gcm,               &
17                                 ptimestep,local_time,IRtoVIScoef,      &
18                                 dustliftday)
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        USE dust_param_mod, only: odpref, t_scenario_sol
26        USE read_dust_scenario_mod, only: read_dust_scenario
27        use write_output_mod, only: write_output
28       
29        IMPLICIT NONE
30       
31        include "callkeys.h"
32       
33        INTEGER, INTENT(in) :: ngrid
34        INTEGER, INTENT(in) :: nlayer
35        REAL, INTENT(in) :: zday ! date at lon=0, in fraction of sols
36        REAL, INTENT(in) :: pplev(ngrid,nlayer+1) ! pressure (Pa)
37        REAL, INTENT(in) :: tau_pref_gcm(ngrid) ! Visible dust opacity column
38                            ! at 610Pa as computed in the GCM
39        REAL, INTENT(in) :: ptimestep
40        REAL, INTENT(in) :: local_time(ngrid)
41        REAL, INTENT(in) :: IRtoVIScoef(ngrid) ! conversion coefficient to apply on
42                                               ! scenario absorption IR (9.3um) CDOD
43                                               ! = tau_pref_gcm_VIS / tau_pref_gcm_IR
44        REAL, INTENT(out) :: dustliftday(ngrid) ! Dust injection rate (s-1)
45       
46        INTEGER :: ig, l
47        INTEGER, SAVE :: nb_daystep ! nomber of step a day
48        REAL :: tau_pref_target(ngrid) ! dust opacity column at odpref=610 Pa
49                ! as extracted from dust scenario
50        REAL :: zday_scenario
51        REAL,ALLOCATABLE,SAVE :: local_time_prev(:)
52       
53        LOGICAL, SAVE :: firstcall=.TRUE. ! signals first call to physics
54       
55!$OMP THREADPRIVATE(nb_daystep,local_time_prev,firstcall)
56       
57       
58        IF(firstcall)THEN
59                ALLOCATE(local_time_prev(ngrid))
60                DO ig=1,ngrid
61                   local_time_prev(ig)=modulo(1.+(zday-ptimestep/daysec)&
62                                      -INT(zday-ptimestep/daysec)       &
63                                      +(longitude_deg(ig)/15)/24,1.)
64                ENDDO
65                nb_daystep=(daysec/dtphys)
66                ! Local time in sol fraction
67                ti_injection_sol=ti_injection/24.
68                tf_injection_sol=tf_injection/24.
69                firstcall=.FALSE.
70        ENDIF
71       
72        ! 1. Obtain tau_pref_target from dust scenario at zday+1
73        if (iaervar.eq.1) then
74          tau_pref_target = tauvis
75        else
76          zday_scenario=zday-modulo(zday,1.) ! integer value of the day: the scenario opacity is measured at 14:00
77          zday_scenario=zday_scenario+1      ! opacity of the dust scenario is read the day after
78          call read_dust_scenario(ngrid,nlayer,zday_scenario,pplev,     &
79                                         IRtoVIScoef,                   &
80                                         tau_pref_target)
81        endif
82       ! for diagnostics
83        call write_output("tau_pref_target", &
84                          "target visible dust opacity column at 610Pa", &
85                          "",tau_pref_target(:))
86
87        ! 2. Compute dtau() and dustliftday()
88        DO ig=1,ngrid
89         IF ((local_time(ig).ge.t_scenario_sol).and.                    &
90                 (local_time_prev(ig).lt.(t_scenario_sol)))THEN
91                 dtau(ig)=tau_pref_target(ig)-tau_pref_gcm(ig)
92         ENDIF
93
94        ! Use dtau (when positive) to compute dustliftday
95         IF (dtau(ig).LT.0) THEN
96             dustliftday(ig)=0.
97         ELSE
98             dustliftday(ig)=coeff_injection*                           &
99                        (dtau(ig)*pplev(ig,1)/odpref)                   &
100                        /(daysec*(tf_injection_sol-ti_injection_sol))
101         ENDIF
102        ENDDO ! of DO ig=1,ngrid
103
104       ! for diagnostics
105        call write_output("dtau","opacity difference wrt scenario",&
106                          "",dtau(:))
107        call write_output("dustliftday","dust injection rate",     &
108                          "s-1",dustliftday(:))
109         
110        ! 4. Save local time
111        local_time_prev(1:ngrid)=local_time(1:ngrid)
112       
113        end subroutine compute_dtau
114
115!=======================================================================
116! Initialization of the module variables
117
118        subroutine ini_compute_dtau_mod(ngrid)
119       
120          implicit none
121       
122          integer, intent(in) :: ngrid
123       
124          allocate(dtau(ngrid))
125       
126        end subroutine ini_compute_dtau_mod
127       
128        subroutine end_compute_dtau_mod
129       
130          implicit none
131       
132          if (allocated(dtau)) deallocate(dtau)
133
134          end subroutine end_compute_dtau_mod       
135
136       END MODULE compute_dtau_mod
Note: See TracBrowser for help on using the repository browser.