from math import log, exp
def _euclidean_dist(vec1, vec2):
euc_sum = 0
for elem1, elem2 in zip(vec1, vec2):
euc_sum += abs(elem1 - elem2)
return euc_sum
[docs]def vapor_pressure(temperature, rel_humidity):
return rel_humidity * 611 * exp(2.5e6 / 461.5 * (1 / 273.15 - 1 / temperature))
[docs]def mixed_saturated_parcel(args, initial_args):
"""
Solve iteratively to find the temperature and vapor pressure of a supersaturated parcel after condensation.
Args:
args (list): A list containing temperature in K, pressure in Pa, and vapor pressure in Pa for this iteration.
initial_args (list): A list containing temperature in K, pressure in Pa, and vapor pressure in Pa from the
initial iteration.
Returns:
The next guess at the temperature, pressure, and vapor pressure for a supersaturated parcel after condensation.
"""
mixed_temperature, pressure, mixed_vapor_pres = args
init_mixed_temperature, init_pressure, init_mixed_vapor_pres = initial_args
mixed_saturation_temp = init_mixed_temperature + (2.5e6 * 0.622) / (pressure * 1005.) * \
(init_mixed_vapor_pres - vapor_pressure(mixed_temperature, 1.))
mixed_saturation_vapor_pres = vapor_pressure(mixed_temperature, 1.)
solution = [(mixed_saturation_temp + mixed_temperature) / 2, pressure,
(mixed_saturation_vapor_pres + mixed_vapor_pres) / 2]
return solution
[docs]def condensation_temp(args, initial_args):
"""
Solve iteratively to find the temperature and pressure at which a parcel will condense if lifted adiabatically.
Args:
args (list): A list containing temperature in K, pressure in kPa, and mixing ratio in g/g for this iteration.
initial_args (list): A list containing temperature in K, pressure in kPa, and mixing ratio in g/g from the
initial iteration.
Returns:
The next guess at the temperature, pressure, and mixing ratio for the condensation temperature of a parcel
lifted adiabatically.
"""
A = 2.53e8
B = 5.42e3
temperature, pressure, mixing_ratio = args
init_temperature, init_pressure, init_mixing_ratio = initial_args
condens_pres = init_pressure * (temperature / init_temperature) ** (7. / 2.)
condens_temp = B / log((0.622 * A) / (mixing_ratio * init_pressure) * (init_temperature / temperature) ** (7. / 2.))
return [condens_temp, condens_pres, mixing_ratio]
[docs]def wet_bulb_temp(args, initial_args):
"""
Solve iteratively to find the wet-bulb temperature for a parcel
Args:
args (list): A list containing temperature in K, pressure in kPa, and mixing ratio in g/g for this iteration.
initial_args (list): A list containing temperature in K, pressure in kPa, and mixing ratio in g/g from the
initial iteration.
Returns:
The next guess at temperature, pressure, and mixing ratio for the wet-bulb temperature of a parcel.
"""
A = 2.53e8
B = 5.42e3
temperature, pressure, mixing_ratio = args
init_temperature, init_pressure, init_mixing_ratio = initial_args
wet_bulb = init_temperature - (2.5e6 / 1005.) * (0.622 / pressure * A * exp(-B / temperature) - mixing_ratio)
return [wet_bulb, pressure, mixing_ratio]
[docs]def solve_iteratively(func, args, first_guess=None, tolerance=1e-5):
"""
Solve a function iteratively.
Args:
func (function): The function to solve. It should take two lists of arguments, the first for this iteration, the
second from the initial iteration.
args (list): Initial set of arguments, representing the initial state of the parcel.
first_guess (list): A first guess at the solution. Optional, defaults to the initial state if not given.
tolerance (float): Tolerance on the change between successive solutions. Optional, defaults to 1e-5 if not
given.
Returns:
An iterative solution to the provided function.
"""
if first_guess is None:
first_guess = args
counter = 0
solution = func(first_guess, args)
last_solution = [s + tolerance * 2 for s in solution]
while _euclidean_dist(solution, last_solution) > tolerance:
last_solution = solution
solution = func(last_solution, args)
counter += 1
return solution
[docs]def iterative_main():
print "Condensation Temperature:"
solution = solve_iteratively(condensation_temp, [299.87, 95.25, 17.005e-3])
print "Temperature:", solution[0] - 273.15, "C"
print "Pressure:", solution[1] * 10, "hPa"
print "Mixing Ratio:", solution[2] * 1000, "g/kg"
print
print "Wet Bulb Temperature:"
solution = solve_iteratively(wet_bulb_temp, [268.15, 80., 2.13e-3])
print "Temperature:", solution[0] - 273.15, "C"
print "Pressure:", solution[1] * 10, "hPa"
print "Mixing Ratio:", solution[2] * 1000, "g/kg"
print
print "Wet Bulb Potential Temperature:"
solution = solve_iteratively(wet_bulb_temp, [285.7, 100., 2.13e-3], [280., 100., 2.13e-3])
print "Temperature:", solution[0], "K"
print "Pressure:", solution[1] * 10, "hPa"
print "Mixing Ratio:", solution[2] * 1000, "g/kg"
print
print
print "Mixed Saturation Temperature:"
mixed_temperature = (303.15 + 275.15) / 2
mixed_vapor_pressure = (vapor_pressure(303.15, .9) + vapor_pressure(275.15, .8)) / 2
solution = solve_iteratively(mixed_saturated_parcel, [mixed_temperature, 100000., mixed_vapor_pressure])
print "Temperature:", solution[0], "K"
print "Pressure:", solution[1] / 100, "hPa"
print "Vapor Pressure:", solution[2] / 100, "hPa"
return
if __name__ == "__main__":
iterative_main()