source: trunk/mars/libf/phymars/initracer.F @ 77

Last change on this file since 77 was 38, checked in by emillour, 14 years ago

Ajout du modè Martien (mon LMDZ.MARS.BETA, du 28/01/2011) dans le rértoire mars, pour pouvoir suivre plus facilement les modifs.
EM

File size: 22.1 KB
Line 
1      SUBROUTINE initracer(qsurf,co2ice)
2
3       IMPLICIT NONE
4c=======================================================================
5c   subject:
6c   --------
7c   Initialization related to tracer
8c   (transported dust, water, chemical species, ice...)
9c
10c   Name of the tracer
11c
12c   Test of dimension :
13c   Initialize COMMON tracer in tracer.h, using tracer names provided
14c   by the dynamics in "advtrac.h"
15c
16c   Old conventions: (not used any more)
17c
18c   If water=T : q(iq=nqmx) is the water mass mixing ratio
19c     and q(iq=nqmx-1) is the ice mass mixing ratio
20
21c   If there is transported dust, it uses iq=1 to iq=dustbin
22c   If there is no transported dust : dustbin=0
23c   If doubleq=T : q(iq=1) is the dust mass mixing ratio
24c                  q(iq=2) is the dust number mixing ratio
25
26c   If (photochem.or.thermochem) there is "ncomp" chemical species (ncomp
27c   is set in aeronomars/chimiedata.h) using the ncomp iq values starting at
28c      iq=nqchem_min = dustbin+1   (nqchem_min is defined in inifis.F)
29c   
30c
31c   author: F.Forget
32c   ------
33c    Modifs: Franck Montmessin, Sebastien Lebonnois (june 2003)
34c            Ehouarn Millour (oct. 2008) identify tracers by their names
35c=======================================================================
36
37
38#include "dimensions.h"
39#include "dimphys.h"
40#include "comcstfi.h"
41#include "callkeys.h"
42#include "tracer.h"
43#include "advtrac.h"
44#include "comgeomfi.h"
45#include "watercap.h"
46#include "chimiedata.h"
47
48
49      real qsurf(ngridmx,nqmx)       ! tracer on surface (e.g.  kg.m-2)
50      real co2ice(ngridmx)           ! co2 ice mass on surface (e.g.  kg.m-2)
51      integer iq,ig,count
52      real r0_lift , reff_lift, nueff_lift
53c     Ratio of small over large dust particles (used when both
54c       doubleq and the submicron mode are active); In Montmessin
55c       et al. (2002), a value of 25 has been deduced;
56      real, parameter :: popratio = 25.
57      logical :: oldnames ! =.true. if old tracer naming convention (q01,...)
58      character(len=20) :: txt ! to store some text
59
60c-----------------------------------------------------------------------
61c  radius(nqmx)      ! aerosol particle radius (m)
62c  rho_q(nqmx)       ! tracer densities (kg.m-3)
63c  alpha_lift(nqmx)  ! saltation vertical flux/horiz flux ratio (m-1)
64c  alpha_devil(nqmx) ! lifting coeeficient by dust devil
65c  rho_dust          ! Mars dust density
66c  rho_ice           ! Water ice density
67c  nuice_ref         ! Effective variance nueff of the
68c                    !   water-ice size distributions
69c  doubleq           ! if method with mass (iq=1) and number(iq=2) mixing ratio
70c  varian            ! Characteristic variance of log-normal distribution
71c-----------------------------------------------------------------------
72
73      integer :: nqchem_start
74
75! Initialization: get tracer names from the dynamics and check if we are
76!                 using 'old' tracer convention ('q01',q02',...)
77!                 or new convention (full tracer names)
78      ! check if tracers have 'old' names
79      oldnames=.false.
80      count=0
81      do iq=1,nqmx
82        txt=" "
83        write(txt,'(a1,i2.2)') 'q',iq
84        if (txt.eq.tnom(iq)) then
85          count=count+1
86        endif
87      enddo ! of do iq=1,nqmx
88     
89      if (count.eq.nqmx) then
90        write(*,*) "initracer: tracers seem to follow old naming ",
91     &             "convention (q01,q02,...)"
92        write(*,*) "   => will work for now ... "
93        write(*,*) "      but you should run newstart to rename them"
94        oldnames=.true.
95      endif
96
97      ! copy/set tracer names
98      if (oldnames) then ! old names (derived from iq & option values)
99        ! 1. dust:
100        if (dustbin.ne.0) then ! transported dust
101          do iq=1,dustbin
102            txt=" "
103            write(txt,'(a4,i2.2)') 'dust',iq
104            noms(iq)=txt
105            mmol(iq)=100.
106          enddo
107          ! special case if doubleq
108          if (doubleq) then
109            if (dustbin.ne.2) then
110              write(*,*) 'initracer: error expected dustbin=2'
111            else
112              noms(1)='dust_mass'   ! dust mass mixing ratio
113              noms(2)='dust_number' ! dust number mixing ratio
114            endif
115          endif
116        endif
117        ! 2. water & ice
118        if (water) then
119          noms(nqmx)='h2o_vap'
120          mmol(nqmx)=18.
121          noms(nqmx-1)='h2o_ice'
122          mmol(nqmx-1)=18.
123        endif
124        ! 3. Chemistry
125        if (photochem .or. callthermos) then
126          if (doubleq) then
127            nqchem_start=3
128          else
129            nqchem_start=dustbin+1
130          end if
131          do iq=nqchem_start,nqchem_start+ncomp-1
132            noms(iq)=nomchem(iq-nqchem_start+1)
133            mmol(iq)=mmolchem(iq-nqchem_start+1)
134          enddo
135        endif ! of if (photochem .or. callthermos)
136        ! 4. Other tracers ????
137        if ((dustbin.eq.0).and.(.not.water)) then
138          noms(1)='co2'
139          mmol(1)=44
140          if (nqmx.eq.2) then
141            noms(nqmx)='Ar_N2'
142            mmol(nqmx)=30
143          endif
144        endif
145      else
146        ! copy tracer names from dynamics
147        do iq=1,nqmx
148          noms(iq)=tnom(iq)
149        enddo
150        ! mmol(:) array is initialized later (see below)
151      endif ! of if (oldnames)
152
153      ! special modification when using 'old' tracers:
154      ! qsurf(nqmx) was h2o ice whereas q(nqmx) was water vapour
155      ! and (if iceparty) q(nqmx-1) was null whereas q(nqmx-1) was water ice
156      if (oldnames.and.water) then
157        write(*,*)'initracer: moving surface water ice to index ',nqmx-1
158        qsurf(1:ngridmx,nqmx-1)=qsurf(1:ngridmx,nqmx)
159        qsurf(1:ngridmx,nqmx)=0
160      endif
161     
162c------------------------------------------------------------
163c     Test Dimensions tracers
164c------------------------------------------------------------
165
166! Ehouarn: testing number of tracers is obsolete...
167!      if(photochem.or.thermochem) then
168!          if (water) then
169!              if ((dustbin+ncomp+2).ne.nqmx) then
170!                 print*,'initracer: tracer dimension problem:'
171!                 print*,"(dustbin+ncomp+2).ne.nqmx"
172!                 print*,"ncomp: ",ncomp
173!                 print*,"dustbin: ",dustbin
174!                 print*,"nqmx: ",nqmx
175!                 print*,'Change ncomp in chimiedata.h'
176!               endif
177!          else
178!              if ((dustbin+ncomp+1).ne.nqmx) then
179!                 print*,'initracer: tracer dimension problem:'
180!                 print*,"(dustbin+ncomp+1).ne.nqmx"
181!                 print*,"ncomp: ",ncomp
182!                 print*,"dustbin: ",dustbin
183!                 print*,"nqmx: ",nqmx
184!                 print*,'Change ncomp in chimiedata.h'
185!                 STOP
186!               endif
187!            endif
188!      endif
189
190c------------------------------------------------------------
191c         NAME and molar mass of the tracer
192c------------------------------------------------------------
193
194   
195! Identify tracers by their names: (and set corresponding values of mmol)
196      ! 0. initialize tracer indexes to zero:
197      do iq=1,nqmx
198        igcm_dustbin(iq)=0
199      enddo
200      igcm_dust_mass=0
201      igcm_dust_number=0
202      igcm_dust_submicron=0
203      igcm_h2o_vap=0
204      igcm_h2o_ice=0
205      igcm_co2=0
206      igcm_co=0
207      igcm_o=0
208      igcm_o1d=0
209      igcm_o2=0
210      igcm_o3=0
211      igcm_h=0
212      igcm_h2=0
213      igcm_oh=0
214      igcm_ho2=0
215      igcm_h2o2=0
216      igcm_n2=0
217      igcm_ar=0
218      igcm_ar_n2=0
219
220      ! 1. find dust tracers
221      count=0
222      if (dustbin.gt.0) then
223        do iq=1,nqmx
224          txt=" "
225          write(txt,'(a4,i2.2)')'dust',count+1
226          if (noms(iq).eq.txt) then
227            count=count+1
228            igcm_dustbin(count)=iq
229            mmol(iq)=100.
230          endif
231        enddo !do iq=1,nqmx
232      endif ! of if (dustbin.gt.0)
233      if (doubleq) then
234        do iq=1,nqmx
235          if (noms(iq).eq."dust_mass") then
236            igcm_dust_mass=iq
237            count=count+1
238          endif
239          if (noms(iq).eq."dust_number") then
240            igcm_dust_number=iq
241            count=count+1
242          endif
243        enddo
244      endif ! of if (doubleq)
245      if (submicron) then
246        do iq=1,nqmx
247          if (noms(iq).eq."dust_submicron") then
248            igcm_dust_submicron=iq
249            mmol(iq)=100.
250            count=count+1
251          endif
252        enddo
253      endif ! of if (submicron)
254      ! 2. find chemistry and water tracers
255      do iq=1,nqmx
256        if (noms(iq).eq."co2") then
257          igcm_co2=iq
258          mmol(igcm_co2)=44.
259          count=count+1
260        endif
261        if (noms(iq).eq."co") then
262          igcm_co=iq
263          mmol(igcm_co)=28.
264          count=count+1
265        endif
266        if (noms(iq).eq."o") then
267          igcm_o=iq
268          mmol(igcm_o)=16.
269          count=count+1
270        endif
271        if (noms(iq).eq."o1d") then
272          igcm_o1d=iq
273          mmol(igcm_o1d)=16.
274          count=count+1
275        endif
276        if (noms(iq).eq."o2") then
277          igcm_o2=iq
278          mmol(igcm_o2)=32.
279          count=count+1
280        endif
281        if (noms(iq).eq."o3") then
282          igcm_o3=iq
283          mmol(igcm_o3)=48.
284          count=count+1
285        endif
286        if (noms(iq).eq."h") then
287          igcm_h=iq
288          mmol(igcm_h)=1.
289          count=count+1
290        endif
291        if (noms(iq).eq."h2") then
292          igcm_h2=iq
293          mmol(igcm_h2)=2.
294          count=count+1
295        endif
296        if (noms(iq).eq."oh") then
297          igcm_oh=iq
298          mmol(igcm_oh)=17.
299          count=count+1
300        endif
301        if (noms(iq).eq."ho2") then
302          igcm_ho2=iq
303          mmol(igcm_ho2)=33.
304          count=count+1
305        endif
306        if (noms(iq).eq."h2o2") then
307          igcm_h2o2=iq
308          mmol(igcm_h2o2)=34.
309          count=count+1
310        endif
311        if (noms(iq).eq."n2") then
312          igcm_n2=iq
313          mmol(igcm_n2)=28.
314          count=count+1
315        endif
316        if (noms(iq).eq."ar") then
317          igcm_ar=iq
318          mmol(igcm_ar)=40.
319          count=count+1
320        endif
321        if (noms(iq).eq."h2o_vap") then
322          igcm_h2o_vap=iq
323          mmol(igcm_h2o_vap)=18.
324          count=count+1
325        endif
326        if (noms(iq).eq."h2o_ice") then
327          igcm_h2o_ice=iq
328          mmol(igcm_h2o_ice)=18.
329          count=count+1
330        endif
331        ! Other stuff: e.g. for simulations using co2 + neutral gaz
332        if (noms(iq).eq."Ar_N2") then
333          igcm_ar_n2=iq
334          mmol(igcm_ar_n2)=30.
335          count=count+1
336        endif
337      enddo ! of do iq=1,nqmx
338!      count=count+nbqchem
339     
340      ! check that we identified all tracers:
341      if (count.ne.nqmx) then
342        write(*,*) "initracer: found only ",count," tracers"
343        write(*,*) "               expected ",nqmx
344        do iq=1,count
345          write(*,*)'      ',iq,' ',trim(noms(iq))
346        enddo
347        stop
348      else
349        write(*,*) "initracer: found all expected tracers, namely:"
350        do iq=1,nqmx
351          write(*,*)'      ',iq,' ',trim(noms(iq))
352        enddo
353      endif
354
355      ! if water cycle but iceparty=.false., there will nevertheless be
356      ! water ice at the surface (iceparty is not used anymore, but this
357      ! part is still relevant, as we want to stay compatible with the
358      ! older versions).
359      if (water.and.(igcm_h2o_ice.eq.0)) then
360        igcm_h2o_ice=igcm_h2o_vap ! so that qsurf(i_h2o_ice) is identified
361                                  ! even though there is no q(i_h2o_ice)
362      else
363       ! surface ice qsurf(i_h2o_ice) was loaded twice by phyetat0,
364       ! as qsurf(i_h2o_vap) & as qsurf(i_h2o_ice), so to be clean:
365       if (igcm_h2o_vap.ne.0) then
366         qsurf(1:ngridmx,igcm_h2o_vap)=0
367       endif
368      endif
369
370c------------------------------------------------------------
371c     Initialisation tracers ....
372c------------------------------------------------------------
373      call zerophys(nqmx,rho_q)
374
375      rho_dust=2500.  ! Mars dust density (kg.m-3)
376      rho_ice=920.    ! Water ice density (kg.m-3)
377      nuice_ref=0.1   ! Effective variance nueff of the
378                      ! water-ice size distributions
379
380      if (doubleq) then
381c       "doubleq" technique
382c       -------------------
383c      (transport of mass and number mixing ratio)
384c       iq=1: Q mass mixing ratio, iq=2: N number mixing ratio
385
386        if( (nqmx.lt.2).or.(water.and.(nqmx.lt.4)) ) then
387          write(*,*)'initracer: nqmx is too low : nqmx=', nqmx
388          write(*,*)'water= ',water,' doubleq= ',doubleq   
389        end if
390
391        nueff_lift = 0.5
392        varian=sqrt(log(1.+nueff_lift))
393
394        rho_q(igcm_dust_mass)=rho_dust
395        rho_q(igcm_dust_number)=rho_dust
396
397c       Intermediate calcul for computing geometric mean radius r0
398c       as a function of mass and number mixing ratio Q and N
399c       (r0 = (r3n_q * Q/ N)^(1/3))
400        r3n_q = exp(-4.5*varian**2)*(3./4.)/(pi*rho_dust)
401
402c       Intermediate calcul for computing effective radius reff
403c       from geometric mean radius r0
404c       (reff = ref_r0 * r0)
405        ref_r0 = exp(2.5*varian**2)
406       
407c       lifted dust :
408c       '''''''''''
409        reff_lift = 3.0e-6 !3.e-6 !Effective radius of lifted dust (m)
410        alpha_devil(igcm_dust_mass)=9.e-9   !  dust devil lift mass coeff
411c       alpha_lift(igcm_dust_mass)=3.0e-15  !  Lifted mass coeff
412        alpha_lift(igcm_dust_mass)=1.e-6 !1.e-6 !Lifted mass coeff
413
414        r0_lift = reff_lift/ref_r0
415        alpha_devil(igcm_dust_number)=r3n_q*
416     &                        alpha_devil(igcm_dust_mass)/r0_lift**3
417        alpha_lift(igcm_dust_number)=r3n_q*
418     &                        alpha_lift(igcm_dust_mass)/r0_lift**3
419
420        radius(igcm_dust_mass) = reff_lift
421        radius(igcm_dust_number) = reff_lift
422
423        write(*,*) "initracer: doubleq_param reff_lift:", reff_lift
424        write(*,*) "initracer: doubleq_param nueff_lift:", nueff_lift
425        write(*,*) "initracer: doubleq_param alpha_lift:",
426     &    alpha_lift(igcm_dust_mass)
427
428      else
429
430       if (dustbin.gt.1) then
431        print*,'initracer: STOP!',
432     $   ' properties of dust need to be set in initracer !!!'
433        stop
434
435       else if (dustbin.eq.1) then
436
437c       This will be used for 1 dust particle size:
438c       ------------------------------------------
439        radius(igcm_dustbin(1))=3.e-6
440        alpha_lift(igcm_dustbin(1))=0.0e-6
441        alpha_devil(igcm_dustbin(1))=7.65e-9
442        rho_q(igcm_dustbin(1))=rho_dust
443
444       endif
445      end if    ! (doubleq)
446
447c     Submicron dust mode:
448c     --------------------
449
450      if (submicron) then
451        radius(igcm_dust_submicron)=0.1e-6
452        rho_q(igcm_dust_submicron)=rho_dust
453        if (doubleq) then
454c         If doubleq is also active, we use the population ratio:
455          alpha_lift(igcm_dust_submicron) =
456     &      alpha_lift(igcm_dust_number)*popratio*
457     &      rho_q(igcm_dust_submicron)*4./3.*pi*
458     &      radius(igcm_dust_submicron)**3.
459          alpha_devil(igcm_dust_submicron)=1.e-30
460        else
461          alpha_lift(igcm_dust_submicron)=1e-6
462          alpha_devil(igcm_dust_submicron)=1.e-30
463        endif ! (doubleq)
464      end if  ! (submicron)
465
466c     Initialization for photochemistry:
467c     ---------------------------------
468      if (photochem) then
469      ! initialize chemistry+water (water will be correctly initialized below)
470      ! by initializing everything which is not dust ...
471        do iq=1,nqmx
472          txt=noms(iq)
473          if (txt(1:4).ne."dust") then
474            radius(iq)=0.
475            alpha_lift(iq) =0.
476            alpha_devil(iq)=0.
477          endif
478        enddo ! do iq=1,nqmx
479      endif
480
481c     Initialization for water vapor
482c     ------------------------------
483      if(water) then
484         radius(igcm_h2o_vap)=0.
485         alpha_lift(igcm_h2o_vap) =0.
486         alpha_devil(igcm_h2o_vap)=0.
487
488c       "Dryness coefficient" controlling the evaporation and
489c        sublimation from the ground water ice (close to 1)
490c        HERE, the goal is to correct for the fact
491c        that the simulated permanent water ice polar caps
492c        is larger than the actual cap and the atmospheric
493c        opacity not always realistic.
494
495         do ig=1,ngridmx
496           if (ngridmx.ne.1) watercaptag(ig)=.false.
497           dryness(ig) = 1.
498         enddo
499
500         IF (caps) THEN
501c Perennial H20 north cap defined by watercaptag=true (allows surface to be
502c hollowed by sublimation in vdifc).
503         do ig=1,ngridmx
504           if (lati(ig)*180./pi.gt.84) then
505             if (ngridmx.ne.1) watercaptag(ig)=.true.
506             dryness(ig) = 1.
507c Use the following cap definition for high spatial resolution (latitudinal bin <= 5 deg)
508c             if (lati(ig)*180./pi.lt.85.and.long(ig).ge.0) then
509c               if (ngridmx.ne.1) watercaptag(ig)=.true.
510c               dryness(ig) = 1.
511c             endif
512c             if (lati(ig)*180./pi.ge.85) then
513c               if (ngridmx.ne.1) watercaptag(ig)=.true.
514c               dryness(ig) = 1.
515c             endif
516           endif  ! (lati>80 deg)
517         end do ! (ngridmx)
518        ENDIF ! (caps)
519
520         if(water.and.(nqmx.ge.2)) then
521           radius(igcm_h2o_ice)=3.e-6
522           rho_q(igcm_h2o_ice)=rho_ice
523           alpha_lift(igcm_h2o_ice) =0.
524           alpha_devil(igcm_h2o_ice)=0.
525         elseif(water.and.(nqmx.lt.2)) then
526            write(*,*) 'nqmx is too low : nqmx=', nqmx
527            write(*,*) 'water= ',water
528         endif
529
530      end if  ! (water)
531
532c     Output for records:
533c     ~~~~~~~~~~~~~~~~~~
534      write(*,*)
535      Write(*,*) '******** initracer : dust transport parameters :'
536      write(*,*) 'alpha_lift = ', alpha_lift
537      write(*,*) 'alpha_devil = ', alpha_devil
538      write(*,*) 'radius  = ', radius
539      if(doubleq) then
540        write(*,*) 'reff_lift (um) =  ', reff_lift
541        write(*,*) 'size distribution variance  = ', varian
542        write(*,*) 'r3n_q , ref_r0 : ', r3n_q , ref_r0
543      end if
544
545!
546!     some extra (possibly redundant) sanity checks for tracers:
547!     ---------------------------------------------------------
548
549       if (doubleq) then
550       ! verify that we indeed have dust_mass and dust_number tracers
551         if (igcm_dust_mass.eq.0) then
552           write(*,*) "initracer: error !!"
553           write(*,*) "  cannot use doubleq option without ",
554     &                "a dust_mass tracer !"
555           stop
556         endif
557         if (igcm_dust_number.eq.0) then
558           write(*,*) "initracer: error !!"
559           write(*,*) "  cannot use doubleq option without ",
560     &                "a dust_number tracer !"
561           stop
562         endif
563       endif
564
565       if ((.not.doubleq).and.(dustbin.gt.0)) then
566       ! verify that we indeed have 'dustbin' dust tracers
567         count=0
568         do iq=1,dustbin
569           if (igcm_dustbin(iq).ne.0) then
570             count=count+1
571           endif
572         enddo
573         if (count.ne.dustbin) then
574           write(*,*) "initracer: error !!"
575           write(*,*) "  dusbin is set to ",dustbin,
576     &                " but we only have the following dust tracers:"
577           do iq=1,count
578             write(*,*)"   ",trim(noms(igcm_dustbin(iq)))
579           enddo
580           stop
581         endif
582       endif
583
584       if (water) then
585       ! verify that we indeed have h2o_vap and h2o_ice tracers
586         if (igcm_h2o_vap.eq.0) then
587           write(*,*) "initracer: error !!"
588           write(*,*) "  cannot use water option without ",
589     &                "an h2o_vap tracer !"
590           stop
591         endif
592         if (igcm_h2o_ice.eq.0) then
593           write(*,*) "initracer: error !!"
594           write(*,*) "  cannot use water option without ",
595     &                "an h2o_ice tracer !"
596           stop
597         endif
598       endif
599
600       if (photochem .or. callthermos) then
601       ! verify that we indeed have the chemistry tracers
602         if (igcm_co2.eq.0) then
603           write(*,*) "initracer: error !!"
604           write(*,*) "  cannot use chemistry option without ",
605     &                "a co2 tracer !"
606         stop
607         endif
608         if (igcm_co.eq.0) then
609           write(*,*) "initracer: error !!"
610           write(*,*) "  cannot use chemistry option without ",
611     &                "a co tracer !"
612         stop
613         endif
614         if (igcm_o.eq.0) then
615           write(*,*) "initracer: error !!"
616           write(*,*) "  cannot use chemistry option without ",
617     &                "a o tracer !"
618         stop
619         endif
620         if (igcm_o1d.eq.0) then
621           write(*,*) "initracer: error !!"
622           write(*,*) "  cannot use chemistry option without ",
623     &                "a o1d tracer !"
624         stop
625         endif
626         if (igcm_o2.eq.0) then
627           write(*,*) "initracer: error !!"
628           write(*,*) "  cannot use chemistry option without ",
629     &                "an o2 tracer !"
630         stop
631         endif
632         if (igcm_o3.eq.0) then
633           write(*,*) "initracer: error !!"
634           write(*,*) "  cannot use chemistry option without ",
635     &                "an o3 tracer !"
636         stop
637         endif
638         if (igcm_h.eq.0) then
639           write(*,*) "initracer: error !!"
640           write(*,*) "  cannot use chemistry option without ",
641     &                "a h tracer !"
642         stop
643         endif
644         if (igcm_h2.eq.0) then
645           write(*,*) "initracer: error !!"
646           write(*,*) "  cannot use chemistry option without ",
647     &                "a h2 tracer !"
648         stop
649         endif
650         if (igcm_oh.eq.0) then
651           write(*,*) "initracer: error !!"
652           write(*,*) "  cannot use chemistry option without ",
653     &                "an oh tracer !"
654         stop
655         endif
656         if (igcm_ho2.eq.0) then
657           write(*,*) "initracer: error !!"
658           write(*,*) "  cannot use chemistry option without ",
659     &                "a ho2 tracer !"
660         stop
661         endif
662         if (igcm_h2o2.eq.0) then
663           write(*,*) "initracer: error !!"
664           write(*,*) "  cannot use chemistry option without ",
665     &                "a h2o2 tracer !"
666         stop
667         endif
668         if (igcm_n2.eq.0) then
669           write(*,*) "initracer: error !!"
670           write(*,*) "  cannot use chemistry option without ",
671     &                "a n2 tracer !"
672         stop
673         endif
674         if (igcm_ar.eq.0) then
675           write(*,*) "initracer: error !!"
676           write(*,*) "  cannot use chemistry option without ",
677     &                "an ar tracer !"
678         stop
679         endif
680       endif ! of if (photochem .or. callthermos)
681
682      end
Note: See TracBrowser for help on using the repository browser.