MAIC-2  Revision 19
 All Classes Files Functions Variables
evaporation.f90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Module : e v a p o r a t i o n
4 !
5 !> @file
6 !!
7 !! Computation of the evaporation rate
8 !! (buoyant-diffusion approach by Ingersoll).
9 !!
10 !! @section Copyright
11 !!
12 !! Copyright 2010-2013 Ralf Greve, Bjoern Grieger, Oliver J. Stenzel
13 !!
14 !! @section License
15 !!
16 !! This file is part of MAIC-2.
17 !!
18 !! MAIC-2 is free software: you can redistribute it and/or modify
19 !! it under the terms of the GNU General Public License as published by
20 !! the Free Software Foundation, either version 3 of the License, or
21 !! (at your option) any later version.
22 !!
23 !! MAIC-2 is distributed in the hope that it will be useful,
24 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
25 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
26 !! GNU General Public License for more details.
27 !!
28 !! You should have received a copy of the GNU General Public License
29 !! along with MAIC-2. If not, see <http://www.gnu.org/licenses/>.
30 !<
31 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
32 
33 !-------------------------------------------------------------------------------
34 !> Computation of the evaporation rate
35 !! (buoyant-diffusion approach by Ingersoll).
36 !<------------------------------------------------------------------------------
38 
39 use maic2_types
40 
41 implicit none
42 real(dp), private :: evap_fact, gamma_reg, &
43  r_univ, mol_w, mol_c, diff_w_c, visc_c, g
44 
45 contains
46 
47 !-------------------------------------------------------------------------------
48 !> Setting of parameters.
49 !<------------------------------------------------------------------------------
50  subroutine setevappar_bd(evap_fact_para, gamma_reg_para, &
51  r_univ_para, mol_w_para, mol_c_para, &
52  diff_w_c_para, visc_c_para, g_para)
53 
54  implicit none
55 
56  real(dp), optional :: evap_fact_para, gamma_reg_para, &
57  r_univ_para, mol_w_para, mol_c_para, &
58  diff_w_c_para, visc_c_para, g_para
59 
60  if ( present(evap_fact_para) ) then
61  evap_fact = evap_fact_para ! evaporation factor
62  else
63  evap_fact = 1.0_dp
64  end if
65 
66  if ( present(gamma_reg_para) ) then
67  gamma_reg = gamma_reg_para ! regolith-insolation coefficient [in m]
68  else
69  gamma_reg = 1.11e+11_dp
70  end if
71 
72  if ( present(r_univ_para) ) then
73  r_univ = r_univ_para ! universal gas constant [J/(mol*K)]
74  else
75  r_univ = 8.314_dp
76  end if
77 
78  if ( present(mol_w_para) ) then
79  mol_w = mol_w_para ! molar mass of water [kg/mol]
80  else
81  mol_w = 1.802e-02_dp
82  end if
83 
84  if ( present(mol_c_para) ) then
85  mol_c = mol_c_para ! molar mass of CO2 [kg/mol]
86  else
87  mol_c = 4.401e-02_dp
88  end if
89 
90  if ( present(diff_w_c_para) ) then
91  diff_w_c = diff_w_c_para ! diffusion coefficient of water in CO2 [m2/s]
92  else
93  diff_w_c = 1.4e-03_dp
94  end if
95 
96  if ( present(visc_c_para) ) then
97  visc_c = visc_c_para ! kinematic viscosity of CO2 [m2/s]
98  else
99  visc_c = 6.93e-04_dp
100  end if
101 
102  if ( present(g_para) ) then
103  g = g_para ! gravity acceleration [m/s2]
104  else
105  g = 3.7_dp
106  end if
107 
108  end subroutine setevappar_bd
109 
110 !-------------------------------------------------------------------------------
111 !> Computation of evaporation.
112 !<------------------------------------------------------------------------------
113  subroutine getevap_bd(temp, temp_amp, p, H, evap)
114 
115  implicit none
116 
117  real(dp) :: temp(:), temp_amp(:), p(:), h(:), evap(:)
118 
119  integer(i4b) :: i, hr, n
120  real(dp) :: inv_visc_c_2, one_third, one_eighth
121  real(dp) :: p_sat, delta_eta, rho, &
122  delta_rho_rho_1, delta_rho_rho_2, delta_rho_rho, &
123  inv_length, inv_gamma_reg, &
124  temp_hr, evap_hr
125  real(dp) :: sat_pressure
126  real(dp), dimension(8) :: cos_factor
127  real(dp), parameter :: pi = 3.141592653589793_dp
128 
129  n = size(temp)
130 
131  inv_visc_c_2 = 1.0_dp/visc_c**2
132  inv_gamma_reg = 1.0_dp/gamma_reg
133  one_third = 1.0_dp/3.0_dp
134  one_eighth = 1.0_dp/8.0_dp
135  cos_factor = (/(cos(2.0_dp*pi*one_eighth*real(hr,dp)),hr=1,8)/)
136 
137  do i = 1, n
138 
139  evap(i) = 0.0_dp
140 
141  do hr=1, 8 ! counter for one Martian day (sampled in 3-hour intervals)
142 
143  temp_hr = temp(i) - temp_amp(i)*cos_factor(hr)
144  sat_pressure = p_sat(temp_hr)
145  delta_eta = mol_w*sat_pressure/(mol_c*p(i))
146  rho = mol_c*p(i)/(r_univ*temp_hr)
147 
148  delta_rho_rho_1 = (mol_c-mol_w)*sat_pressure
149  delta_rho_rho_2 = mol_c*p(i)-(mol_c-mol_w)*sat_pressure
150 
151  if (delta_rho_rho_1 <= delta_rho_rho_2) then
152  ! physically reasonable case
153  delta_rho_rho = delta_rho_rho_1/delta_rho_rho_2
154  else ! physically unreasonable case
155  delta_rho_rho = 1.0_dp
156  end if
157 
158  inv_length = (delta_rho_rho*g*inv_visc_c_2)**one_third
159 
160  evap_hr = evap_fact * 0.17_dp*delta_eta*rho*diff_w_c*inv_length
161 
162  evap(i) = evap(i) + one_eighth*evap_hr
163 
164  end do
165 
166  if (h(i) < 0.0_dp) evap(i) = evap(i) * exp(inv_gamma_reg*h(i))
167 
168  end do
169 
170  end subroutine getevap_bd
171 
172 end module evaporation
173 !