00001
00002
00003
00004
00005 !> @file
00006 !!
00007 !! Diffusive transport of water in the Martian atmosphere.
00008 !!
00009 !! @section Copyright
00010 !!
00011 !! Copyright 2010, 2011 Ralf Greve, Bjoern Grieger, Oliver J. Stenzel
00012 !!
00013 !! @section License
00014 !!
00015 !! This file is part of MAIC-2.
00016 !!
00017 !! MAIC-2 is free software: you can redistribute it and/or modify
00018 !! it under the terms of the GNU General Public License as published by
00019 !! the Free Software Foundation, either version 3 of the License, or
00020 !! (at your option) any later version.
00021 !!
00022 !! MAIC-2 is distributed in the hope that it will be useful,
00023 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
00024 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
00025 !! GNU General Public License for more details.
00026 !!
00027 !! You should have received a copy of the GNU General Public License
00028 !! along with MAIC-2. If not, see <http://www.gnu.org/licenses/>.
00029 !<
00030 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
00031
00032
00033 !> Diffusive transport of water in the Martian atmosphere.
00034 !<------------------------------------------------------------------------------
00035 subroutine diff_trans(dtime)
00036
00037 use maic2_types
00038 use maic2_variables
00039
00040 implicit none
00041
00042 real(dp), intent(in) :: dtime
00043
00044 integer(i4b) :: l
00045 real(dp), dimension(0:LMAX) :: sle_a0, sle_a1, sle_a2, sle_b
00046 real(dp) :: transport(LMAX)
00047 real(dp) :: water_mean
00048
00049 #if SOLV_DIFF==0 /* Scheme by Bjoern */
00050
00051
00052
00053 do l=1, LMAX
00054 transport(l) = dtime * DIFF_WATER_MAIC * 2_dp*acos(-1_dp)*cos_phi_cb1(l) &
00055 * ( water(l-1) - water(l) ) * dphi_inv(l)
00056 end do
00057
00058
00059
00060 l=0
00061 water_new(l) = water(l) + dtime * ( evap(l) - cond(l) ) &
00062 + ( 0_dp - transport(l+1) ) &
00063 / ( 2_dp * acos(-1_dp) * R**2 &
00064 * ( sin_phi_cb2(l) - sin_phi_cb1(l) ) &
00065 )
00066
00067 do l=1, LMAX-1
00068 water_new(l) = water(l) + dtime * ( evap(l) - cond(l) ) &
00069 + ( transport(l) - transport(l+1) ) &
00070 / ( 2_dp * acos(-1_dp) * R**2 &
00071 * ( sin_phi_cb2(l) - sin_phi_cb1(l) ) &
00072 )
00073 end do
00074
00075 l=LMAX
00076 water_new(l) = water(l) + dtime * ( evap(l) - cond(l) ) &
00077 + ( transport(l) - 0_dp ) &
00078 / ( 2_dp * acos(-1_dp) * R**2 &
00079 * ( sin_phi_cb2(l) - sin_phi_cb1(l) ) &
00080 )
00081
00082 #elif SOLV_DIFF==1 /* Explicit scheme (Euler forward) */
00083
00084
00085
00086 l=0
00087 water_new(l) = water(l) + dtime &
00088 * ( evap(l) - cond(l) &
00089 + diff_aux(l) &
00090 * ( cos_phi_cb2(l)*(water(l+1)-water(l)) &
00091 *dphi_inv(l+1) ) &
00092 )
00093
00094 do l=1, LMAX-1
00095 water_new(l) = water(l) + dtime &
00096 * ( evap(l) - cond(l) &
00097 + diff_aux(l) &
00098 * ( cos_phi_cb2(l)*(water(l+1)-water(l)) &
00099 *dphi_inv(l+1) &
00100 -cos_phi_cb1(l)*(water(l)-water(l-1)) &
00101 *dphi_inv(l) ) &
00102 )
00103 end do
00104
00105 l=LMAX
00106 water_new(l) = water(l) + dtime &
00107 * ( evap(l) - cond(l) &
00108 + diff_aux(l) &
00109 * ( -cos_phi_cb1(l)*(water(l)-water(l-1)) &
00110 *dphi_inv(l) ) &
00111 )
00112
00113 #elif SOLV_DIFF==2 /* Implicit scheme (Euler backward) */
00114
00115
00116
00117 l=0
00118 sle_a1(l) = 1.0_dp + dtime * diff_aux(l)*cos_phi_cb2(l)*dphi_inv(l+1)
00119 sle_a2(l) = - dtime * diff_aux(l)*cos_phi_cb2(l)*dphi_inv(l+1)
00120 sle_b(l) = water(l) + dtime * (evap(l)-cond(l))
00121
00122 do l=1, LMAX-1
00123 sle_a0(l) = - dtime * diff_aux(l)*cos_phi_cb1(l)*dphi_inv(l)
00124 sle_a1(l) = 1.0_dp + dtime * diff_aux(l) &
00125 * ( cos_phi_cb1(l)*dphi_inv(l) &
00126 +cos_phi_cb2(l)*dphi_inv(l+1) )
00127 sle_a2(l) = - dtime * diff_aux(l)*cos_phi_cb2(l)*dphi_inv(l+1)
00128 sle_b(l) = water(l) + dtime * (evap(l)-cond(l))
00129 end do
00130
00131 l=LMAX
00132 sle_a0(l) = - dtime * diff_aux(l)*cos_phi_cb1(l)*dphi_inv(l)
00133 sle_a1(l) = 1.0_dp + dtime * diff_aux(l)*cos_phi_cb1(l)*dphi_inv(l)
00134 sle_b(l) = water(l) + dtime * (evap(l)-cond(l))
00135
00136
00137
00138 call tri_sle(sle_a0, sle_a1, sle_a2, water_new, sle_b, LMAX)
00139
00140 #elif SOLV_DIFF==3 /* Instantaneous mixing (infinite diffusivity) */
00141
00142
00143
00144 do l=0, LMAX
00145 water_new(l) = water(l) + dtime * ( evap(l) - cond(l) )
00146 end do
00147
00148
00149
00150 water_mean = 0.0_dp
00151
00152 do l=0, LMAX
00153 water_mean = water_mean + 0.5_dp*water_new(l)*(sin_phi_cb2(l)-sin_phi_cb1(l))
00154 end do
00155
00156 water_new = water_mean
00157
00158 #else /* Wrong value of parameter SOLV_DIFF */
00159
00160 stop ' Wrong value of parameter SOLV_DIFF (must be 0, 1, 2 or 3)!'
00161
00162 #endif
00163
00164 end subroutine diff_trans
00165