MAIC-2  Revision 19
 All Classes Files Functions Variables
diff_trans.F90
Go to the documentation of this file.
1 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 !
3 ! Subroutine : d i f f _ t r a n s
4 !
5 !> @file
6 !!
7 !! Diffusive transport of water in the Martian atmosphere.
8 !!
9 !! @section Copyright
10 !!
11 !! Copyright 2010-2013 Ralf Greve, Bjoern Grieger, Oliver J. Stenzel
12 !!
13 !! @section License
14 !!
15 !! This file is part of MAIC-2.
16 !!
17 !! MAIC-2 is free software: you can redistribute it and/or modify
18 !! it under the terms of the GNU General Public License as published by
19 !! the Free Software Foundation, either version 3 of the License, or
20 !! (at your option) any later version.
21 !!
22 !! MAIC-2 is distributed in the hope that it will be useful,
23 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
24 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 !! GNU General Public License for more details.
26 !!
27 !! You should have received a copy of the GNU General Public License
28 !! along with MAIC-2. If not, see <http://www.gnu.org/licenses/>.
29 !<
30 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
31 
32 !-------------------------------------------------------------------------------
33 !> Diffusive transport of water in the Martian atmosphere.
34 !<------------------------------------------------------------------------------
35 subroutine diff_trans(dtime)
36 
37 use maic2_types
39 
40 implicit none
41 
42 real(dp), intent(in) :: dtime
43 
44 integer(i4b) :: l
45 real(dp), dimension(0:LMAX) :: sle_a0, sle_a1, sle_a2, sle_b
46 real(dp) :: ratio_water_np_sp
47 real(dp) :: water_mean, water_mean_normalized
48 
49 #if SOLV_DIFF==0
50 
51 stop ' diff_trans: Option SOLV_DIFF==0 not available any more!'
52 
53 #elif SOLV_DIFF==1 /* Explicit scheme (Euler forward) */
54 
55 !-------- Direct solution of the Euler forward scheme --------
56 
57 l=0
58 water_new(l) = water(l) + dtime &
59  * ( evap(l) - cond(l) &
60  + diff_aux(l) &
61  * ( cos_phi_cb2(l)*(water(l+1)-water(l)) &
62  *dphi_inv(l+1) ) &
63  )
64 
65 do l=1, lmax-1
66  water_new(l) = water(l) + dtime &
67  * ( evap(l) - cond(l) &
68  + diff_aux(l) &
69  * ( cos_phi_cb2(l)*(water(l+1)-water(l)) &
70  *dphi_inv(l+1) &
71  -cos_phi_cb1(l)*(water(l)-water(l-1)) &
72  *dphi_inv(l) ) &
73  )
74 end do
75 
76 l=lmax
77 water_new(l) = water(l) + dtime &
78  * ( evap(l) - cond(l) &
79  + diff_aux(l) &
80  * ( -cos_phi_cb1(l)*(water(l)-water(l-1)) &
81  *dphi_inv(l) ) &
82  )
83 
84 #elif SOLV_DIFF==2 /* Implicit scheme (Euler backward) */
85 
86 !-------- Setting of matrix and vector components --------
87 
88 l=0
89 sle_a1(l) = 1.0_dp + dtime * diff_aux(l)*cos_phi_cb2(l)*dphi_inv(l+1)
90 sle_a2(l) = - dtime * diff_aux(l)*cos_phi_cb2(l)*dphi_inv(l+1)
91 sle_b(l) = water(l) + dtime * (evap(l)-cond(l))
92 
93 do l=1, lmax-1
94  sle_a0(l) = - dtime * diff_aux(l)*cos_phi_cb1(l)*dphi_inv(l)
95  sle_a1(l) = 1.0_dp + dtime * diff_aux(l) &
96  * ( cos_phi_cb1(l)*dphi_inv(l) &
97  +cos_phi_cb2(l)*dphi_inv(l+1) )
98  sle_a2(l) = - dtime * diff_aux(l)*cos_phi_cb2(l)*dphi_inv(l+1)
99  sle_b(l) = water(l) + dtime * (evap(l)-cond(l))
100 end do
101 
102 l=lmax
103 sle_a0(l) = - dtime * diff_aux(l)*cos_phi_cb1(l)*dphi_inv(l)
104 sle_a1(l) = 1.0_dp + dtime * diff_aux(l)*cos_phi_cb1(l)*dphi_inv(l)
105 sle_b(l) = water(l) + dtime * (evap(l)-cond(l))
106 
107 !-------- Solution of the system of linear equations --------
108 
109 call tri_sle(sle_a0, sle_a1, sle_a2, water_new, sle_b, lmax)
110 
111 #elif SOLV_DIFF==3 /* Instantaneous mixing (optional north-south gradient) */
112 
113 !-------- Ratio of the atmospheric water content at the poles
114 ! (north pole relative to south pole) --------
115 
116 #if defined(RATIO_WATER_NP_SP)
117 ratio_water_np_sp = ratio_water_np_sp ! prescribed north-south gradient
118 #else
119 ratio_water_np_sp = 1.0_dp ! no north-south gradient,
120  ! constant water content everywhere on the planet
121 #endif
122 
123 !-------- Predictor step without transport --------
124 
125 do l=0, lmax
126  water_new(l) = water(l) + dtime * ( evap(l) - cond(l) )
127 end do
128 
129 !-------- Mixing --------
130 
131 ! ------ Mean water content from predictor
132 
133 water_mean = 0.0_dp
134 
135 do l=0, lmax
136  water_mean = water_mean &
137  + 0.5_dp*water_new(l)*(sin_phi_cb2(l)-sin_phi_cb1(l))
138 end do
139 
140 ! ------ Water content from linear distribution
141 ! with normalized value 1 at the south pole
142 
143 do l=0, lmax
144  water_new(l) = 0.5_dp * (ratio_water_np_sp + 1.0_dp) &
145  + pi_inv * (ratio_water_np_sp - 1.0_dp) * phi_node(l)
146 end do
147 
148 ! ------ Correction such that the correct mean water content results
149 
150 water_mean_normalized = 0.0_dp
151 
152 do l=0, lmax
153  water_mean_normalized = water_mean_normalized &
154  + 0.5_dp*water_new(l)*(sin_phi_cb2(l)-sin_phi_cb1(l))
155 end do
156 
157 water_new = water_new * (water_mean/water_mean_normalized)
158 
159 #else /* Wrong value of parameter SOLV_DIFF */
160 
161 stop ' Wrong value of parameter SOLV_DIFF (must be 1, 2 or 3)!'
162 
163 #endif
164 
165 end subroutine diff_trans
166 !