source: codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/initracer.F

Last change on this file was 227, checked in by milmd, 10 years ago

Last LMDZ version (1315) with OpenMP directives and other stuff

File size: 10.8 KB
Line 
1      SUBROUTINE initracer(ngrid,nq,nametrac)
2
3      use surfdat_h
4      USE comgeomfi_h
5      USE tracer_h
6      IMPLICIT NONE
7c=======================================================================
8c   subject:
9c   --------
10c   Initialization related to tracer 
11c   (transported dust, water, chemical species, ice...)
12c
13c   Name of the tracer
14c
15c   Test of dimension :
16c   Initialize COMMON tracer in tracer.h, using tracer names provided
17c   by the argument nametrac
18c
19c   author: F.Forget
20c   ------
21c            Ehouarn Millour (oct. 2008) identify tracers by their names
22c=======================================================================
23
24
25#include "dimensions.h"
26!#include "dimphys.h"
27#include "comcstfi.h"
28#include "callkeys.h"
29
30      integer :: ngrid,nq
31
32!      real qsurf(ngrid,nq)       ! tracer on surface (e.g.  kg.m-2)
33!      real co2ice(ngrid)           ! co2 ice mass on surface (e.g.  kg.m-2)
34      character(len=20) :: txt ! to store some text
35      integer iq,ig,count
36      real r0_lift , reff_lift
37!      logical :: oldnames ! =.true. if old tracer naming convention (q01,...)
38
39      character*20 nametrac(nq)   ! name of the tracer from dynamics
40
41
42c-----------------------------------------------------------------------
43c  radius(nq)      ! aerosol particle radius (m)
44c  rho_q(nq)       ! tracer densities (kg.m-3)
45c  qext(nq)        ! Single Scat. Extinction coeff at 0.67 um
46c  alpha_lift(nq)  ! saltation vertical flux/horiz flux ratio (m-1)
47c  alpha_devil(nq) ! lifting coeeficient by dust devil
48c  rho_dust          ! Mars dust density
49c  rho_ice           ! Water ice density
50c  doubleq           ! if method with mass (iq=1) and number(iq=2) mixing ratio
51c  varian            ! Characteristic variance of log-normal distribution
52c-----------------------------------------------------------------------
53
54       !! we allocate once for all arrays in common in tracer_h.F90
55       !! (supposedly those are not used before call to initracer)
56       IF (.NOT.ALLOCATED(noms)) ALLOCATE(noms(nq))
57       ALLOCATE(mmol(nq))
58       ALLOCATE(radius(nq))
59       ALLOCATE(rho_q(nq))
60       ALLOCATE(qext(nq))
61       ALLOCATE(alpha_lift(nq))
62       ALLOCATE(alpha_devil(nq))
63       ALLOCATE(qextrhor(nq))
64       ALLOCATE(igcm_dustbin(nq))
65       !! initialization
66       alpha_lift(:)=0.
67       alpha_devil(:)=0.
68
69! Initialization: get tracer names from the dynamics and check if we are
70!                 using 'old' tracer convention ('q01',q02',...)
71!                 or new convention (full tracer names)
72      ! check if tracers have 'old' names
73
74! copy tracer names from dynamics
75        do iq=1,nq
76          noms(iq)=nametrac(iq)
77        enddo
78
79
80! Identify tracers by their names: (and set corresponding values of mmol)
81      ! 0. initialize tracer indexes to zero:
82      ! NB: igcm_* indexes are commons in 'tracer.h'
83      do iq=1,nq
84        igcm_dustbin(iq)=0
85      enddo
86      igcm_dust_mass=0
87      igcm_dust_number=0
88      igcm_h2o_vap=0
89      igcm_h2o_ice=0
90      igcm_co2=0
91      igcm_co=0
92      igcm_o=0
93      igcm_o1d=0
94      igcm_o2=0
95      igcm_o3=0
96      igcm_h=0
97      igcm_h2=0
98      igcm_oh=0
99      igcm_ho2=0
100      igcm_h2o2=0
101      igcm_n2=0
102      igcm_ar=0
103      igcm_ar_n2=0
104      igcm_co2_ice=0
105
106      write(*,*) 'initracer: noms() ', noms
107
108
109      !print*,'Setting dustbin = 0 in initracer.F'
110      !dustbin=0
111
112      ! 1. find dust tracers
113      count=0
114!      if (dustbin.gt.0) then
115!        do iq=1,nq
116!          txt=" "
117!          write(txt,'(a4,i2.2)')'dust',count+1   
118!          if (noms(iq).eq.txt) then
119!            count=count+1
120!            igcm_dustbin(count)=iq
121!            mmol(iq)=100.
122!          endif
123!        enddo !do iq=1,nq
124!      endif ! of if (dustbin.gt.0)
125
126
127!      if (doubleq) then
128!        do iq=1,nq
129!          if (noms(iq).eq."dust_mass") then
130!            igcm_dust_mass=iq
131!            count=count+1
132!          endif
133!          if (noms(iq).eq."dust_number") then
134!            igcm_dust_number=iq
135!            count=count+1
136!          endif
137!        enddo
138!      endif ! of if (doubleq)
139      ! 2. find chemistry and water tracers
140      do iq=1,nq
141        if (noms(iq).eq."co2") then
142          igcm_co2=iq
143          mmol(igcm_co2)=44.
144          count=count+1
145!          write(*,*) 'co2: count=',count
146        endif
147        if (noms(iq).eq."co2_ice") then
148          igcm_co2_ice=iq
149          mmol(igcm_co2_ice)=44.
150          count=count+1
151!          write(*,*) 'co2_ice: count=',count
152        endif
153        if (noms(iq).eq."h2o_vap") then
154          igcm_h2o_vap=iq
155          mmol(igcm_h2o_vap)=18.
156          count=count+1
157!          write(*,*) 'h2o_vap: count=',count
158        endif
159        if (noms(iq).eq."h2o_ice") then
160          igcm_h2o_ice=iq
161          mmol(igcm_h2o_ice)=18.
162          count=count+1
163!          write(*,*) 'h2o_ice: count=',count
164        endif
165      enddo ! of do iq=1,nq
166     
167      ! check that we identified all tracers:
168      if (count.ne.nq) then
169        write(*,*) "initracer: found only ",count," tracers"
170        write(*,*) "               expected ",nq
171        do iq=1,count
172          write(*,*)'      ',iq,' ',trim(noms(iq))
173        enddo
174!        stop
175      else
176        write(*,*) "initracer: found all expected tracers, namely:"
177        do iq=1,nq
178          write(*,*)'      ',iq,' ',trim(noms(iq))
179        enddo
180      endif
181
182
183c------------------------------------------------------------
184c     Initialisation tracers ....
185c------------------------------------------------------------
186      call zerophys(nq,rho_q)
187
188      rho_dust=2500.  ! Mars dust density (kg.m-3)
189      rho_ice=920.    ! Water ice density (kg.m-3)
190      rho_co2=1620.   ! CO2 ice density (kg.m-3)
191
192
193
194c$$$      if (doubleq) then
195c$$$c       "doubleq" technique 
196c$$$c       -------------------
197c$$$c      (transport of mass and number mixing ratio)
198c$$$c       iq=1: Q mass mixing ratio, iq=2: N number mixing ratio
199c$$$
200c$$$        if( (nq.lt.2).or.(water.and.(nq.lt.3)) ) then
201c$$$          write(*,*)'initracer: nq is too low : nq=', nq
202c$$$          write(*,*)'water= ',water,' doubleq= ',doubleq   
203c$$$        end if
204c$$$
205c$$$        varian=0.637    ! Characteristic variance   
206c$$$        qext(igcm_dust_mass)=3.04   ! reference extinction at 0.67 um for ref dust
207c$$$        qext(igcm_dust_number)=3.04 ! reference extinction at 0.67 um for ref dust
208c$$$        rho_q(igcm_dust_mass)=rho_dust
209c$$$        rho_q(igcm_dust_number)=rho_dust
210c$$$
211c$$$c       Intermediate calcul for computing geometric mean radius r0
212c$$$c       as a function of mass and number mixing ratio Q and N
213c$$$c       (r0 = (r3n_q * Q/ N)
214c$$$        r3n_q = exp(-4.5*varian**2)*(3./4.)/(pi*rho_dust)
215c$$$
216c$$$c       Intermediate calcul for computing effective radius reff
217c$$$c       from geometric mean radius r0
218c$$$c       (reff = ref_r0 * r0)
219c$$$        ref_r0 = exp(2.5*varian**2)
220c$$$       
221c$$$c       lifted dust :
222c$$$c       '''''''''''
223c$$$        reff_lift = 3.e-6      !  Effective radius of lifted dust (m)
224c$$$        alpha_devil(igcm_dust_mass)=9.e-9   !  dust devil lift mass coeff
225c$$$        alpha_lift(igcm_dust_mass)=3.0e-15  !  Lifted mass coeff
226c$$$
227c$$$        r0_lift = reff_lift/ref_r0
228c$$$        alpha_devil(igcm_dust_number)=r3n_q*
229c$$$     &                        alpha_devil(igcm_dust_mass)/r0_lift**3
230c$$$        alpha_lift(igcm_dust_number)=r3n_q*
231c$$$     &                        alpha_lift(igcm_dust_mass)/r0_lift**3
232c$$$
233c$$$c       Not used:
234c$$$        radius(igcm_dust_mass) = 0.
235c$$$        radius(igcm_dust_number) = 0.
236c$$$
237c$$$      else
238
239
240c$$$       if (dustbin.gt.1) then
241c$$$        print*,'ATTENTION:',
242c$$$     $   ' properties of dust need input in initracer !!!'
243c$$$        stop
244c$$$
245c$$$       else if (dustbin.eq.1) then
246c$$$
247c$$$c       This will be used for 1 dust particle size:
248c$$$c       ------------------------------------------
249c$$$        radius(igcm_dustbin(1))=3.e-6
250c$$$        Qext(igcm_dustbin(1))=3.04
251c$$$        alpha_lift(igcm_dustbin(1))=0.0e-6
252c$$$        alpha_devil(igcm_dustbin(1))=7.65e-9
253c$$$        qextrhor(igcm_dustbin(1))=(3./4.)*Qext(igcm_dustbin(1))
254c$$$     &                         /(rho_dust*radius(igcm_dustbin(1)))
255c$$$        rho_q(igcm_dustbin(1))=rho_dust
256c$$$
257c$$$       endif
258c$$$!      end if    ! (doubleq)
259
260c     Initialization for water vapor
261c     ------------------------------
262      if(water) then
263         radius(igcm_h2o_vap)=0.
264         Qext(igcm_h2o_vap)=0.
265         alpha_lift(igcm_h2o_vap) =0.
266         alpha_devil(igcm_h2o_vap)=0.
267         qextrhor(igcm_h2o_vap)= 0.
268
269c       "Dryness coefficient" controlling the evaporation and
270c        sublimation from the ground water ice (close to 1)
271c        HERE, the goal is to correct for the fact
272c        that the simulated permanent water ice polar caps
273c        is larger than the actual cap and the atmospheric
274c        opacity not always realistic.
275
276
277!         if(ngrid.eq.1)
278
279
280!     to be modified for BC+ version?
281
282         !! this is defined in surfdat_h.F90
283         IF (.not.ALLOCATED(dryness)) ALLOCATE(dryness(ngrid))
284         IF (.not.ALLOCATED(watercaptag)) ALLOCATE(watercaptag(ngrid))
285
286         do ig=1,ngrid
287           if (ngrid.ne.1) watercaptag(ig)=.false.
288           dryness(ig) = 1.
289         enddo
290
291
292
293
294!         IF (caps) THEN
295c Perennial H20 north cap defined by watercaptag=true (allows surface to be
296c hollowed by sublimation in vdifc).
297!         do ig=1,ngrid
298!           if (lati(ig)*180./pi.gt.84) then
299!             if (ngrid.ne.1) watercaptag(ig)=.true.
300!             dryness(ig) = 1.
301c Use the following cap definition for high spatial resolution (latitudinal bin <= 5 deg)
302c             if (lati(ig)*180./pi.lt.85.and.long(ig).ge.0) then
303c               if (ngrid.ne.1) watercaptag(ig)=.true.
304c               dryness(ig) = 1.
305c             endif
306c             if (lati(ig)*180./pi.ge.85) then
307c               if (ngrid.ne.1) watercaptag(ig)=.true.
308c               dryness(ig) = 1.
309c             endif
310!           endif  ! (lati>80 deg)
311!         end do ! (ngrid)
312!        ENDIF ! (caps)
313
314!         if(iceparty.and.(nq.ge.2)) then
315
316           radius(igcm_h2o_ice)=3.e-6
317           rho_q(igcm_h2o_ice)=rho_ice
318           Qext(igcm_h2o_ice)=0.
319!           alpha_lift(igcm_h2o_ice) =0.
320!           alpha_devil(igcm_h2o_ice)=0.
321           qextrhor(igcm_h2o_ice)= (3./4.)*Qext(igcm_h2o_ice) 
322     $       / (rho_ice*radius(igcm_h2o_ice))
323
324
325
326!         elseif(iceparty.and.(nq.lt.2)) then
327!            write(*,*) 'nq is too low : nq=', nq
328!            write(*,*) 'water= ',water,' iceparty= ',iceparty
329!         endif
330
331      end if  ! (water)
332
333c     Output for records:
334c     ~~~~~~~~~~~~~~~~~~
335      write(*,*)
336      Write(*,*) '******** initracer : dust transport parameters :'
337      write(*,*) 'alpha_lift = ', alpha_lift
338      write(*,*) 'alpha_devil = ', alpha_devil
339      write(*,*) 'radius  = ', radius
340!      if(doubleq) then
341!        write(*,*) 'reff_lift (um) =  ', reff_lift
342!        write(*,*) 'size distribution variance  = ', varian
343!        write(*,*) 'r3n_q , ref_r0 : ', r3n_q , ref_r0
344!      end if
345      write(*,*) 'Qext  = ', qext 
346      write(*,*)
347
348      end
Note: See TracBrowser for help on using the repository browser.