source: trunk/LMDZ.MARS/libf/aeronomars/hrtherm.F @ 937

Last change on this file since 937 was 705, checked in by emillour, 13 years ago

Mars GCM:

  • Added possibility to run with a varying EUV cycle following real one. The flag solvarmod=1 triggers this behaviour, with companion flag solvaryear=## , where ## is the Mars Year (from 23 to 30). Setting solvarmod=0 reverts to 'old' behaviour, where there is a constant EUV forcing throughout the run, set by the "solarcondate" flag.
  • Needs corresponding input data files ("param_v6" subdirectory of "EUV" subdirectory in "datadir").
  • Added files jthermcalc_e107.F and param_read_e107.F in "aeronomars", modified files euvheat.F90, hrtherm.F, chemthermos.F90, param_v4.h and param_read.F in "aeronomars" and files inifis.F, physiq.F and callkeys.h in "phymars".

FGG

File size: 3.6 KB
Line 
1c**********************************************************************
2
3      subroutine hrtherm(ig,euvmod,rm,nespeuv,tx,iz,zenit,zday,jtot)
4
5
6c     feb 2002        fgg           first version
7c     nov 2002        fgg           second version
8
9c**********************************************************************
10
11      implicit none
12
13c     common variables and constants
14
15
16      include 'dimensions.h'
17      include 'dimphys.h'
18      include 'param.h'
19      include 'param_v4.h'
20      include "callkeys.h"
21
22
23c    local parameters and variables
24
25      real       xabsi(nabs,nlayermx)                   !densities
26      real       jergs(ninter,nabs,nlayermx)
27     
28      integer    i,j,k,indexint          !indexes
29      character  dn
30
31
32c     input and output variables
33
34      integer    ig  ,euvmod
35      integer    nespeuv
36      real       rm(nlayermx,nespeuv)              !density matrix (cm^-3)
37      real       jtot(nlayermx)                    !output: heating rate(erg/s)
38      real       tx(nlayermx)                      !temperature
39      real       zenit
40      real       iz(nlayermx)
41      real       zday
42
43      ! tracer indexes for the EUV heating:
44!!! ATTENTION. These values have to be identical to those in chemthermos.F90
45!!! If the values are changed there, the same has to be done here  !!!
46      integer,parameter :: i_co2=1
47      integer,parameter :: i_o2=2
48      integer,parameter :: i_o=3
49      integer,parameter :: i_co=4
50      integer,parameter :: i_h=5
51      integer,parameter :: i_h2=8
52      integer,parameter :: i_h2o=9
53      integer,parameter :: i_h2o2=10
54      integer,parameter :: i_o3=12
55      integer,parameter :: i_n2=13
56      integer,parameter :: i_n=14
57      integer,parameter :: i_no=15
58      integer,parameter :: i_no2=17
59
60c*************************PROGRAM STARTS*******************************
61
62      !If nighttime, photoabsorption coefficient set to 0
63      if(zenit.gt.140.) then
64         dn='n'
65         else
66         dn='d'
67      end if
68      if(dn.eq.'n') then
69        do i=1,nlayermx                                   
70              jtot(i)=0.
71        enddo       
72        return
73      endif
74
75      !initializations
76      jergs(:,:,:)=0.
77      xabsi(:,:)=0.
78      jtot(:)=0.
79      !All number densities to a single array, xabsi(species,layer)
80      do i=1,nlayermx
81         xabsi(1,i)  = rm(i,i_co2)
82         xabsi(2,i)  = rm(i,i_o2)
83         xabsi(3,i)  = rm(i,i_o)
84         xabsi(4,i)  = rm(i,i_h2o)
85         xabsi(5,i)  = rm(i,i_h2)
86         xabsi(6,i)  = rm(i,i_h2o2)
87         !Only if O3, N or ion chemistry requested
88         if(euvmod.ge.1) then
89            xabsi(7,i)  = rm(i,i_o3)
90         endif
91         !Only if N or ion chemistry requested
92         if(euvmod.ge.2) then
93            xabsi(8,i)  = rm(i,i_n2)
94            xabsi(9,i)  = rm(i,i_n)
95            xabsi(10,i) = rm(i,i_no)
96            xabsi(13,i) = rm(i,i_no2)
97         endif
98         xabsi(11,i) = rm(i,i_co)
99         xabsi(12,i) = rm(i,i_h)
100      end do
101
102      !Calculation of photoabsortion coefficient
103      if(solvarmod.eq.0) then
104         call jthermcalc(ig,euvmod,rm,nespeuv,tx,iz,zenit)
105      else if (solvarmod.eq.1) then
106         call jthermcalc_e107(ig,euvmod,rm,nespeuv,tx,iz,zenit,zday)
107         do indexint=1,ninter
108            fluxtop(indexint)=1.
109         enddo
110      endif
111
112      !Total photoabsorption coefficient
113      do i=1,nlayermx
114         jtot(i)=0.
115        do j=1,nabs
116          do indexint=1,ninter
117            jergs(indexint,j,i) = jfotsout(indexint,j,i)
118     $              * xabsi (j,i) * fluxtop(indexint) 
119     $              / (0.5e9 * freccen(indexint))
120            jtot(i)=jtot(i)+jergs(indexint,j,i)
121          end do
122        end do
123      end do
124
125      return
126
127      end
128
Note: See TracBrowser for help on using the repository browser.