Goal: In this example, check whether the Meteorology component produces output for land_surface_net-longwave-radiation__energy_flux (Qn_LW internally) when the model state is updated, given scalar inputs for temperature and emissivity.
Begin with an import and some magic:
In [ ]:
%matplotlib inline
import numpy as np
Import the BMI Meteorology component, create an instance, and initialize the component:
In [ ]:
from topoflow.components.met_base import met_component
m = met_component()
m.initialize('./input/meteorology.cfg')
There's no input value for Qn_LW; rather, it's calculated from temperature and emissivity values. It's initialized to an np.zeros array with the dimensions of the grid (line 951 of met_base.py):
In [ ]:
Q = m.get_value('land_surface_net-longwave-radiation__energy_flux') # `Qn_LW` internally
print type(Q)
print Q.size
print Q.shape
Q
Advance the model by one time step:
In [ ]:
m.update()
print '\nCurrent time: {} s'.format(m.get_current_time())
Check the calculated value of Qn_LW:
In [ ]:
Q = m.get_value('land_surface_net-longwave-radiation__energy_flux')
print type(Q)
print Q.size
print Q.shape
Q
The initial Qn_LW array has been converted into a scalar.
The Qn_LW calculation starts on line 1631 of met_base.py:
#--------------------------------
# Compute Qn_LW for this time
#--------------------------------
T_air_K = self.T_air + self.C_to_K
T_surf_K = self.T_surf + self.C_to_K
LW_in = self.em_air * self.sigma * (T_air_K)** 4.0
LW_out = self.em_surf * self.sigma * (T_surf_K)** 4.0
LW_out = LW_out + ((1.0 - self.em_surf) * LW_in)
self.Qn_LW = (LW_in - LW_out) # [W m-2]
The temperature and emissivity values are set to scalars in the cfg file, so Qn_LW is also a scalar.
Result: A scalar value of Qn_LW is calculated from scalar values of temperature and emissivity, which seems legit.
Import the Babel-wrapped Meteorology component, create an instance, and initialize the component:
In [ ]:
from cmt.components import Meteorology
bm = Meteorology()
bm.initialize('./input/meteorology.cfg')
As above, Qn_LW is initialized to an np.zeros array with the dimensions of the grid. However, the Bocca implementation flattens the array:
In [ ]:
bQ = bm.get_value('land_surface_net-longwave-radiation__energy_flux') # `Qn_LW` internally
print type(bQ)
print bQ.size
print bQ.shape
bQ
Get the model time step (which has been broadcast over the grid):
In [ ]:
time_step = bm.get_value('model__time_step') # `dt` internally
time_step = time_step.max()
print time_step
Advance the model by one time step:
In [ ]:
bm.update(time_step)
print '\nCurrent time: {} s'.format(bm.get_current_time())
Check the calculated value of Qn_LW:
In [ ]:
Q = bm.get_value('land_surface_net-longwave-radiation__energy_flux')
print type(Q)
print Q.size
print Q.shape
Q
The calculated value of Qn_LW matches that of the BMI case above, but it has been broadcast over the grid. This is what we expect.
Advance the model by a second time step:
In [ ]:
bm.update(2*time_step)
print '\nCurrent time: {} s'.format(bm.get_current_time())
Again, check the calculated value of Qn_LW:
In [ ]:
Q = bm.get_value('land_surface_net-longwave-radiation__energy_flux')
print type(Q)
print Q.size
print Q.shape
Q
There are no changes to the value of Qn_LW, which is expected, because there are no changes to the temperature or emissivity values.
Result: An array of Qn_LW values, matching the shape of its associated grid, are calculated from scalar values of temperature and emissivity, which is what we expect!