source: trunk/LMDZ.MARS/libf/aeronomars/moldiffcoeff_red.F @ 1119

Last change on this file since 1119 was 1047, checked in by emillour, 11 years ago

Mars GCM:

  • IMPORTANT CHANGE: Removed all reference/use of ngridmx (dimphys.h) in routines (necessary prerequisite to using parallel dynamics); in most cases this just means adding 'ngrid' as routine argument, and making local saved variables allocatable (and allocated at first call). In the process, had to convert many *.h files to equivalent modules: yomaer.h => yomaer_h.F90 , surfdat.h => surfdat_h.F90 , comsaison.h => comsaison_h.F90 , yomlw.h => yomlw_h.F90 , comdiurn.h => comdiurn_h.F90 , dimradmars.h => dimradmars_mod.F90 , comgeomfi.h => comgeomfi_h.F90, comsoil.h => comsoil_h.F90 , slope.h => slope_mod.F90
  • Also updated EOF routines, everything is now in eofdump_mod.F90
  • Removed unused routine lectfux.F (in dyn3d)

EM
File size: 8.7 KB
Line 
1      subroutine moldiffcoeff_red(dij,indic,gcmind,ncompdiff2)
2
3       use tracer_mod, only: nqmx, noms, mmol
4       IMPLICIT NONE
5c=======================================================================
6c   subject:
7c   --------
8c   Computing molecular diffusion coefficients
9c   following Nair 94 (pg 131)
10c   author:  MAC 2002
11c   ------
12c
13c=======================================================================
14#include "dimensions.h"
15#include "dimphys.h"
16#include "callkeys.h"
17!#include "comdiurn.h"
18#include "chimiedata.h"
19!#include "tracer.h"
20!#include "conc.h"
21#include "diffusion.h"
22
23c-----------------------------------------------------------------------
24c    Input/Output
25c    ------------
26c       integer,parameter :: ncompmoldiff = 12
27        integer ncompdiff2
28        integer gcmind(ncompdiff2)
29      real dij(ncompdiff2,ncompdiff2)
30      integer indic(nqmx)
31
32c    Local variables:
33c    ---------------
34      INTEGER nq, n, nn, i,iq
35cccccccccccccccccccccccccccccccccccccccccccccccccccccccc
36c     tracer numbering in the molecular diffusion
37cccccccccccccccccccccccccccccccccccccccccccccccccccccccc
38
39      real :: dijh2co,dijh2n2,dijh2co2,dijh2o2,dijho,dijref
40!       integer :: i_h2,i_h,i_o
41        integer :: g_h2,g_h,g_o
42!      integer,parameter :: i_o   = 1
43!      integer,parameter :: i_n2   = 2
44!      integer,parameter :: i_co   = 3
45!      integer,parameter :: i_ar  = 4
46!      integer,parameter :: i_h2   = 5
47!      integer,parameter :: i_h    = 6
48!      integer,parameter :: i_o2   = 7
49!      integer,parameter :: i_oh  = 8
50!      integer,parameter :: i_ho2  = 9
51!      integer,parameter :: i_h2o = 10
52!      integer,parameter :: i_h2o2  = 11
53!      integer,parameter :: i_o1d   = 12
54!      integer,parameter :: i_o3   = 13
55!      integer,parameter :: i_n    = 13
56!      integer,parameter :: i_no   = 14
57!      integer,parameter :: i_no2  = 15
58!      integer,parameter :: i_n2d  = 17
59!      integer,parameter :: i_oplus = 18
60!      integer,parameter :: i_co2    = 16
61!      integer,parameter :: i_oplus = 17
62!      integer,parameter :: i_hplus = 18
63
64! Tracer indexes in the GCM:
65!      integer,save :: g_co2=0
66!      integer,save :: g_co=0
67!      integer,save :: g_o=0
68!      integer,save :: g_o1d=0
69!      integer,save :: g_o2=0
70!      integer,save :: g_o3=0
71!      integer,save :: g_h=0
72!      integer,save :: g_h2=0
73!      integer,save :: g_oh=0
74!      integer,save :: g_ho2=0
75!      integer,save :: g_h2o2=0
76!      integer,save :: g_n2=0
77!      integer,save :: g_ar=0
78!      integer,save :: g_h2o=0
79!      integer,save :: g_n=0
80!      integer,save :: g_no=0
81!      integer,save :: g_no2=0
82!      integer,save :: g_n2d=0
83!      integer,save :: g_oplus=0
84!      integer,save :: g_hplus=0
85
86!      integer,save :: gcmind(ncompdiff)
87
88      real dnh
89      logical,save :: firstcall=.true.
90      logical,parameter :: outputcoeffs=.false. ! to output 'coeffs.dat' file,
91                                                ! set outputcoeffs=.true.
92
93! Initializations at first call (and some sanity checks)
94      if (firstcall) then
95        ! identify the indexes of the tracers we'll need
96!        g_co2=igcm_co2
97!        if (g_co2.eq.0) then
98!          write(*,*) "moldiffcoeff: Error; no CO2 tracer !!!"
99!          stop
100!        endif
101!        g_n2=igcm_n2
102!        if (g_n2.eq.0) then
103!          write(*,*) "moldiffcoeff: Error; no N2 tracer !!!"
104!          stop
105!        endif
106!        g_ar=igcm_ar
107!        if (g_ar.eq.0) then
108!          write(*,*) "moldiffcoeff: Error; no Ar tracer !!!"
109!          stop
110!        endif       
111!        g_h2=igcm_h2
112!        if (g_h2.eq.0) then
113!          write(*,*) "moldiffcoeff: Error; no H2 tracer !!!"
114!          stop
115!        endif
116!        g_h=igcm_h
117!        if (g_h.eq.0) then
118!          write(*,*) "moldiffcoeff: Error; no H tracer !!!"
119!          stop
120!        endif
121!        g_co=igcm_co
122!        if (g_co.eq.0) then
123!          write(*,*) "moldiffcoeff: Error; no CO tracer !!!"
124!          stop
125!        endif
126!        g_o2=igcm_o2
127!        if (g_o2.eq.0) then
128!          write(*,*) "moldiffcoeff: Error; no O2 tracer !!!"
129!          stop
130!        endif
131!        g_oh=igcm_oh
132!        if (g_oh.eq.0) then
133!          write(*,*) "moldiffcoeff: Error; no OH tracer !!!"
134!          stop
135!        endif
136!        g_ho2=igcm_ho2
137!        if (g_ho2.eq.0) then
138!          write(*,*) "moldiffcoeff: Error; no HO2 tracer !!!"
139!          stop
140!        endif
141!        g_h2o=igcm_h2o_vap
142!        if (g_h2o.eq.0) then
143!          write(*,*) "moldiffcoeff: Error; no H2O tracer !!!"
144!          stop
145!        endif
146!        g_h2o2=igcm_h2o2
147!        if (g_h2o2.eq.0) then
148!          write(*,*) "moldiffcoeff: Error; no H2O2 tracer !!!"
149!          stop
150!        endif
151!        g_o1d=igcm_o1d
152!        if (g_h.eq.0) then
153!          write(*,*) "moldiffcoeff: Error; no O1d tracer !!!"
154!          stop
155!        endif
156!        g_o3=igcm_o3
157!        if (g_o3.eq.0) then
158!          write(*,*) "moldiffcoeff: Error; no O3 tracer !!!"
159!          stop
160!        endif
161!        g_n=igcm_n
162!        if (g_n.eq.0) then
163!          write(*,*) "moldiffcoeff: Error; no N tracer !!!"
164!          stop
165!        endif
166!        g_no=igcm_no
167!        if (g_no.eq.0) then
168!          write(*,*) "moldiffcoeff: Error; no NO tracer !!!"
169!          stop
170!        endif
171!        g_no2=igcm_no2
172!        if (g_no2.eq.0) then
173!          write(*,*) "moldiffcoeff: Error; no NO2 tracer !!!"
174!          stop
175!        endif
176!        g_n2d=igcm_n2d
177!        if (g_n2d.eq.0) then
178!          write(*,*) "moldiffcoeff: Error; no N2(D) tracer !!!"
179!          stop
180!        endif
181!        g_oplus=igcm_oplus
182!        if (g_oplus .eq. 0) then
183!        write(*,*) "moldiffcoeff: Error; no Oplus tracer !!!"
184!        stop
185!        endif
186!       g_hplus=igcm_hplus
187!        if (g_hplus .eq. 0) then
188!        write(*,*) "moldiffcoeff: Error; no Hplus tracer !!!"
189!        stop
190!        endif
191!        g_o=igcm_o
192!        if (g_o.eq.0) then
193!          write(*,*) "moldiffcoeff: Error; no O tracer !!!"
194!          stop
195!        endif
196       
197c
198cccccccccccccccccccccccccccccccccccccccccccccccccccccccc
199c    fill array to relate local indexes to gcm indexes
200cccccccccccccccccccccccccccccccccccccccccccccccccccccccc
201
202!        gcmind(i_co2)  =   g_co2
203!        gcmind(i_n2)  =   g_n2
204!        gcmind(i_ar)  =   g_ar
205!        gcmind(i_h2) =   g_h2
206!        gcmind(i_h)  =   g_h
207!        gcmind(i_co)   =   g_co
208!        gcmind(i_o2) =   g_o2
209!        gcmind(i_oh)=   g_oh
210!        gcmind(i_ho2)  =   g_ho2
211!        gcmind(i_h2o) = g_h2o
212!        gcmind(i_h2o2)= g_h2o2
213!        gcmind(i_o1d) = g_o1d
214!        gcmind(i_o3) = g_o3
215!        gcmind(i_n)= g_n
216!        gcmind(i_no) = g_no
217!        gcmind(i_no2) = g_no2
218!        gcmind(i_n2d) = g_n2d
219!        gcmind(i_oplus) =  g_oplus
220!        gcmind(i_hplus) = g_hplus
221!        gcmind(i_o)   =   g_o
222
223c
224cccccccccccccccccccccccccccccccccccccccccccccccccccccccc
225        firstcall= .false.
226      endif ! of if (firstcall)
227
228        dijh2co = 0.0000651
229        dijh2n2 = 0.0000674
230        dijh2o2 = 0.0000697
231        dijh2co2 = 0.0000550
232        dijho = 0.000114
233
234!      dij(i_h2,i_co)   = 0.0000651
235!      dij(i_h2,i_n2)   = 0.0000674
236!      dij(i_h2,i_o2)   = 0.0000697
237!      dij(i_h2,i_co2)  = 0.0000550
238!      dij(i_h2,i_h2)   = 0.0
239!      dij(i_h2,i_h)    = 0.0
240!      dij(i_h2,i_h2o)  = 0.0   !0003
241!      dij(i_h2,i_h2o2) = 0.0   !0003
242!      dij(i_h2,i_o3)   = 0.0   !0003
243!      dij(i_h2,i_o)    = 0.0
244!      dij(i_h2,i_ar)   = 0.0
245!      dij(i_h2,i_n)    = 0.0
246
247!c      dij(i_h,i_o)     = 0.0000144
248!      dij(i_h,i_o)     = 0.000114
249
250! find h2, h and o index in gcm
251! these species are used to define the diffusion coefficients
252
253        do n=1,nqmx
254        if (noms(n) .eq. 'h2') g_h2=n
255        if (noms(n) .eq. 'h') g_h=n
256        if (noms(n) .eq. 'o') g_o=n
257        enddo
258        print*,'moldiffcoeff_red: gh2',g_h2,g_h,g_o
259
260       print*,'moldiffcoeff_red: COEFF CALC'
261
262      do n=1,ncompdiff2
263        dijref=0.
264        if (noms(gcmind(n)) .eq. 'co') dijref=dijh2co
265        if (noms(gcmind(n)) .eq. 'n2') dijref=dijh2n2
266        if (noms(gcmind(n)) .eq. 'o2') dijref=dijh2o2
267        if (noms(gcmind(n)) .eq. 'co2') dijref=dijh2co2
268!       print*,'test',n,dijref
269        if (dijref .gt. 0.0) then
270          do nn=n,ncompdiff2
271            dij(nn,n)=dijref
272     &                  *sqrt(mmol(g_h2)/mmol(gcmind(nn)))
273            if(n.eq.nn) dij(nn,n)=1.0
274            dij(n,nn)=dij(nn,n)
275          enddo
276        endif
277        if (dijref .eq. 0.0) then
278        dijref=dijho
279          dnh=dijref*sqrt(mmol(g_o)/mmol(gcmind(n)))
280          do nn=n,ncompdiff2
281            dij(nn,n)=dnh*sqrt(mmol(g_h)/mmol(gcmind(nn)))
282            if(n.eq.nn) dij(nn,n)=1.0
283            dij(n,nn)=dij(nn,n)
284          enddo
285        endif
286      enddo
287
288      if (outputcoeffs) then
289       ! output coefficients in 'coeffs.dat' file
290       open(56,file='coeffs.dat',status='unknown')
291       do n=1,ncompdiff2
292        do nn=n,ncompdiff2
293          write(56,*) n,nn,dij(n,nn)    !*1.e5/1.381e-23/(273**1.75)
294        enddo
295       enddo
296       close(56)
297      endif
298
299
300      return   
301      end
302     
Note: See TracBrowser for help on using the repository browser.