32 - 能量平衡与叶温¶
热红外模拟需要地物的温度。温度可以直接给定（第 26 章），也可以由能量平衡根据光照、气温、湿度和风速求出。本章介绍后者。
求解什么¶
对每个参与求解的叶片或地面单元，LESS 求出使能量收支平衡的温度：
\[R_n = H + LE + G\]
- R_n：净辐射，等于吸收的短波和长波减去自身的长波发射；
- H：感热通量，由表面与空气的温差和交换阻力决定；
- LE：潜热通量，即蒸腾或蒸发，由气孔阻力和水汽压差决定；
- G：土壤热通量，叶片为 0。
长波辐射在叶片、地面和天空之间的交换在三维场景中计算，因此一片叶子的温度会影响它周围的叶子和地面。
所需属性¶
参与能量平衡的组件需要光学、热和生物物理三类属性：
import less
scene = less.Scene()
scene.size = 6.0
scene.terrain = less.Terrain(property=less.Lambertian(reflectance=0.15))
scene.terrain.set_property(less.ThermalProperty(temperature="air", emissivity=0.95))
scene.terrain.set_property(less.BiophysicalProperty(two_sided=False, stomata_side="none",
ground_heat_fraction=0.35))
scene.illumination = less.Illumination(
source=less.Sun(zenith=30, azimuth=150),
atmosphere=less.SimpleSpectralAtmosphere(turbidity=2.0))
plant = less.Object("plant")
plant.add_component("leaves", less.Mesh(less.asset_path("tutorial_tree.obj"), group="leaves"))
plant.set_property("leaves", less.Prospect(N=1.5, cab=45, car=10, cw=0.012, cm=0.006))
plant.set_property("leaves", less.ThermalProperty(temperature="air", emissivity=0.98))
plant.set_property("leaves", less.BiophysicalProperty(
leaf_width=0.05, stomata_side="bottom", stomatal_resistance=150.0))
scene.add(plant, positions=[(3, 3, 0)])
ThermalProperty.temperature 是求解的初始温度："air" 表示从气温开始，也可以给出开尔文温度。
BiophysicalProperty 描述表面如何与环境交换能量：
| 参数 | 含义 |
|---|---|
two_sided |
两面都参与辐射交换（叶片为真，地面总是单面） |
stomata_side |
气孔所在面："top"、"bottom"、"both" 或 "none"（不蒸腾） |
stomatal_resistance |
气孔阻力或土壤表面阻力（s m⁻¹） |
leaf_width |
叶片特征宽度（m），用于边界层阻力 r_a = 200 √(d/u) |
aerodynamic_resistance |
直接给出交换阻力（s m⁻¹），此时不用叶宽公式 |
ground_heat_fraction |
G 占净辐射的比例 |
ground_heat_flux |
常数 G（W m⁻²），与上一项二选一 |
solve |
为假时温度保持热属性中的值，不参与求解 |
常见设置：
| 对象 | two_sided |
stomata_side |
G |
|---|---|---|---|
| 叶片 | 是 | bottom |
0 |
| 土壤 | 否 | none 或 top |
比例或常数 |
| 树干 | 是 | none |
0 |
| 水面 | 否 | top |
0 |
微气象¶
climate = less.Microclimate(
air_temperature=25.0, # °C
humidity=60.0, # 相对湿度，%
wind_speed=2.0, # m s⁻¹
pressure=101325.0, # Pa
ca=400.0, # CO₂，ppm
)
downwelling_longwave（W m⁻²）给出天空向下的长波辐射，不给出时由气温和湿度估算。
求解¶
energy = scene.simulate(less.EnergyBalanceProcess(quality="preview"), microclimate=climate)
result = energy.result
print(result.summary())
quality 可以是 "preview"、"standard" 或 "high"，也可以直接设置 n_samples_per_face（每个三角形每面的长波采样数）、max_iterations、convergence_tol（K）和 energy_tolerance（W m⁻²）。
求解后，热红外成像会自动使用求得的温度：
thermal = scene.simulate(less.ThermalImager(
less.Orthographic(image_size=64), bands=[10600.0], quality=32), microclimate=climate)
thermal.to_brightness_temperature(unit="K").save("leaf_temperature.png")
结果¶
import numpy as np
area = np.asarray(result.primitive_areas)
temperature = np.asarray(result.temperatures) - 273.15
print("温度范围（°C）:", temperature.min(), temperature.max())
for name in ("R_n", "H", "LE", "G"):
print(name, "总量（W）:", float(np.sum(np.asarray(getattr(result, name)) * area)))
| 结果 | 含义 |
|---|---|
temperatures |
每个元素的温度（K） |
R_n、H、LE、G |
每单位面积的通量（W m⁻²） |
S_abs、L_abs、L_emit |
吸收的短波、吸收的长波、发射的长波 |
primitive_areas |
每个元素的面积（m²） |
temperature_sunlit、temperature_shaded、sunlit_fraction |
受光、阴影部分的温度及受光比例 |
residual、converged、n_iterations |
能量收支残差、是否收敛、迭代次数 |
检查能量收支是否闭合：
active = result.active_mask
print("最大残差（W m⁻²）:", abs(result.residual[active]).max())
print("收敛:", result.converged, "迭代次数:", result.n_iterations)
不收敛时会抛出 less.ConvergenceError，可以检查微气象、初始温度，或增加 max_iterations。
统计植被¶
统计植被的用法相同，结果按组件给出：组件内部分别求受光叶片和阴影叶片的温度，再按受光比例汇总。
canopy = less.Object("canopy")
canopy.add_component("crown", less.TurbidBoundary(less.asset_path("tutorial_tree.obj"), group="crown"))
canopy.set_property("crown", less.OpticalVegetation(
leaf=less.Prospect(cab=45), leaf_area_density=1.5, leaf_angle_distribution="spherical"))
canopy.set_property("crown", less.ThermalProperty(temperature="air", emissivity=0.98))
canopy.set_property("crown", less.BiophysicalProperty(leaf_width=0.05, stomata_side="bottom",
stomatal_resistance=150.0))
scene.add(canopy, positions=[(3, 3, 0)])
mixed = scene.simulate(less.EnergyBalanceProcess(quality="preview"), microclimate=climate)
print(mixed.result.summary())
同一场景可以同时包含三角网格叶片、统计植被和地形，它们在一次求解中相互交换辐射。
冠层地表层与大气稳定度¶
上面的叶片用叶宽计算边界层阻力。研究整片冠层与大气的交换时，可以使用 SCOPE 中的冠层地表层模型：冠层通过冠层阻力、土壤通过土壤阻力与参考高度处的空气交换感热和潜热，并根据整个场景的感热通量计算 Monin–Obukhov 长度，修正大气稳定度对阻力的影响。
field = less.Scene()
field.size = (1.0, 1.0)
field.repetitive = True
field.terrain = less.Terrain(property=less.Lambertian(reflectance=0.15))
field.terrain.set_property(less.ThermalProperty(temperature="air", emissivity=0.94))
field.terrain.set_property(less.BiophysicalProperty(two_sided=False, stomata_side="top",
stomatal_resistance=500.0,
ground_heat_fraction=0.35))
field.illumination = less.Illumination(
source=less.Sun(zenith=30, azimuth=180, irradiance=0.4), atmosphere=less.NoAtmosphere())
from pathlib import Path
Path("canopy_box.obj").write_text(
"g canopy\n"
+ "".join(f"v {x} {y} {z}\n" for x, y, z in [(0, 0, 0), (-1, 0, 0), (-1, 2, 0), (0, 2, 0),
(0, 0, -1), (-1, 0, -1), (-1, 2, -1), (0, 2, -1)])
+ "f 1 4 3\nf 1 3 2\nf 5 6 7\nf 5 7 8\nf 1 5 8\nf 1 8 4\n"
+ "f 2 3 7\nf 2 7 6\nf 1 2 6\nf 1 6 5\nf 4 8 7\nf 4 7 3\n")
crop = less.Object("crop")
crop.add_component("leaves", less.TurbidBoundary("canopy_box.obj", group="canopy"))
crop.set_property("leaves", less.OpticalVegetation(
leaf=less.Prospect(cab=45), leaf_area_density=1.5, leaf_angle_distribution="spherical"))
crop.set_property("leaves", less.ThermalProperty(temperature="air", emissivity=0.98))
crop.set_property("leaves", less.BiophysicalProperty(leaf_width=0.1, stomata_side="bottom",
stomatal_resistance=600.0))
field.add(crop, positions=[(0.0, 0.0, 0.0)])
layer = less.SurfaceLayer(canopy_height=2.0, leaf_area_index=3.0, measurement_height=5.0)
calm = less.Microclimate(air_temperature=25.0, humidity=50.0, wind_speed=0.5,
downwelling_longwave=350.0, surface_layer=layer)
coupled = field.simulate(less.EnergyBalanceProcess(quality="preview", max_iterations=100),
microclimate=calm)
print(coupled.metadata["surface_layer"])
SurfaceLayer 的主要参数：
| 参数 | 含义 |
|---|---|
canopy_height、leaf_area_index |
冠层高度（m）和叶面积指数 |
measurement_height |
风速和气温的测量高度（m），高于冠层 |
roughness_length、displacement |
粗糙度长度和零平面位移（m），不给出时按 SCOPE 由冠层参数计算 |
drag_coefficient、soil_boundary_resistance、within_canopy_resistance |
叶片拖曳系数、土壤边界层阻力、冠层内阻力 |
stability |
是否考虑大气稳定度，默认是 |
使用地表层时，叶片和土壤不要再单独给出 aerodynamic_resistance。coupled.metadata["surface_layer"] 给出最终的 Obukhov 长度、摩擦速度、冠层阻力、土壤阻力和场景的感热通量。白天冠层比空气热时大气不稳定，阻力变小、叶温降低；风速很小时这一影响可达 1–2 °C。由于 Obukhov 长度要与温度一起迭代，max_iterations 宜设为 60–100。
何时需要重新求解¶
修改光照、微气象、属性、几何或实例后，原来的结果失效，需要重新求解。把过期的结果交给传感器时，LESS 会报错。