source: trunk/LMDZ.GENERIC/libf/phystd/inifis.F @ 164

Last change on this file since 164 was 135, checked in by aslmd, 14 years ago

CHANGEMENT ARBORESCENCE ETAPE 2 -- NON COMPLET
File size: 26.1 KB
RevLine 
[135]1      SUBROUTINE inifis(ngrid,nlayer,
2     $           day_ini,pdaysec,ptimestep,
3     $           plat,plon,parea,
4     $           prad,pg,pr,pcpp)
5!
6!=======================================================================
7!
8!   purpose:
9!   -------
10!
11!   Initialisation for the physical parametrisations of the LMD
12!   martian atmospheric general circulation modele.
13!
14!   author: Frederic Hourdin 15 / 10 /93
15!   -------
16!   modified: Sebastien Lebonnois 11/06/2003 (new callphys.def)
17!             Ehouarn Millour (oct. 2008) tracers are now identified
18!              by their names and may not be contiguously
19!              stored in the q(:,:,:,:) array
20!             E.M. (june 2009) use getin routine to load parameters
21!
22!
23!   arguments:
24!   ----------
25!
26!   input:
27!   ------
28!
29!    ngrid                 Size of the horizontal grid.
30!                          All internal loops are performed on that grid.
31!    nlayer                Number of vertical layers.
32!    pdayref               Day of reference for the simulation
33!    pday                  Number of days counted from the North. Spring
34!                          equinoxe.
35!
36!=======================================================================
37!
38!-----------------------------------------------------------------------
39!   declarations:
40!   -------------
41! to use  'getin'
42      USE ioipsl_getincom
43      IMPLICIT NONE
44#include "dimensions.h"
45#include "dimphys.h"
46#include "planete.h"
47#include "comcstfi.h"
48#include "comsaison.h"
49#include "comdiurn.h"
50#include "comgeomfi.h"
51#include "callkeys.h"
52#include "surfdat.h"
53
54
55      REAL prad,pg,pr,pcpp,pdaysec,ptimestep
56 
57      INTEGER ngrid,nlayer
58      REAL plat(ngrid),plon(ngrid),parea(ngridmx)
59      integer day_ini
60      INTEGER ig,ierr
61 
62      EXTERNAL iniorbit,orbite
63      EXTERNAL SSUM
64      REAL SSUM
65 
66      CHARACTER ch1*12
67      CHARACTER ch80*80
68
69      logical chem, h2o
70      logical :: parameter, doubleq=.false.
71
72      rad=prad
73      daysec=pdaysec
74      dtphys=ptimestep
75      cpp=pcpp
76      g=pg
77      r=pr
78      rcp=r/cpp
79
80
81      avocado = 6.02214179e23   ! added by RW
82
83
84! --------------------------------------------------------
85!     The usual Tests
86!     --------------
87      IF (nlayer.NE.nlayermx) THEN
88         PRINT*,'STOP in inifis'
89         PRINT*,'Probleme de dimensions :'
90         PRINT*,'nlayer     = ',nlayer
91         PRINT*,'nlayermx   = ',nlayermx
92         STOP
93      ENDIF
94
95      IF (ngrid.NE.ngridmx) THEN
96         PRINT*,'STOP in inifis'
97         PRINT*,'Probleme de dimensions :'
98         PRINT*,'ngrid     = ',ngrid
99         PRINT*,'ngridmx   = ',ngridmx
100         STOP
101      ENDIF
102
103! --------------------------------------------------------------
104!  Reading the "callphys.def" file controlling some key options
105! --------------------------------------------------------------
106     
107      ! check that 'callphys.def' file is around
108      OPEN(99,file='callphys.def',status='old',form='formatted'
109     &     ,iostat=ierr)
110      CLOSE(99)
111     
112      IF(ierr.EQ.0) THEN
113         PRINT*
114         PRINT*
115         PRINT*,'--------------------------------------------'
116         PRINT*,' inifis: Parametres pour la physique (callphys.def)'
117         PRINT*,'--------------------------------------------'
118
119
120         write(*,*) "Run with or without tracer transport ?"
121         tracer=.false. ! default value
122         call getin("tracer",tracer)
123         write(*,*) " tracer = ",tracer
124
125         write(*,*) "Diurnal cycle ?"
126         write(*,*) "(if diurnal=false, diurnal averaged solar heating)"
127         diurnal=.true. ! default value
128         call getin("diurnal",diurnal)
129         write(*,*) " diurnal = ",diurnal
130
131         write(*,*) "Seasonal cycle ?"
132         write(*,*) "(if season=false, Ls stays constant, to value ",
133     &   "set in 'start'"
134         season=.true. ! default value
135         call getin("season",season)
136         write(*,*) " season = ",season
137
138         write(*,*) "Tidally resonant rotation ?"
139         tlocked=.false. ! default value
140         call getin("tlocked",tlocked)
141         write(*,*) "tlocked = ",tlocked
142
143! Test of incompatibility:
144! if tlocked, then diurnal should be false
145         if (tlocked.and.diurnal) then
146           print*,'If diurnal=true, we should turn off tlocked.'
147           stop
148         endif
149
150         write(*,*) "Tidal resonance ratio ?"
151         nres=0          ! default value
152         call getin("nres",nres)
153         write(*,*) "nres = ",nres
154
155         write(*,*) "Write some extra output to the screen ?"
156         lwrite=.false. ! default value
157         call getin("lwrite",lwrite)
158         write(*,*) " lwrite = ",lwrite
159
160         write(*,*) "Save statistics in file stats.nc ?"
161         callstats=.true. ! default value
162         call getin("callstats",callstats)
163         write(*,*) " callstats = ",callstats
164
165         write(*,*) "Test energy conservation of model physics ?"
166         enertest=.false. ! default value
167         call getin("enertest",enertest)
168         write(*,*) " enertest = ",enertest
169
170         write(*,*) "Save EOF profiles in file 'profiles' for ",
171     &              "Climate Database?"
172         calleofdump=.false. ! default value
173         call getin("calleofdump",calleofdump)
174         write(*,*) " calleofdump = ",calleofdump
175
176         write(*,*) "Dust scenario:"
177         iaervar=3 ! default value
178         call getin("iaervar",iaervar)
179         write(*,*) " iaervar = ",iaervar
180
181         write(*,*) "Dust vertical distribution:"
182         write(*,*) "(=1 Dust opt.deph read in startfi;",
183     & " =2 Viking scenario; =3 MGS scenario,",
184     & " =4 Mars Year 24 from TES assimilation)"
185         iddist=3 ! default value
186         call getin("iddist",iddist)
187         write(*,*) " iddist = ",iddist
188
189         write(*,*) "Dust top altitude (km). (Matters only if iddist=1)"
190         topdustref= 90.0 ! default value
191         call getin("topdustref",topdustref)
192         write(*,*) " topdustref = ",topdustref
193
194
195         write(*,*) "call radiative transfer ?"
196         callrad=.true. ! default value
197         call getin("callrad",callrad)
198         write(*,*) " callrad = ",callrad
199
200         write(*,*) "call correlated-k radiative transfer ?"
201         corrk=.true. ! default value
202         call getin("corrk",corrk)
203         write(*,*) " corrk = ",corrk
204
205!         write(*,*) "call NLTE radiative schemes ?",
206!     &              "(matters only if callrad=T)"
207!         callnlte=.false. ! default value
208!         call getin("callnlte",callnlte)
209!         write(*,*) " callnlte = ",callnlte
210         
211         write(*,*) "call gaseous absorption in the visible bands?",
212     &              "(matters only if callrad=T)"
213         callgasvis=.false. ! default value
214         call getin("callgasvis",callgasvis)
215         write(*,*) " callgasvis = ",callgasvis
216         
217         write(*,*) "call turbulent vertical diffusion ?"
218         calldifv=.true. ! default value
219         call getin("calldifv",calldifv)
220         write(*,*) " calldifv = ",calldifv
221
222         write(*,*) "call convective adjustment ?"
223         calladj=.true. ! default value
224         call getin("calladj",calladj)
225         write(*,*) " calladj = ",calladj
226       
227         write(*,*)"Gas is nonideal CO2 ?"
228         nonideal=.false.
229         call getin("nonideal",nonideal)
230         write(*,*)" nonideal = ",nonideal
231
232! Test of incompatibility
233         if (enertest.and.nonideal) then
234            print*,'Energy conservation calculations currently
235     &           assume ideal gas!'
236            call abort
237         endif
238
239         write(*,*) "call CO2 condensation ?"
240         co2cond=.true. ! default value
241         call getin("co2cond",co2cond)
242         write(*,*) " co2cond = ",co2cond
243
244         write(*,*)"call thermal conduction in the soil ?"
245         callsoil=.true. ! default value
246         call getin("callsoil",callsoil)
247         write(*,*) " callsoil = ",callsoil
248         
249!         write(*,*)"call Lott's gravity wave/subgrid topography ",
250!     &             "scheme ?"
251!         calllott=.true. ! default value
252!         call getin("calllott",calllott)
253!         write(*,*)" calllott = ",calllott
254
255         write(*,*)"rad.transfer is computed every iradia",
256     &             " physical timestep"
257         iradia=1 ! default value
258         call getin("iradia",iradia)
259         write(*,*)" iradia = ",iradia
260       
261         write(*,*)"Rayleigh scattering ?"
262         rayleigh=.false.
263         call getin("rayleigh",rayleigh)
264         write(*,*)" rayleigh = ",rayleigh
265
266         write(*,*)"Parametrized Earth-like ozone absorption ?"
267         ozone=.false.     
268         call getin("ozone",ozone)
269         write(*,*) " ozone = ",ozone
270
271         write(*,*)"Output mean OLR in 1D?"
272         meanOLR=.false.
273         call getin("meanOLR",meanOLR)
274         write(*,*)" meanOLR = ",meanOLR
275
276         write(*,*)"Output spectral OLR in 3D?"
277         specOLR=.false.
278         call getin("specOLR",specOLR)
279         write(*,*)" specOLR = ",specOLR
280
281         write(*,*)"Default planetary temperature?"
282         tplanet=215.0
283         call getin("tplanet",tplanet)
284         write(*,*)" tplanet = ",tplanet
285
286         write(*,*)"Which star?"
287         startype=1 ! default value = Sol
288         call getin("startype",startype)
289         write(*,*)" startype = ",startype
290
291         write(*,*)"Value of stellar flux at 1 AU?"
292         Fat1AU=1356.0 ! default value = Sol today
293         call getin("Fat1AU",Fat1AU)
294         write(*,*)" Fat1AU = ",Fat1AU
295
296         write(*,*)"Set temperature to 1 K above CO2 condensation?"
297         nearco2cond=.false.
298         call getin("nearco2cond",nearco2cond)
299         write(*,*)" nearco2cond = ",nearco2cond
300
301! TRACERS:
302
303         write(*,*)"Transported dust ? (if >0, use 'dustbin' dust bins)"
304         dustbin=0 ! default value
305         call getin("dustbin",dustbin)
306         write(*,*)" dustbin = ",dustbin
307
308         write(*,*)"Radiatively active aerosols?"
309         aerofixed=.true. ! default value
310         call getin("aerofixed",aerofixed)
311         write(*,*)" aerofixed = ",aerofixed
312
313         write(*,*)"Number mixing ratio of CO2 ice particles:"
314         Nmix_co2=100000. ! default value
315         call getin("Nmix_co2",Nmix_co2)
316         write(*,*)" Nmix_co2 = ",Nmix_co2
317
318         write(*,*)"Number mixing ratio of H2O ice particles:"
319         Nmix_h2o=10000000. ! default value
320         call getin("Nmix_h2o",Nmix_h2o)
321         write(*,*)" Nmix_h2o = ",Nmix_h2o
322
323         write(*,*)"Is the variable gas species radiatively active?"
324         varactive=.false.
325         call getin("varactive",varactive)
326         write(*,*)" varactive = ",varactive
327
328         write(*,*)"Is the variable gas species distribution set?"
329         varfixed=.false.
330         call getin("varfixed",varfixed)
331         write(*,*)" varfixed = ",varfixed
332
333         write(*,*)"What is the saturation % of the variable species?"
334         satval=0.8
335         call getin("satval",satval)
336         write(*,*)" satval = ",satval
337
338!         write(*,*)"Radiatively active dust ? (matters if dustbin>0)"
339!         active=.false. ! default value
340!         call getin("active",active)
341!         write(*,*)" active = ",active
342
343! Test of incompatibility:
344! if no tracers, then aerofixed should be true
345         if ((.not.tracer).and.(.not.aerofixed)) then
346           print*,'if tracers are off, aerofixed must be ON!'
347           stop
348         endif
349
350! Test of incompatibility:
351! if varactive, then varfixed should be false
352         if (varactive.and.varfixed) then
353           print*,'if varactive, varfixed must be OFF!'
354           stop
355         endif
356
357! Test of incompatibility:
358! if active is used, then dustbin should be > 0
359
360!         if (active.and.(dustbin.lt.1)) then
361!           print*,'if active is used, then dustbin should > 0'
362!           stop
363!         endif
364
365!         write(*,*)"use mass and number mixing ratios to predict",
366!     &             " dust size ?"
367!         doubleq=.false. ! default value
368!         call getin("doubleq",doubleq)
369!         write(*,*)" doubleq = ",doubleq
370
371! Test of incompatibility:
372! if doubleq is used, then dustbin should be 1
373
374!         if (doubleq.and.(dustbin.ne.1)) then
375!           print*,'if doubleq is used, then dustbin should be 1'
376!           stop
377!         endif
378
379!         write(*,*)"dust lifted by GCM surface winds ?"
380!         lifting=.false. ! default value
381!         call getin("lifting",lifting)
382!         write(*,*)" lifting = ",lifting
383
384! Test of incompatibility:
385! if lifting is used, then dustbin should be > 0
386
387!         if (lifting.and.(dustbin.lt.1)) then
388!           print*,'if lifting is used, then dustbin should > 0'
389!           stop
390!         endif
391
392!         write(*,*)" dust lifted by dust devils ?"
393!         callddevil=.false. !default value
394!         call getin("callddevil",callddevil)
395!         write(*,*)" callddevil = ",callddevil
396         
397
398! Test of incompatibility:
399! if dustdevil is used, then dustbin should be > 0
400
401!         if (callddevil.and.(dustbin.lt.1)) then
402!           print*,'if dustdevil is used, then dustbin should > 0'
403!           stop
404!         endif
405
406!         write(*,*)"Dust scavenging by CO2 snowfall ?"
407!         scavenging=.false. ! default value
408!         call getin("scavenging",scavenging)
409!         write(*,*)" scavenging = ",scavenging
410         
411
412! Test of incompatibility:
413! if scavenging is used, then dustbin should be > 0
414
415!         if (scavenging.and.(dustbin.lt.1)) then
416!           print*,'if scavenging is used, then dustbin should > 0'
417!           stop
418!         endif
419
420         write(*,*) "Gravitationnal sedimentation ?"
421         sedimentation=.true. ! default value
422         call getin("sedimentation",sedimentation)
423         write(*,*) " sedimentation = ",sedimentation
424
425!         write(*,*) "includes water ice",
426!     &              "(if true, 'water' must also be .true.)"
427!         iceparty=.false. ! default value
428!         call getin("iceparty",iceparty)
429!         write(*,*) " iceparty = ",iceparty
430
431!         write(*,*) "Radiatively active transported atmospheric ",
432!     &              "water ice ?"
433!         activice=.false. ! default value
434!         call getin("activice",activice)
435!         write(*,*) " activice = ",activice
436
437
438! Test of incompatibility:
439! if activice is used, then iceparty should be used too
440
441!         if (activice.and..not.iceparty) then
442!           print*,'if activice is used, iceparty should be used too'
443!           stop
444!         endif
445
446         write(*,*) "Compute water cycle ?"
447         water=.false. ! default value
448         call getin("water",water)
449         write(*,*) " water = ",water
450         
451         write(*,*) "Include water condensation ?"
452         watercond=.false. ! default value
453         call getin("watercond",watercond)
454         write(*,*) " watercond = ",watercond
455
456         write(*,*) "Include water precipitation ?"
457         waterrain=.false. ! default value
458         call getin("waterrain",waterrain)
459         write(*,*) " waterrain = ",waterrain
460
461         write(*,*) "Precipitation threshold ?"
462         rainthreshold=0.011 ! default value (Emmanuel 1997)
463         call getin("rainthreshold",rainthreshold)
464         write(*,*) " rainthreshold = ",rainthreshold
465
466! Test of incompatibility:
467! if watercond is used, then water should be used too
468
469         if (watercond.and.(.not.watercond)) then
470           print*,'if watercond is used, water should be used too'
471           stop
472         endif
473
474         write(*,*) "photochemistry: include chemical species"
475         photochem=.false. ! default value
476         call getin("photochem",photochem)
477         write(*,*) " photochem = ",photochem
478
479
480! THERMOSPHERE
481
482         write(*,*) "call thermosphere ?"
483         callthermos=.false. ! default value
484         call getin("callthermos",callthermos)
485         write(*,*) " callthermos = ",callthermos
486         
487         write(*,*) " water included without cycle ",
488     &              "(only if water=.false.)"
489         thermoswater=.false. ! default value
490         call getin("thermoswater",thermoswater)
491         write(*,*) " thermoswater = ",thermoswater
492
493         write(*,*) "call thermal conduction ?",
494     &    " (only if callthermos=.true.)"
495         callconduct=.false. ! default value
496         call getin("callconduct",callconduct)
497         write(*,*) " callconduct = ",callconduct
498
499         write(*,*) "call EUV heating ?",
500     &   " (only if callthermos=.true.)"
501         calleuv=.false.  ! default value
502         call getin("calleuv",calleuv)
503         write(*,*) " calleuv = ",calleuv
504
505         write(*,*) "call molecular viscosity ?",
506     &   " (only if callthermos=.true.)"
507         callmolvis=.false. ! default value
508         call getin("callmolvis",callmolvis)
509         write(*,*) " callmolvis = ",callmolvis
510
511         write(*,*) "call molecular diffusion ?",
512     &   " (only if callthermos=.true.)"
513         callmoldiff=.false. ! default value
514         call getin("callmoldiff",callmoldiff)
515         write(*,*) " callmoldiff = ",callmoldiff
516         
517         write(*,*) "call thermospheric photochemistry ?",
518     &   " (only if callthermos=.true.)"
519         thermochem=.false. ! default value
520         call getin("thermochem",thermochem)
521         write(*,*) " thermochem = ",thermochem
522
523         write(*,*) "date for solar flux calculation:",
524     &   " (1985 < date < 2002)"
525         write(*,*) "(Solar min=1996.4 ave=1993.4 max=1990.6)"
526         solarcondate=1993.4 ! default value
527         call getin("solarcondate",solarcondate)
528         write(*,*) " solarcondate = ",solarcondate
529         
530
531         if (.not.callthermos) then
532           if (thermoswater) then
533             print*,'if thermoswater is set, callthermos must be true'
534             stop
535           endif         
536           if (callconduct) then
537             print*,'if callconduct is set, callthermos must be true'
538             stop
539           endif       
540           if (calleuv) then
541             print*,'if calleuv is set, callthermos must be true'
542             stop
543           endif         
544           if (callmolvis) then
545             print*,'if callmolvis is set, callthermos must be true'
546             stop
547           endif       
548           if (callmoldiff) then
549             print*,'if callmoldiff is set, callthermos must be true'
550             stop
551           endif         
552           if (thermochem) then
553             print*,'if thermochem is set, callthermos must be true'
554             stop
555           endif         
556        endif
557
558! Test of incompatibility:
559! if photochem is used, then water should be used too
560
561         if (photochem.and..not.water) then
562           print*,'if photochem is used, water should be used too'
563           stop
564         endif
565
566! if callthermos is used, then thermoswater should be used too
567! (if water not used already)
568
569         if (callthermos .and. .not.water) then
570           if (callthermos .and. .not.thermoswater) then
571             print*,'if callthermos is used, water or thermoswater
572     &               should be used too'
573             stop
574           endif
575         endif
576
577         PRINT*,'--------------------------------------------'
578         PRINT*
579         PRINT*
580      ELSE
581         write(*,*)
582         write(*,*) 'Cannot read file callphys.def. Is it here ?'
583         stop
584      ENDIF
585
5868000  FORMAT(t5,a12,l8)
5878001  FORMAT(t5,a12,i8)
588
589      PRINT*
590      PRINT*,'inifis: daysec',daysec
591      PRINT*
592      PRINT*,'inifis: The radiative transfer is computed:'
593      PRINT*,'           each ',iradia,' physical time-step'
594      PRINT*,'        or each ',iradia*dtphys,' seconds'
595      PRINT*
596! --------------------------------------------------------------
597!  Managing the Longwave radiative transfer
598! --------------------------------------------------------------
599
600!     In most cases, the run just use the following values :
601!     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
602      callemis=.true.     
603!     ilwd=10*int(daysec/dtphys) ! bug before 22/10/01       
604      ilwd=10*int(daysec/dtphys)
605      ilwn=2               
606      linear=.true.       
607      ncouche=3
608      alphan=0.4
609      ilwb=2
610      semi=0
611
612c$$$!     BUT people working hard on the LW may want to read them in 'radia.def'
613c$$$!     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
614c$$$
615c$$$      OPEN(99,file='radia.def',status='old',form='formatted'
616c$$$     .     ,iostat=ierr)
617c$$$      IF(ierr.EQ.0) THEN
618c$$$         write(*,*) 'inifis: Reading radia.def !!!'
619c$$$         READ(99,fmt='(a)') ch1
620c$$$         READ(99,*) callemis
621c$$$         WRITE(*,8000) ch1,callemis
622c$$$
623c$$$         READ(99,fmt='(a)') ch1
624c$$$         READ(99,*) iradia
625c$$$         WRITE(*,8001) ch1,iradia
626c$$$
627c$$$         READ(99,fmt='(a)') ch1
628c$$$         READ(99,*) ilwd
629c$$$         WRITE(*,8001) ch1,ilwd
630c$$$
631c$$$         READ(99,fmt='(a)') ch1
632c$$$         READ(99,*) ilwn
633c$$$         WRITE(*,8001) ch1,ilwn
634c$$$
635c$$$         READ(99,fmt='(a)') ch1
636c$$$         READ(99,*) linear
637c$$$         WRITE(*,8000) ch1,linear
638c$$$
639c$$$         READ(99,fmt='(a)') ch1
640c$$$         READ(99,*) ncouche
641c$$$         WRITE(*,8001) ch1,ncouche
642c$$$
643c$$$         READ(99,fmt='(a)') ch1
644c$$$         READ(99,*) alphan
645c$$$         WRITE(*,*) ch1,alphan
646c$$$
647c$$$         READ(99,fmt='(a)') ch1
648c$$$         READ(99,*) ilwb
649c$$$         WRITE(*,8001) ch1,ilwb
650c$$$
651c$$$
652c$$$         READ(99,fmt='(a)') ch1
653c$$$         READ(99,'(l1)') callg2d
654c$$$         WRITE(*,8000) ch1,callg2d
655c$$$
656c$$$         READ(99,fmt='(a)') ch1
657c$$$         READ(99,*) semi
658c$$$         WRITE(*,*) ch1,semi
659c$$$      end if
660c$$$      CLOSE(99)
661
662!-----------------------------------------------------------------------
663!     Some more initialization:
664!     ------------------------
665
666      CALL SCOPY(ngrid,plon,1,long,1)
667      CALL SCOPY(ngrid,plat,1,lati,1)
668      CALL SCOPY(ngrid,parea,1,area,1)
669      totarea=SSUM(ngridmx,area,1)
670
671      DO ig=1,ngrid
672         sinlat(ig)=sin(plat(ig))
673         coslat(ig)=cos(plat(ig))
674         sinlon(ig)=sin(plon(ig))
675         coslon(ig)=cos(plon(ig))
676      ENDDO
677
678      pi=2.*asin(1.) ! NB: pi is a common in comcstfi.h
679
680!     managing the tracers, and tests:
681!     -------------------------------
682
683      if(tracer) then
684
685!          when photochem is used, nqchem_min is the rank
686!          of the first chemical species
687
688! Ehouarn: nqchem_min is now meaningless and no longer used
689!       nqchem_min = 1
690       if (photochem .or. callthermos) then
691         chem = .true.
692       end if
693
694       if (water .or. thermoswater)then
695          h2o = .true.
696       else
697          h2o = .false.
698       endif
699
700!          TESTS
701
702       print*,'inifis: TRACERS:'
703
704c$$$       if ((doubleq).and.(h2o).and.
705c$$$     $     (chem).and.(iceparty)) then
706c$$$         print*,' 2 dust tracers (doubleq)'
707c$$$         print*,' 1 water vapour tracer'
708c$$$         print*,' 1 water ice tracer'
709c$$$         print*,nqmx-4,' chemistry tracers'
710c$$$       endif
711c$$$
712c$$$       if ((doubleq).and.(h2o).and.
713c$$$     $     (chem).and..not.(iceparty)) then
714c$$$         print*,' 2 dust tracers (doubleq)'
715c$$$         print*,' 1 water vapour tracer'
716c$$$         print*,nqmx-3,' chemistry tracers'
717c$$$       endif
718c$$$
719c$$$       if ((doubleq).and.(h2o).and.
720c$$$     $     .not.(chem).and.(iceparty)) then
721c$$$         print*,' 2 dust tracers (doubleq)'
722c$$$         print*,' 1 water vapour tracer'
723c$$$         print*,' 1 water ice tracer'
724c$$$         if (nqmx.ne.4) then
725c$$$           print*,'nqmx should be 4 with these options.'
726c$$$               print*,'(or check callphys.def)'
727c$$$           stop
728c$$$         endif
729c$$$       endif
730c$$$
731c$$$       if ((doubleq).and.(h2o).and.
732c$$$     $     .not.(chem).and..not.(iceparty)) then
733c$$$         print*,' 2 dust tracers (doubleq)'
734c$$$         print*,' 1 water vapour tracer'
735c$$$         if (nqmx.ne.3) then
736c$$$           print*,'nqmx should be 3 with these options...'
737c$$$               print*,'(or check callphys.def)'
738c$$$           stop
739c$$$         endif
740c$$$       endif
741c$$$
742c$$$       if ((doubleq).and..not.(h2o)) then
743c$$$         print*,' 2 dust tracers (doubleq)'
744c$$$         if (nqmx.ne.2) then
745c$$$           print*,'nqmx should be 2 with these options...'
746c$$$               print*,'(or check callphys.def)'
747c$$$           stop
748c$$$         endif
749c$$$       endif
750
751       if (.not.(doubleq).and.(h2o).and.
752!     $     (chem).and.(iceparty)) then
753     $     (chem).and.(watercond)) then
754         if (dustbin.gt.0) then
755           print*,dustbin,' dust bins'
756         endif
757         print*,nqmx-2-dustbin,' chemistry tracers'
758         print*,' 1 water vapour tracer'
759         print*,' 1 water ice tracer'
760       endif
761
762       if (.not.(doubleq).and.(h2o).and.
763!     $     (chem).and..not.(iceparty)) then
764     $     (chem).and..not.(watercond)) then
765         if (dustbin.gt.0) then
766           print*,dustbin,' dust bins'
767         endif
768         print*,nqmx-1-dustbin,' chemistry tracers'
769         print*,' 1 water vapour tracer'
770       endif
771
772       if (.not.(doubleq).and.(h2o).and.
773!     $     .not.(chem).and.(iceparty)) then
774     $     .not.(chem).and.(watercond)) then
775         if (dustbin.gt.0) then
776           print*,dustbin,' dust bins'
777         endif
778         print*,' 1 water vapour tracer'
779         print*,' 1 water ice tracer'
780
781         if (nqmx.gt.(dustbin+2)) then
782           print*,'nqmx should be ',(dustbin+2),
783     $            ' with these options...'
784                   print*,'(or check callphys.def)'
785         endif
786         if (nqmx.lt.(dustbin+2)) then
787            print*,dustbin,' dust bins, but this should be ok I think.'
788!           stop
789         endif
790       endif
791
792       if (.not.(doubleq).and.(h2o).and.
793!     $     .not.(chem).and..not.(iceparty)) then
794     $     .not.(chem).and..not.(watercond)) then
795         if (dustbin.gt.0) then
796           print*,dustbin,' dust bins'
797         endif
798         print*,' 1 water vapour tracer'
799         if (nqmx.gt.(dustbin+1)) then
800           print*,'nqmx should be ',(dustbin+1),
801     $            ' with these options...'
802                   print*,'(or check callphys.def)'
803         if (nqmx.lt.(dustbin+1)) then
804           stop
805         endif
806         endif
807       endif
808
809!       if (.not.(doubleq).and..not.(h2o)) then
810!         if (dustbin.gt.0) then
811!           print*,dustbin,' dust bins'
812!           if (nqmx.ne.dustbin) then
813!             print*,'nqmx should be ',dustbin,
814!     $              ' with these options...'
815!             print*,'(or check callphys.def)'
816!             stop
817!           endif
818!         else
819           print*,'dustbin=',dustbin,
820     $            ': tracer should be F with these options...'
821     $           ,'UNLESS you just want to move tracers around '
822!         endif
823!       endif
824
825      endif ! of if (tracer)
826
827      RETURN
828      END
Note: See TracBrowser for help on using the repository browser.