source: trunk/LMDZ.MARS/libf/aeronomars/inichim_newstart.F @ 503

Last change on this file since 503 was 460, checked in by emillour, 14 years ago

Mars GCM: more updates for photochemistry:
-aeronomars/inichim_newstart.F : add methane + bug correction in "oldnames" case
-aeronomars/read_phototable.F : adapt routine so that photolysis table file

name may be set in the def file with
phototable = filename
(default is "jmars.20111014")

EM

File size: 19.6 KB
RevLine 
[38]1      subroutine inichim_newstart(pq,ps,sig,nq,latfi,lonfi,airefi,
2     $ thermo,qsurf)
3
4      implicit none
5
6c=======================================================================
7c
8c   subject:
9c   --------
10c
11c  Initialization of the composition for use in the building of a new start file
12c  used by program newstart.F
13c  also used by program testphys1d.F
14c
15c   Author:   Sebastien Lebonnois (08/11/2002)
16c   -------
[460]17c   Revised 07/2003 by Monica Angelats-i-Coll to use input files
18c   Modified 10/2008 Identify tracers by their names Ehouarn Millour
19c   Modified 11/2011 Addition of methane (Franck Lefevre)
[38]20c
21c   Arguments:
22c   ----------
23c
24c    pq(iip1,jjp1,llm,nqmx)  Advected fields, ie chemical species here
25c    sig = sigma                vertical coordinate (interface of layers)
26c    nq: number of tracers to initialize (used to evaluate if only
27c        chemistry species are to be initialized or if water vapour
28c        and water ice should also be initialized.
29c        NB: number of chemistry tracers is ncomp (set in chimiedata.h)
30c=======================================================================
31
32c    Declarations :
33c    --------------
34
35#include "dimensions.h"
36#include "dimphys.h"
37#include "paramet.h"
38#include "chimiedata.h"
39#include "tracer.h"
40#include "comcstfi.h"
41#include "comdiurn.h"
42#include "callkeys.h"
43#include "temps.h"
44#include "datafile.h"
45
46c Arguments :
47c -----------
48
49      real    ps(iip1,jjp1)
50      real    pq(iip1,jjp1,llm,nqmx)     ! Advected fields, ie chemical species
51      real    sig(llm+1)                 ! vertical coordinate (interface of layers)
52      integer nq            !  =nqmx
53                            ! or =nqmx-1 if H2O kept
54                            ! or =nqmx-2 if H2O and ice kept
55      REAL latfi(ngridmx),lonfi(ngridmx),airefi(ngridmx)
56      integer thermo ! flag for thermosphere initialisation only
57      real :: qsurf(ngridmx,nqmx) ! surface tracers
58
59c Local variables :
60c -----------------
61
62      integer  iq,i,j,l,n, nspecini,inil,iqmax,iiq
63      INTEGER aa(nqmx)
64      REAL qtot(iip1,jjp1,llm)
65      real densconc,zgcm,zfile(252),vmrint,nt, mmean
66      real nxf(252),zfilemin(252),sigfile(252)
67      real vmrf(252,14),nf(252,14)
68      real tfile(252),zzfile(252)
69      real ni(14),niold(14)
70      real mmolf(14)    ! molecular mass in amu (species in files)
71      data mmolf/44.,40.,28.,32.,28.,16.,2.,   ! majors
72     &           1.,17.,33.,18.,34.,16.,48./   ! minors
73      character*3 tmp ! temporary variable
74      integer ierr,lnblnk
75      external lnblnk
76
77      logical :: oldnames ! =.true. if old tracer naming convention (q01,...)
78      character(len=20) :: txt ! to store some text
79      integer :: nqchem_start,count
80      integer :: nqchem(nqmx) ! local chemistry tracer index array
81      integer :: nbqchem ! total number of chemical species (water included)
82
83c-----------------------------------------------------------------------
84c    Input/Output
85c    ------------
86! read in 'callphys.def'
87       call inichim_readcallphys()
88
89      ! check if tracers have 'old' names
90      oldnames=.false.
91      count=0
92      do iq=1,nqmx
93        txt=" "
94        write(txt,'(a1,i2.2)') 'q',iq
95        if (txt.eq.noms(iq)) then
96          count=count+1
97        endif
98      enddo ! of do iq=1,nqmx
99
100      if (count.eq.nqmx) then
101        write(*,*) "inichim_newstart: tracers seem to follow old ",
102     &             "naming convention (q01,q02,...)"
103        write(*,*) "   => will work for now ... "
104        write(*,*) "      but you should rename them"
105        oldnames=.true.
106      endif
107
108      ! copy/set tracer names
109      if (oldnames) then ! old names (derived from iq & option values)
110        ! 1. dust:
111        if (dustbin.ne.0) then ! transported dust
112          do iq=1,dustbin
113            txt=" "
114            write(txt,'(a4,i2.2)') 'dust',iq
115            noms(iq)=txt
116            mmol(iq)=100.
117          enddo
118          ! special case if doubleq
119          if (doubleq) then
120            if (dustbin.ne.2) then
121              write(*,*) 'initracer: error expected dustbin=2'
122            else
[164]123!              noms(1)='dust_mass'   ! dust mass mixing ratio
124!              noms(2)='dust_number' ! dust number mixing ratio
125! same as above, but avoid explict possible out-of-bounds on array
126               noms(1)='dust_mass'   ! dust mass mixing ratio
127               do iq=2,2
128               noms(iq)='dust_number' ! dust number mixing ratio
129               enddo
[38]130            endif
131          endif
132        endif
133        ! 2. water & ice
134        if (water) then
135          noms(nqmx)='h2o_vap'
136          mmol(nqmx)=18.
[164]137!            noms(nqmx-1)='h2o_ice'
138!            mmol(nqmx-1)=18.
139          !"loop" to avoid potential out-of-bounds in array
140          do iq=nqmx-1,nqmx-1
141            noms(iq)='h2o_ice'
142            mmol(iq)=18.
143          enddo
[38]144        endif
145        ! 3. Chemistry
146        if (photochem .or. callthermos) then
147          if (doubleq) then
148            nqchem_start=3
149          else
150            nqchem_start=dustbin+1
151          end if
152        endif ! of if (photochem .or. callthermos)
[460]153        do iq=nqchem_start,nqchem_start+ncomp-2
[38]154          noms(iq)=nomchem(iq-nqchem_start+1)
155          mmol(iq)=mmolchem(iq-nqchem_start+1)
156        enddo
157        ! 4. Other tracers ????
158        if ((dustbin.eq.0).and.(.not.water)) then
159          noms(1)='co2'
160          mmol(1)=44
161          if (nqmx.eq.2) then
162            noms(nqmx)='Ar_N2'
163            mmol(nqmx)=30
164          endif
165        endif
166      else
167        ! noms(:) already contain tracer names
168        ! mmol(:) array is initialized later (see below)
169      endif ! of if (oldnames)
170
171      ! special modification when using 'old' tracers:
172      ! qsurf(nqmx) was h2o ice whereas q(nqmx) was water vapour
173      ! and (if iceparty) q(nqmx-1) was null whereas q(nqmx-1) was water ice
174      if (oldnames.and.water) then
175        write(*,*)'inichim_newstart: moving surface water ice to index '
176     &            ,nqmx-1
[164]177        ! "loop" to avoid potential out-of-bounds on arrays
178        do iq=nqmx-1,nqmx-1
179        qsurf(1:ngridmx,iq)=qsurf(1:ngridmx,iq+1)
180        enddo
[38]181        qsurf(1:ngridmx,nqmx)=0
182      endif
183
184c Dimension test:
185c ---------------
186
187       if (water) then
188!         if ((nqchem_min+ncomp+1).ne.nqmx) then
189         if ((dustbin+ncomp+2).ne.nqmx) then
190!            print*,'********* Dimension problem! ********'
191!            print*,"nqchem_min+ncomp+1).ne.nqmx"
192            print*,'inichim_newstart: tracer dimension problem:'
193            print*,"(dustbin+ncomp+2).ne.nqmx"
194            print*,"ncomp: ",ncomp
195!            print*,"nqchem_min: ",nqchem_min
196            print*,"nqmx: ",nqmx
197            print*,'Change ncomp in chimiedata.h'
198            STOP
199         endif
200       else
201!         if ((nqchem_min+ncomp).ne.nqmx) then
202          if ((dustbin+ncomp+1).ne.nqmx) then
203!            print*,'********* Dimension problem! ********'
204!            print*,"nqchem_min+ncomp).ne.nqmx"
205            print*,'initracer: tracer dimension problem:'
206            print*,"(dustbin+ncomp+1).ne.nqmx"
207            print*,"ncomp: ",ncomp
208!            print*,"nqchem_min: ",nqchem_min
209            print*,"nqmx: ",nqmx
210            print*,'Change ncomp in chimiedata.h'
211            STOP
212         endif
213       endif
214
215c noms and mmol vectors:
216c ----------------------
217!
218!       if (iceparty) then
219!          do iq=nqchem_min,nqmx-2
220!            noms(iq) = nomchem(iq-nqchem_min+1)
221!            mmol(iq) = mmolchem(iq-nqchem_min+1)
222!          enddo
223!          noms(nqmx-1) = 'ice'
224!          mmol(nqmx-1) = 18.
225!          noms(nqmx)   = 'h2o'
226!          mmol(nqmx)   = 18.
227!       else
228!          do iq=nqchem_min,nqmx-1
229!            noms(iq) = nomchem(iq-nqchem_min+1)
230!            mmol(iq) = mmolchem(iq-nqchem_min+1)
231!          enddo
232!          noms(nqmx)   = 'h2o'
233!          mmol(nqmx)   = 18.
234!       endif
235!       if (nqchem_min.gt.1) then
236!          do iq=1,nqchem_min-1
237!            noms(iq) = 'dust'
238!            mmol(iq) = 100.
239!          enddo
240!       endif
241
242
243! Identify tracers by their names: (and set corresponding values of mmol)
244      ! 0. initialize tracer indexes to zero:
245      do iq=1,nqmx
246        igcm_dustbin(iq)=0
247      enddo
248      igcm_dust_mass=0
249      igcm_dust_number=0
250      igcm_h2o_vap=0
251      igcm_h2o_ice=0
252      igcm_co2=0
253      igcm_co=0
254      igcm_o=0
255      igcm_o1d=0
256      igcm_o2=0
257      igcm_o3=0
258      igcm_h=0
259      igcm_h2=0
260      igcm_oh=0
261      igcm_ho2=0
262      igcm_h2o2=0
[460]263      igcm_ch4=0
[38]264      igcm_n2=0
265      igcm_ar=0
266      igcm_ar_n2=0
267
268      ! 1. find dust tracers
269      count=0
270      if (dustbin.gt.0) then
271        do iq=1,nqmx
272          txt=" "
273          write(txt,'(a4,i2.2)')'dust',count+1
274          if (noms(iq).eq.txt) then
275            count=count+1
276            igcm_dustbin(count)=iq
277            mmol(iq)=100.
278          endif
279        enddo !do iq=1,nqmx
280      endif ! of if (dustbin.gt.0)
281      if (doubleq) then
282        do iq=1,nqmx
283          if (noms(iq).eq."dust_mass") then
284            igcm_dust_mass=iq
285            count=count+1
286          endif
287          if (noms(iq).eq."dust_number") then
288            igcm_dust_number=iq
289            count=count+1
290          endif
291        enddo
292      endif ! of if (doubleq)
293      ! 2. find chemistry and water tracers
294      nbqchem=0
295      do iq=1,nqmx
296        if (noms(iq).eq."co2") then
297          igcm_co2=iq
298          mmol(igcm_co2)=44.
299          count=count+1
300          nbqchem=nbqchem+1
301        endif
302        if (noms(iq).eq."co") then
303          igcm_co=iq
304          mmol(igcm_co)=28.
305          count=count+1
306          nbqchem=nbqchem+1
307        endif
308        if (noms(iq).eq."o") then
309          igcm_o=iq
310          mmol(igcm_o)=16.
311          count=count+1
312          nbqchem=nbqchem+1
313        endif
314        if (noms(iq).eq."o1d") then
315          igcm_o1d=iq
316          mmol(igcm_o1d)=16.
317          count=count+1
318          nbqchem=nbqchem+1
319        endif
320        if (noms(iq).eq."o2") then
321          igcm_o2=iq
322          mmol(igcm_o2)=32.
323          count=count+1
324          nbqchem=nbqchem+1
325        endif
326        if (noms(iq).eq."o3") then
327          igcm_o3=iq
328          mmol(igcm_o3)=48.
329          count=count+1
330          nbqchem=nbqchem+1
331        endif
332        if (noms(iq).eq."h") then
333          igcm_h=iq
334          mmol(igcm_h)=1.
335          count=count+1
336          nbqchem=nbqchem+1
337        endif
338        if (noms(iq).eq."h2") then
339          igcm_h2=iq
340          mmol(igcm_h2)=2.
341          count=count+1
342          nbqchem=nbqchem+1
343        endif
344        if (noms(iq).eq."oh") then
345          igcm_oh=iq
346          mmol(igcm_oh)=17.
347          count=count+1
348          nbqchem=nbqchem+1
349        endif
350        if (noms(iq).eq."ho2") then
351          igcm_ho2=iq
352          mmol(igcm_ho2)=33.
353          count=count+1
354          nbqchem=nbqchem+1
355        endif
356        if (noms(iq).eq."h2o2") then
357          igcm_h2o2=iq
358          mmol(igcm_h2o2)=34.
359          count=count+1
360          nbqchem=nbqchem+1
361        endif
[460]362        if (noms(iq).eq."ch4") then
363          igcm_ch4=iq
364          mmol(igcm_ch4)=16.
365          count=count+1
366          nbqchem=nbqchem+1
367        endif
[38]368        if (noms(iq).eq."n2") then
369          igcm_n2=iq
370          mmol(igcm_n2)=28.
371          count=count+1
372          nbqchem=nbqchem+1
373        endif
374        if (noms(iq).eq."ar") then
375          igcm_ar=iq
376          mmol(igcm_ar)=40.
377          count=count+1
378          nbqchem=nbqchem+1
379        endif
380        if (noms(iq).eq."h2o_vap") then
381          igcm_h2o_vap=iq
382          mmol(igcm_h2o_vap)=18.
383          count=count+1
384          nbqchem=nbqchem+1
385        endif
386        if (noms(iq).eq."h2o_ice") then
387          igcm_h2o_ice=iq
388          mmol(igcm_h2o_ice)=18.
389          count=count+1
390          nbqchem=nbqchem+1
391        endif
392        ! Other stuff: e.g. for simulations using co2 + neutral gaz
393        if (noms(iq).eq."Ar_N2") then
394          igcm_ar_n2=iq
395          mmol(igcm_ar_n2)=30.
396          count=count+1
397        endif
398      enddo ! of do iq=1,nqmx
399     
400      ! check that we identified all tracers:
401      if (count.ne.nqmx) then
402        write(*,*) "inichim_newstart: found only ",count," tracers"
403        write(*,*) "               expected ",nqmx
404        do iq=1,count
405          write(*,*)'      ',iq,' ',trim(noms(iq))
406        enddo
407        stop
408      else
409       write(*,*)"inichim_newstart: found all expected tracers, namely:"
410        do iq=1,nqmx
411          write(*,*)'      ',iq,' ',trim(noms(iq))
412        enddo
413      endif
414
415      ! build local chemistry tracer index array nqchem(:)
416      ! as in the old days, water vapour is last and water ice,
417      ! if present, is just before water vapour
[460]418      do iq=1,16 ! loop so as to avoid out-of-bounds on array
[164]419      if (iq==1) nqchem(i)=igcm_co2
420      if (iq==2) nqchem(i)=igcm_co
421      if (iq==3) nqchem(i)=igcm_o
422      if (iq==4) nqchem(i)=igcm_o1d
423      if (iq==5) nqchem(i)=igcm_o2
424      if (iq==6) nqchem(i)=igcm_o3
425      if (iq==7) nqchem(i)=igcm_h
426      if (iq==8) nqchem(i)=igcm_h2
427      if (iq==9) nqchem(i)=igcm_oh
428      if (iq==10) nqchem(i)=igcm_ho2
429      if (iq==11) nqchem(i)=igcm_h2o2
[460]430      if (iq==12) nqchem(i)=igcm_ch4
431      if (iq==13) nqchem(i)=igcm_n2
432      if (iq==14) nqchem(i)=igcm_ar
433      if (iq==15) nqchem(i)=igcm_h2o_ice
434      if (iq==16) nqchem(i)=igcm_h2o_vap
[164]435      enddo
[38]436
437! Load in chemistry data for initialization:
438c----------------------------------------------------------------------
439c Order of Major species in file:
440c
441c    1=CO2
442c    2=AR
443c    3=N2 
444c    4=O2 
445c    5=CO 
446c    6=O   
447c    7=H2
448c
449c Order of Minor species in files:
450c
451c    1=H 
452c    2=OH
453c    3=HO2
454c    4=H2O
455c    5=H2O2
456c    6=O1D
457c    7=O3
458c
459c----------------------------------------------------------------------
460c Major species:
461        aa(igcm_co2) = 1
462        aa(igcm_ar)  = 2
463        aa(igcm_n2)  = 3
464        aa(igcm_o2)  = 4
465        aa(igcm_co)  = 5
466        aa(igcm_o)   = 6
467        aa(igcm_h2)  = 7
468c Minor species:
469        aa(igcm_h)       = 8
470        aa(igcm_oh)      = 9
471        aa(igcm_ho2)     = 10
472        aa(igcm_h2o_vap) = 11
473        aa(igcm_h2o2)    = 12
474        aa(igcm_o1d)     = 13
475        aa(igcm_o3)      = 14
476c----------------------------------------------------------------------
477      open(210, iostat=ierr,file=
478     & datafile(1:lnblnk(datafile))//'/atmosfera_LMD_may.dat')
479      if (ierr.ne.0) then
480        write(*,*)'Error : cannot open file atmosfera_LMD_may.dat '
481        write(*,*)'(in aeronomars/inichim.F)'
482        write(*,*)'It should be in :', datafile(1:lnblnk(datafile)),'/'
483        write(*,*)'1) You can change this directory address in '
484        write(*,*)'   file phymars/datafile.h'
485        write(*,*)'2) If necessary atmosfera_LMD_may.dat (and others)'
486        write(*,*)'   can be obtained online on:'
487        write(*,*)' http://www.lmd.jussieu.fr/~forget/datagcm/datafile'
488         STOP
489      endif
490      open(220, iostat=ierr,file=
491     & datafile(1:lnblnk(datafile))//'/atmosfera_LMD_min.dat')
492      if (ierr.ne.0) then
493        write(*,*)'Error : cannot open file atmosfera_LMD_min.dat '
494        write(*,*)'(in aeronomars/inichim.F)'
495        write(*,*)'It should be in :', datafile(1:lnblnk(datafile)),'/'
496        write(*,*)'1) You can change this directory address in '
497        write(*,*)'   file phymars/datafile.h'
498        write(*,*)'2) If necessary atmosfera_LMD_min.dat (and others)'
499        write(*,*)'   can be obtained online on:'
500        write(*,*)' http://www.lmd.jussieu.fr/~forget/datagcm/datafile'
501         STOP
502      endif
503
504      read(210,*) tmp
505      read(220,*) tmp
506        do l=1,252
507          read(210,112) tfile(l),zfile(l),nxf(l),(vmrf(l,n),n=1,7)
508          zfile(l)=zfile(l)*100.              !pressure (Pa)
509          sigfile(l)=zfile(l)/zfile(1)       
510        enddo
511        do l=1,252
512          read(220,113)zfilemin(l),(vmrf(l,n),n=8,14)
513          zfilemin(l)=zfilemin(l)*1000.       !height (m)
514          do n=1,14
515            nf(l,n)=vmrf(l,n)*nxf(l)
516          enddo
517        enddo
518        close(210)
519        close(220)
520
521c flag thermo set to 1 or 0 in newstart
522c inil=33 for initialisation of thermosphere only       
523c inil=1 for initialisation of all layers       
524        if (thermo.eq.1) then
525        inil=33
526        else
527        inil=1
528        endif
529
530! Initialization of tracers
531
532      do i=1,iip1
533       do j=1,jjp1
534        do l=inil,llm
535
536c          zgcm=sig(l)
537          zgcm=sig(l)*ps(i,j)
538          densconc=0.
539          nt=0.
540
541          do n=1,14
542            call intrplf(zgcm,vmrint,zfile,nf(1,n),252)
543c            call intrplf(zgcm,vmrint,sigfile,nf(1,n),252)
544            ni(n)=vmrint
545
546            densconc=densconc+ni(n)*mmolf(n)
547            nt=nt+ni(n)
548          enddo
549
550          mmean=densconc/nt                     ! in amu
551
552          if (nq.ne.nbqchem) then ! don't initialize water vapour
553            do n=1,nq-1
554            pq(i,j,l,nqchem(n))=
555     &                 ni(aa(nqchem(n)))/nt*mmol(nqchem(n))/mmean
556            if(i.eq.iip1) pq(i,j,l,nqchem(n))=pq(1,j,l,nqchem(n))
557            enddo
558            if (water .and. (nq.gt.nbqchem-2)) then
559              pq(i,j,l,nqchem(nq)) = 0.
560              if(i.eq.iip1) pq(i,j,l,nqchem(nq))=pq(1,j,l,nqchem(nq))
561            else
562              pq(i,j,l,nqchem(nq))=
563     &              ni(aa(nqchem(nq)))/nt*mmol(nqchem(nq))/mmean
564              if(i.eq.iip1) pq(i,j,l,nqchem(nq))=pq(1,j,l,nqchem(nq))
565            endif
566          endif ! of if (nq.ne.nbqchem)
567         
568          if (nq.eq.nbqchem) then ! also initialize water vapour
569            do n=1,nq-2
570              pq(i,j,l,nqchem(n))=
571     &                   ni(aa(nqchem(n)))/nt*mmol(nqchem(n))/mmean
572              if(i.eq.iip1) pq(i,j,l,nqchem(n))=pq(1,j,l,nqchem(n))
573            enddo
574            if (water) then
575              pq(i,j,l,nqchem(nq-1)) = 0.
576              if(i.eq.iip1) then
577                pq(i,j,l,nqchem(nq-1))=pq(1,j,l,nqchem(nq-1))
578              endif
579              pq(i,j,l,nqchem(nq))=
580     &                   ni(aa(nqchem(nq)))/nt*mmol(nqchem(nq))/mmean
581              if(i.eq.iip1) pq(i,j,l,nqchem(nq))=pq(1,j,l,nqchem(nq))
582            else
583              do n=nq-1,nq
584                pq(i,j,l,nqchem(nq))=
585     &                   ni(aa(nqchem(nq)))/nt*mmol(nqchem(nq))/mmean
586                if(i.eq.iip1) pq(i,j,l,nqchem(nq))=pq(1,j,l,nqchem(nq))
587              enddo
588            endif ! of if (water)
589          endif ! of if (nq.eq.nqmx)
590        enddo     !nlayer
591       enddo      !ngrid_j
592      enddo       !ngrid_i
593
594cccccccccccccccccccccccccccccc     
595c  make sure that sum of q = 1     
596c dominent species is = 1 - sum(all other species)     
597cccccccccccccccccccccccccccccc     
598!      iqmax=nqchem_min
599       iqmax=1
600     
601!      if ((nqmx-nqchem_min).gt.10) then
602      if (nbqchem.gt.10) then
603       do l=1,llm
604        do j=1,jjp1
605          do i=1,iip1
606!           do iq=nqchem_min,nqmx
607            do iiq=1,nbqchem
608              iq=nqchem(iiq)
609              if ( (pq(i,j,l,iq).gt.pq(i,j,l,iqmax)).and.
610     .           (noms(iq).ne."ice") )  then
611              iqmax=iq
612              endif
613            enddo
614            pq(i,j,l,iqmax)=1.
615            qtot(i,j,l)=0
616!           do iq=nqchem_min,nqmx
617            do iiq=1,nbqchem
618             iq=nqchem(iiq)
619             if ( (iq.ne.iqmax).and.
620     .           (noms(iq).ne."ice") ) then       
621               pq(i,j,l,iqmax)=pq(i,j,l,iqmax)-pq(i,j,l,iq)       
622             endif
623            enddo !iq
624!            do iq=nqchem_min,nqmx
625            do iiq=1,nbqchem
626             iq=nqchem(iiq)
627              if (noms(iq).ne."ice") then
628                qtot(i,j,l)=qtot(i,j,l)+pq(i,j,l,iq)
629              endif
630c            if (i.eq.1.and.j.eq.1.and.l.eq.1) write(*,*) 'qtot(i,j,l)',
631c     $    qtot(i,j,l)
632            enddo !iq
633            if (i.eq.1.and.j.eq.1) write(*,*) 'inichim_newstart: ',
634     $        'qtot(i,j,l)=',qtot(i,j,l)
635          enddo !i   
636         enddo !j   
637       enddo !l 
638      endif ! of if (nbqchem.gt.10)
639ccccccccccccccccccccccccccccccc
640
64166      format(i2,6(1x,e11.5))
642112     format(7x,f7.3,3x,e12.6,3x,e12.6,7(3x,e12.6))
643113     format(1x,f6.2,7(3x,e12.6))
644
645      end
Note: See TracBrowser for help on using the repository browser.