ASABE Logo

Article Request Page ASABE Journal Article

Dairy Cow Thermal Balance Model During Heat Stress: Part 1. Model Development

Chad R. Nelson1, Kevin A. Janni1,*


Published in Journal of the ASABE 66(2): 441-460 (doi: 10.13031/ja.15190). Copyright 2023 American Society of Agricultural and Biological Engineers.


1Bioproducts and Biosystems Engineering, University of Minnesota, St. Paul, Minnesota, USA.

*Correspondence: kjanni@umn.edu

The authors have paid for open access for this article. This work is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License https://creative commons.org/licenses/by-nc-nd/4.0/

Submitted for review on 16 May 2022 as manuscript number PAFS 15190; approved for publication as a Research Article by Associate Editor Dr. Lingjuan Wang-Li and Community Editor Dr. Jun Zhu of the Plant, Animal, & Facility Systems Community of ASABE on 7 December 2022.

Highlights

Abstract. Dairy cow heat stress impacts cow well-being, reduces milk yield, and leads to economic losses. Understanding heat stress mechanics supports ongoing and future efforts to mitigate heat stress. The purpose of this project was to modify a steady-state heat transfer model developed by McGovern and Bruce (2000) by incorporating work by Berman (2005), McArthur (1987), Turnpenny et al. (2000a,b), Thompson et al. (2014), Gwadera et al. (2017), two new empirical relations for tissue insulation and sweat rate, and a new solution method that allowed for overlapping changes in heat exchange. The modified model describes heat exchange between a lactating cow and the environment through respiration, convection, sweating, and shortwave and longwave radiation. This article describes the process-based model equations, compares results from the two new empirical relations used to published work, and presents the inputs and results for a cow on pasture in sunlight. The modified model, which can be solved with a spreadsheet, provides insight into factors and processes that affect lactating cow heat exchange. A companion paper compares the modified model results with published average measured body temperatures, respiration rates, and skin temperatures and unpublished body temperature data for cows on pasture in the sunshine.

Keywords. Body temperature, Dairy, Heat stress, Lactating cow, Respiration rate, Thermal balance model.

Dairy cows are homeothermic, which means that they make behavioral and physiological adjustments to maintain their core body temperature within a narrow range. The behavioral and physiological adjustments that cows make can modify metabolic heat production, reduce heat gains from the environment, or increase heat losses to the environment. Animal caregivers can better help cows deal with heat stress when they understand and can quantitatively describe heat exchange between lactating cows and the surrounding environment.

Dairy cow heat stress is an active research topic. Several studies have reported on heat stress causes and impacts on dairy cows (West, 2003; Berman, 2005; Gebremedhin et al., 2010). Shoshani and Hetzroni (2013) assessed the impact of naturally ventilated barn design factors on high-producing cows in hot weather using a threshold temperature that was a function of ambient temperature, relative humidity, and wind velocity. Chen et al. (2015) reported physiological responses (i.e., body temperatures, respiration rates, and skin temperatures) to water sprayed on lactating Holstein cows in hot weather. Gorczya and Gebremedhin (2020) used individual cow data collected by Chen et al. (2015) to develop machine learning algorithms to rank heat stressor factors. Wang et al. (2018a, 2019) developed an empirical equivalent temperature index for cattle that calculates a temperature that considers dry-bulb temperature, relative humidity, air velocity, and solar radiation. Li et al. (2020) presented relations between Holstein dairy cows' mean rectal temperature and mean respiration rateand developed an empirical equation for finding the mean respiration rate based on air temperature, relative humidity, wind speed, and milk yield.

Work has also been done to assess the expected impact of heat stress on dairy cows. St-Pierre et al. (2003) and Abuajamieh et al. (2013) reported on economic losses due to heat stress on dairy cows. Additional studies by Key et al. (2014) and Mauger et al. (2015) described expected climate change impacts on U.S. milk production.

Numerous researchers have developed models that describe how dairy and beef cows are expected to balance metabolic heat production, shortwave radiation, respiration, convection, longwave radiation, and evaporation to maintain body temperatures in acceptable ranges (McArthur, 1987; McGovern and Bruce, 2000; Turnpenny et al., 2000a; Thompson et al., 2014). These researchers neglected conductive heat loss to the surrounding environment, which corresponds to standing cows (McArthur, 1987; McGovern and Bruce, 2000; Turnpenny et al., 2000a; Thompson et al., 2014).

Mondaca et al. (2013) developed a transient conjugate heat and mass transfer model of conductive, longwave radiation, convective, and evaporative heat exchange for a cow lying down on a cooling surface. They assumed that the cow had a constant body core temperature, did not include respiration heat exchange, and used a control volume-based approach to model metabolic heat generation and transfer within the cows body to maintain a constant core body temperature (Mondaca et al., 2013). The modeled cows metabolic heat production was not related to either cow mass or daily milk yield (Mondaca et al., 2013). The model was validated using data from a fiberglass calf with a heater placed in the body to mimic actual metabolic heat production (Mondaca et al., 2013).

Wang et al. (2018a, 2019) developed an empirical temperature index based on air temperature, relative humidity, air velocity, and solar radiation for high-producing Holstein cows using a data set with over 845 observations. The resulting Equivalent Temperature Index for Cattle (ETIC) did not consider cow mass or daily milk yield. Wang et al. (2018a) reported that the ETIC was able to account for approximately half of the variance in observed respiration rates.

Wang et al. (2018b) used computational fluid dynamics to assess convective heat transfer from cows at five air velocities either standing or reclining and exposed to either cross, front, or downward airflows. Body heat exchange through conduction, radiation, and evaporation was not modeled (Wang et al., 2018b). Zhou et al. (2019) used computational fluid dynamics and the empirical ETIC developed by Wang et al. (2018a, 2019) to assess baffle locations in cross-ventilated dairy barns. Zhou et al. (2019) modeled heat exchange through a cows skin and hair coat via convection, radiation, and evaporation of sweat assuming a constant body temperature but the study neglected to include respiration heat exchange in the model. Both Wang et al. (2018b) and Zhou et al. (2019) conducted steady-state analyses.

Existing models describing heat exchange between lactating cows and the surrounding environment include numerous assumptions and simplifications. Model assumptions and simplifications are acceptable if noted, and the model results provide useful information.

Engineers use models to help them design systems, processes, and facilities that are effective and economical. Asteady-state process-based model describing heat exchange between a lactating cow and the environment through respiration, convection, sweating, and shortwave and longwave radiation heat exchange that produced useful results and could be solved using a spreadsheet will help engineers assess cooling system design factors (e.g., air velocity and evaporative cooling) for different environmental conditions (e.g., hot and dry, hot and humid, sunny, and shaded).

This project aimed to modify a steady-state heat transfer model that McGovern and Bruce (2000) developed to describe a cows heat exchange during heat stress and that can be solved using a spreadsheet. The McGovern and Bruce (2010) model was modified by incorporating work by Berman (2005), McArthur (1987), Turnpenny et al. (2000a, b), Thompson et al. (2014), and Gwadera et al. (2017). The modified model incorporated new empirical relations for finding tissue insulation and sweat rate and used a solution procedure that allowed all heat exchange methods to adjust concurrently rather than sequentially, as done by McGovern and Bruce (2010).

This is the first of two articles. This article describes the modified models development and the equations used. Select intermediate values are compared with published studies. An example set of inputs and outputs is presented. Asecond article compares the modified model respiration rate, sweating rate, and body temperatures to published data for heat-stressed lactating cows from Gebremedhin et al. (2010) and Chen et al. (2015). The second article also compares unpublished body temperatures of grazing cows in the sun at the University of Minnesota West Central Research and Outreach Center (WCROC) that were part of a shade study (Sharpe et al., 2021). The results of the modified model were compared to averaged data corresponding to steady-state conditions.

Model Development

McGovern and Bruce (2000) developed a steady-state energy balance model that described the heat exchange between a standing cow and the surrounding environment, which described a cows metabolic heat generation (M) and heat exchange to the surrounding environment. The heat exchanges were respiratory latent and sensible heat loss (Er), evaporative heat loss from the skin from sweating (Es), shortwave solar radiation heat gain (Rn), and convective (C) and longwave radiative (LWn) heat exchanges (gain or loss), and body heat storage (Gb) (fig. 1). Conductive heat exchange was neglected.

The heat balance model by McGovern and Bruce (2000) achieved a steady-state heat balance by adjusting tissue insulation, sweat rate, and respiration rate in sequential order. The sequence to increase heat exchange from the cow to the environment was: (1) reduce tissue insulation (Ib) to a minimum to account for vasodilation and increased conduction through the cows skin, (2) increase sweat rate to a maximum physiological rate, and (3) increase respiration rate to a maximum limit. If the metabolic heat could not be dissipated to the cow's surroundings as a result of those changes, the excess heat was stored in the body (an unsteady-state condition) and raised the cow's body temperature.

The sequential initiation of heat loss regulating mechanisms used by McGovern and Bruce (2000) (i.e., reduced tissue insulation, increased sweat rate, and increased respiration) was not used. The modified model allowed overlapping changes in heat exchange in response to heat stress, as described by McArthur (1987) and Thompson et al. (2011). Tissue insulation, sweat rate, and body temperature were calculated using functions (i.e., eqs. 26, 38, and 61) that allowed all methods of heat exchange to be adjusted at the same time in the modified model described here.

The modified model also incorporated modifications and empirical equations for heat production, coat insulation, tissue insulation, body surface area, and respiratory cooling using exhaled air temperature and tidal volume developed by Berman (2003, 2004, 2005).

Figure 1. Schematic diagram of modeled heat flows following McGovern and Bruce (2000).

Variables, inputs, and parameters with their corresponding values used in the modified model are listed and defined in table 1. The value column indicates whether the number was calculated using the equation listed in the equation column, was a model input, or was the constant given.

Model Schematic

Figure 2 is a schematic outlining the thermal processes and equations used in the modified spreadsheet model. The equations describe the heat exchanges shown in figure 1. The thermal balance and heat exchange equations are described in the following sections. A spreadsheet add-in program, Solver, was used to solve the system of equations using the given constants in table 1 and the inputs used. Arrows in figure 2 indicate the direction of information flow to solve the model. Bidirectional arrows indicate where information flows in both directions using Solver. The modified spreadsheet model adjusted the assumed inspired volumetric respiration rate, Vresp, and other variables to find intermediate calculated values and the cow body temperature, Tb, skin temperature, Ts, coat temperature, Tc, and respiration rate, rr, that balanced the heat balance given by equation 4.

Overall Heat Balance and Exchanges

The modified model describes heat exchange through the body core, skin surface, and coat surface using heat balance equations from McGovern and Bruce (2000) (fig. 1). The flux units used were watts per square meter of cow surface area (Wm-2). The body core equation (eq. 1) (McGovern and Bruce, 2000) balances metabolic heat production, M, conductive heat exchange from the body to the skin surface, Qb, sensible and latent heat exchange due to respiration, Qr,n, and body heat storage, Gb:

(1)

where

Qb = Heat flux from body core to skin surface (W m-2)

M = Metabolic heat production flux (W m-2)

Qr,n = Net latent and sensible heat exchange flux by respiration (W m-2)

Gb= Body heat storage flux (W m-2).

The skin surface equation (eq. 2) (McGovern and Bruce, 2000) balances conductive heat transfer from the body to the skin surface, Qb, and the evaporative heat exchange from the skin through the coat, Es, and sensible heat exchange from the skin surface to the coat surface, Qs,c:

(2)

where

Qs,c= Heat flux from skin surface to coat surface (W m-2)

Qb = Heat flux from body core to skin surface (W m-2)

Es = Latent heat loss from skin (W m-2)

The coat surface equation (eq. 3) (McGovern and Bruce, 2000) balances the sensible heat exchange from the coat surface to the surrounding environment, Qc,e, the convective heat exchange between the coat surface and the surrounding air, C, longwave radiation between the coat surface and surrounding surfaces, LWn, and shortwave radiation from the sun if in sunshine or under solar lamps, Rn:

(3)

where

Qc,e = Outbound heat flux from coat surface to the surrounding environment (W m-2)

C = Convective heat flux from coat surface to the surrounding air (W m-2)

Table 1. Variable symbol, description, value or source, units, and equation where calculated.
SymbolDescriptionValue[a]UnitsEquation[b],[c]
a1Parameter used in equation 26Calculatedbpm m2 K W-127
a2Parameter in sweat rate function equation 38CalculatedC39
AAnimal surface areaCalculatedm25
AePercent animal surface exposed to airInputpercent
Ah / AHorizontal area absorbing shortwave radiation per unit area of cowCalculatedm2 m-257 or 58
bParameter used in equation 26Calculatedm2 K W-129
CConvective heat flux from coat surface to airCalculatedW m-243
cbHeat capacity of the cows body3,400J kg-1 K-1
CICoat insulation for Holsteins0.30C m2 d Mcal-1 mm-1
CCfCloud cover factorCalculatedDecimal53 or 54
CpaSpecific heat of airCalculatedJ kg-1 K-1P
CSfClear sky fractionInputDecimal
CvaVolumetric specific heat of airCalculatedJ m-3 K-1P
dtEffective cow trunk diameterCalculatedm6
DDiffusivity of water in air adjusted for temperature and pressureCalculatedm2 s-134
EsLatent heat loss from skinCalculatedW m-240
EvapEvaporation flux for finding ground temperature52W m-2
fLongwave exposure factor1Dimensionless
Fat%Butterfat contentInput%
GbBody heat storage fluxCalculatedW m-223
GrcoatGrashof number at coat and air interfaceCalculatedDimensionless48
hConvective coefficient for finding ground temperature33W m-2 K-1
hcCoat conductanceCalculatedW m-2 C42
hmMid cow heightInputm
HPlactLactation heat productionCalculatedMcal d-19
HPmaintMaintenance heat productionCalculatedMcal d-18
IbTissue resistanceCalculatedm2 K W-126
Ib,maxMaximum tissue resistanceCalculatedm2 K W-125
Ib,minMinimum tissue resistance0.016m2 K W-1
kThermal conductivity of air using average of air and coat surface temperaturesCalculatedW m-1 K-143
ltEffective cow trunk lengthCalculatedm7
lcCoat thickness3.0mm
ldCoat depth = lc / (1000 mm/m)Calculatedm35
?lwCoat depth decrease due to windCalculatedm35
LWnLongwave radiant heat flux exchange between cows coat and surrounding surfacesCalculatedW m-250
LWgGround longwave radiation flux for finding Tground63W m-2
MMetabolic heat production fluxCalculatedW m-24, 13
Mmaint+lactMaintenance and lactation heat production fluxCalculatedW m-210
?MbMetabolic heat production adjustment for body temperature changeCalculatedW m-211
?MrrMetabolic heat production adjustment for respiration rate changeCalculatedW m-212
NuNusselt numberCalculatedDimensionless44, 45
PAtmospheric pressureInputPa
PwPartial vapor pressureCalculatedPaP
PwcoatAverage of Pws,Ts and PwCalculatedPa49
PwsSaturation vapor pressureCalculatedPaP
Pws,TsSaturation vapor press at skin surface temperatureCalculatedPaP
Qlw radLongwave radiation flux between coat and surroundings when outsideCalculatedW m-251
QbHeat flux from body core to skin surfaceCalculatedW m-21
Qc,eHeat flux from coat surface to surrounding environmentCalculatedW m-23
Qs,cHeat flux from skin surface to coat surfaceCalculatedW m-22, 41
Qr,minNet latent and sensible heat exchange by respiration at minimum respiration for finding sweat rateCalculatedW m-214
Qr,nNet latent and sensible heat exchange by respirationCalculatedW m-214
rcoatCoat reflectivity0.3Dimensionless
RecoatReynolds number at the coat and air interfaceCalculatedDimensionless46
rvResistance to vapor transfer in the coatCalculateds m-133
rrRespiration rateInterpolatedbpm17 - 22
rr50%respRespiration rate at 50% of first phase respiratory coolingCalculatedbpm28
rr50%totrespRespiration rate at 50% total respiratory cooling capacityCalculatedbpm37
rrmaxRespiration rate at maximum respiration rate limitCalculatedbpm18
rrsevereRespiration rate at severe heat limitCalculatedbpm19
rrwarmRespiration rate at warm conditions limitCalculatedbpm17
RnRadiant heat flux incoming from sun or solar lampsCalculatedW m-260
RHaAmbient air relative humidityCalculated%P

    [a]Value indicates whether the number was calculated using the equation listed, was a model input or the constant given.

    [b]S = Solar equations from ASHRAE (2021).

    [c]P = Psychrometric equations from ASHRAE (2021).

Table 1 (continued). Variable symbol, description, value or source, units and equation where calculated.
SymbolDescriptionValue[a]UnitsEquation[b],[c]
SnhSolar irradiance on a horizontal surface. Calculated using equation 59 or inputted.Input or
calculated
W m-259
S0Extraterrestrial radiant flux at 21st of the monthInputW m-2S
ssSun shading factorInputDecimal
t1Heat storage time 3h
t2Apparent solar time (0 to 24)Inputh
TaAmbient air dry-bulb temperatureInputC
Ta,KAmbient air dry-bulb temperatureInputK
TbCows body temperatureCalculatedC61
TcCoat surface temperatureCalculatedC41
TcoataveAverage of coat surface and air temperatureAveragedC
TdpAmbient air dew-point temperatureInputC
TgroundGround surface temperatureCalculatedK56
TexExhaled air temperatureCalculatedC15
TkAbsolute temperature (i.e., ambient air, exhaust air or skin)Input or
calculated
K
TMRTSurrounding surfaces mean radiant temperatureInputC
TNNormal cow body temperature38.0C
TsSkin temperatureCalculatedC24
Tsky,KSurrounding sky temperatureCalculatedK55
TsurroundingTemperature of the surroundings which was an average of
the ground temperature and the sky temperature
CalculatedK53,54
Ts,LCTLower critical skin temperature for finding sweat rate functionCalculatedC24
Ts, 50%Qr,nSkin temperature at 50% respiratory coolingCalculatedC24
Ambient air virtual temperatureCalculatedK16
Exhaled air virtual temperature at saturated vapor pressureCalculatedK16
Air virtual temperature used in equation 33CalculatedK16
Skin virtual temperature used in equation 33CalculatedK16
ucowAir speed at cow levelInput or
calculated
m s-1
uWWind speedInputm s-1
VrespVolumetric respiration rateVariableL min-1
Vt,maxTidal volume at maximum respirationCalculatedL breath-121
Vt,severeMaximal tidal volume based on body massCalculatedL breath-120
Vt,warmTidal volume at warm conditionsCalculatedL breath-122
Vt,50%respTidal volume at 50% total respiratory cooling capacityCalculatedL breath-136
WCow massInputkg
WaAmbient air saturation humidity ratioCalculateddecimalP
WsSaturation humidity ratio at exhaled air temperatureCalculateddecimalP
xAnimal size ratio = (2lt) / dtCalculatedDimensionless
YmilkMilk yieldInputkg day-1
aaAtmospheric dispersion factor0.1Dimensionless
Solar altitude CalculateddegS
dSolar declination angleCalculateddegS
ecoatEmissivity of cows coat0.98Dimensionless
egroundGround emissivity0.9Dimensionless
eskyCloudy sky emissivityCalculatedDimensionless52
?cowAzimuth of cow based on orientationInputdeg
?sSolar azimuthCalculateddegS
?TAngular difference between the solar azimuth and azimuth of the cowCalculateddegS
?Psychrometric constant CalculatedPa K-132
?aveAverage psychrometric constant value using ambient and skin temperaturesAveragedPa K-1
Sweat rate as function of skin temperature, TsCalculatedg m-2 h-138
Minimum sweat rate14.4g m-2 h-1
Maximum potential environmental sweat rateCalculatedg m-2 h-131
Maximum physiological sweat rate660g m-2 h-1
Sweat rate at 50% total respiratory cooling capacityCalculatedg m-2 h-130
?Latent heat of vaporizationCalculatedJ kg-1P
?aveAverage latent heat of vaporization using ambient and skin temperaturesAveragedJ kg-1P
[a] Value indicates whether the number was calculated using the equation listed, was a model input or the constant given.
[b]S = Solar equations from ASHRAE (2021).
[c]P = Psychrometric equations from ASHRAE (2021).
Figure 2. Schematic of model thermal processes, equations, inputs, and constants used in the spreadsheet model.

LWn = Longwave radiation heat flux exchange with surrounding surfaces (W m-2)

Rn= Shortwave radiation heat flux from the sun or solar lamps (W m-2)

When fluxes Qs,c and Qc,e are equal, equations 1, 2, and 3 can be combined and reorganized to describe the overall heat balance of a cow (eq. 4) with the surrounding environment. Equation 4 neglects conductive heat exchange:

(4)

Calculated Cow Dimensions

The cow dimensions needed in the modified model were based on an input body mass. A total surface area equation (eq. 5) used by Thompson et al. (2014) from Brody (1945) was used:

(5)

Table 1 (continued). Variable symbol, description, value or source, units and equation where calculated.
SymbolDescriptionValue[a]UnitsEquation[b],[c]
?sLatent heat of vaporization at skin temperature CalculatedJ kg-1P
OCloud cover factor (okta 0 to 8)Input
fAlbedo0.25Dimensionless
?aAir densityCalculatedkg m-3P
vcoat airKinematic viscosity at the coat and air interfaceCalculatedm2 s-147
sStefan-Boltzmann constant5.67 x 10-8W m-2 K-4

    [a] Value indicates whether the number was calculated using the equation listed, was a model input or a constant.

    [b]S = Solar equations from ASHRAE (2021).

    [c]P = Psychrometric equations from ASHRAE (2021).

where

A = Calculated cow surface area (m2)

W = Cow mass (kg).

Expressions from McGovern and Bruce (2000) were used to find an effective cow trunk diameter, dt (eq. 6), and lengthlt (eq. 7), assuming that the total cow surface area could be represented by a cylinder:

(6)

(7)

where

dt = Effective cow trunk diameter (m)

lt = Effective cow length (m).

Exposed Surface Area

The heat exchange flux terms used in the modified model were given as W m-2 of cow surface area. The percent of cow surface area assumed to be available for convective, longwave, and skin evaporative heat exchange was adjustable, Ae. An Ae value of 100% represented a standing cow with 100% of her surface area available for heat exchange.

Metabolic Heat Production

Heat production was calculated using procedures developed by Berman (2005). Maintenance heat production was found using equation 8:

(8)

where

HPmaint = Maintenance heat production (Mcal d-1).

Lactation heat production was found using milk yield and percent fat in the milk:

(9)

where

HPlact Lactation heat production (Mcal d-1)

Ymilk = Milk yield (kg d-1)

Fat% = Milk fat (%).

The maintenance and lactation heat production flux per m2 skin surface was found using equation 10:

(10)

where

Mmaint+lact = Maintenance and lactation heat production flux (W m-2)

Unit conversion

The modified model used equations from McGovern and Bruce (2000) to increase metabolic heat production for both increases in body temperature and respiration rate. The body temperature adjustment equation was originally from Blaxter (1967). The increase in metabolic heat production per unit m2 skin area from an increase in body temperature, Tb-TN, was found using equation 11:

(11)

where

?Mb = Adjustment in metabolic heat production due to change in body temperature (W m-2)

Tb = Cows body temperature (C)

TN = Normal cows body temperature (C).

The increase in metabolic heat production due to the respiration rate relation from McGovern and Bruce (2000) was calculated using equation 12:

(12)

where

?Mrr = Adjustment in metabolic heat production due to changing respiration rate (W m-2)

rr = Respiration rate in breaths per minute (bpm).

Total metabolic heat production flux was the sum of the maintenance and milk production heat production flux and the metabolic heat production flux adjustments for increased body temperature, ?Mb, and respiration rate, ?Mrr:

(13)

Cows commonly reduce dry matter intake and daily milk production during heat stress (St-Pierre et al., 2003; Abuajamieh et al., 2013). This process-based model does not adjust metabolic heat production and assumes no change in dry matter intake or milk production by the cow in response to heat stress.

Psychrometric and Solar Values

The modified model uses numerous psychrometric values (e.g., air density, specific volume, enthalpy, and partial pressure) for given or calculated dry-bulb and dew-point temperatures at different locations given in table 1. These psychrometric values were calculated using equations from ASHRAE (2021) for given or calculated dry-bulb and dew-point temperatures and a given atmospheric pressure.

Similarly, several solar values that depend on latitude, day of the year, and apparent solar hour (i.e., solar altitude and solar azimuth) were calculated using equations from ASHRAE (2021). These solar values are used to describe the shortwave radiation flux from the sun.

Respiration Heat Exchange

The respiratory heat flux was calculated using virtual temperatures and an equation from McGovern and Bruce (2000). The air density, ?a, specific heat, Cpa, and latent heat of vaporization, ?, values used in equation 14 were average values at the ambient, Ta, and exhaled, Tex, air temperatures:

(14)

where

Vrep = Volumetric respiration rate, (L min-1)

?a = Average air density of ambient air and exhaled air conditions, (kg m-3)

Cpa = Average specific heat of ambient air and exhaled air, (J kg-1 K-1)

= Exhaled air virtual temperature at saturated vapor pressure, (K)

= Ambient air virtual temperature (K)

? = Average latent heat of vaporization using ambient air and exhaled air temperatures (J kg-1)

Ws = Saturation humidity ratio at exhaled air temperature (decimal)

Wa = Ambient air humidity ratio (decimal)

The exhaled respiration air temperature was calculated using equation 15 from Berman (2005):

(15)

where

Ta = Ambient air temperature (C)

RHa = Ambient air relative humidity (percent)

Virtual temperatures were adjusted using either the air saturation vapor pressure, Pws, or partial vapor pressure, Pw, for either ambient air or exhaled air conditions. Virtual temperatures were calculated for corresponding temperatures using equation 16 from McGovern and Bruce (2000):

(16)

where

T* = Virtual temperature (K)

Tk = Corresponding temperature. Ambient air or exhaled air temperature for heat exchange by respiration (K)

Pw = Partial vapor pressure for ambient air or saturation vapor pressure (Pws) for exhaled air (Pa)

P = Atmospheric pressure (Pa).

Respiration Rate, Tidal Volume, and Inspired Volumetric Rate Relations

Respiration rate in breaths per minute (bpm) times the tidal volume (L per breath) gives the inspired volumetric rate by respiration in L min-1. Respiration rate relations by McGovern and Bruce (2000) and the maximal tidal volume relationship by Berman (2005) based on cow mass were used to describe the change in inspired volumetric rate as the cow went from warm conditions to maximum respiration to severe heat.

Respiration Rate Relations

Three respiration rate relations (eqs. 17, 18, and 19) from McGovern and Bruce (2000), depend on cow mass, W:

Warm conditions

(17)

Maximum respiration rate

(18)

Severe heat

(19)

where

rrwarm = Respiration rate at warm conditions limit (bpm)

rrmax = Respiration rate at maximum respiration rate limit (bpm)

rrsevere = Respiration rate at severe heat limit (bpm).

Tidal Volume Relations

Three maximal tidal volumes, like the respiration rate relations (eqs. 17, 18 and 19), depend on cow mass. The severe heat tidal volume, Vt,severe, relation for Holstein cows is from Berman (2005) (eq. 20). The Vt,severe result is used to find the tidal volumes at maximum respiration (eq. 21) (McGovern and Bruce, 2000). In turn, the Vt,max is used to calculate the tidal volume at warm conditions (eq. 22) given by McGovern and Bruce (2000):

Severe heat tidal volume:

(20)

Tidal volume at maximum respiration:

(21)

Tidal volume at warm conditions:

(22)

where

Vt,severe = Maximal tidal volume based on body mass from Berman (2005) (L breath-1)

Vt,max = Tidal volume at maximum respiration (L breath-1)

Vt,warm = Tidal volume at warm conditions (L breath-1).

Intermediate respiration rates and tidal volumes were found assuming linear relations between either the warm condition limit and maximum respiration limit or between the maximum respiration limit and severe heat limit. Respiration rates were described as rapid shallow panting when they were between the warm condition limit and maximum respiration limit and as deep slower panting when between the maximum respiration limit and severe heat limit (McGovern and Bruce, 2000).

The maximum respiration rate found using equation 18 creates an upper respiration rate limit for the modified model based on cow mass. For example, maximum respiration rates decrease from 138bpm for a 400kg cow to 120 bpm for a 650kg cow.

Body Heat Storage

McGovern and Bruce (2000) included body heat storage to account for the energy used to raise a cows body temperature above a specified reference value, TN. The heat storage flux, Gb, depended on both the magnitude of the body temperature rise and the time over which the temperature rise occurred. The body temperature and heat storage equation in McGovern and Bruce (2000) was rearranged to find body heat storage given the body temperature rise above the reference value (eq. 23):

(23)

where

Gb = heat storage flux (W m-2) over time, t (h)

cb =Heat capacity of the cows body (3,400 J kg-1 C-1)

t1 = Heat storage time (h).

The heat storage time, t1, impacts heat storage flux and other model results. Increasing heat storage time reduces the heat storage flux value and causes the respiration volume and rate to increase, which causes the body temperature to increase too. At sufficiently long heat storage times, the respiration rate exceeds the maximum allowed (rrmax =120 bpm for a 650 kg cow or 138 bpm for a 400 kg cow), and the modified model cannot balance the system of equations. Decreasing the heat storage time increases the heat storage flux and causes the respiration volume and rate to decrease, which causes a decrease in body temperature. Gebremedhin et al. (2010) reported large body temperature rises within 3h when solar loading was 850 W m-2. The heat storage time used was 3 h.

Sensible Heat Flux from Body Core to Skin Surface

The sensible heat flux, Qb, found using equation 1, can also describe heat flow from the cows core through her skin using a conductive heat transfer relation using the temperature difference between the skin surface temperature, Ts, and the body temperature, Tb, and an average tissue insulation value, Ib, for the entire modeled cow (McGovern and Bruce, 2000). The conductive heat relation can be used to calculate the skin temperature using equation 24:

(24)

Figure 3. Tissue insulation as function of respiration rate for 600 kg cow.

where

Ts = Skin temperature (C)

Ib = Tissue insulation (m2 K W-1).

The tissue insulation, Ib, is a proxy for the cows peripheral blood flow, which is under the cows physiological control. Tissue insulation is reduced when peripheral blood flow is increased through vasodilation to increase sensible heat transfer between the body core and the skin surface. Asair dry-bulb temperature rises, skin temperature rises resulting in a reduced temperature gradient between body core and the skin surface (i.e., Tb-Ts). When the temperature gradient is reduced, the effective tissue insulation must also be reduced, for metabolic heat to be transferred to the skin surface and dissipated to the surrounding environment.

McArthur (1987) related decreased tissue insulation to increased skin temperature. McGovern and Bruce (2000) found skin temperature using tissue insulation as an independent variable. Iterative techniques would have been needed to solve these equations, while possible, it was deemed to be too complex an approach.

Berman (2005) related linear tissue insulation reductions to air temperature increases when above a lower critical temperature. This linear tissue insulation relation was not selected for use.

For this model, tissue insulation was found using a non-linear function dependent on respiration rate (fig. 3). The non-linear inverse relation was used between the specified maximum tissue resistance value at the warm conditions respiration rate, rrwarm, (i.e., approximately 12bpm for a 600 kg cow [McGovern and Bruce, 2000]) and the minimum tissue resistance at 50% of the rapid shallow panting phase respiration rate (i.e., approximately 68 bpm for a 600 kg cow [McGovern and Bruce, 2000]). The 50% respiration rate in the rapid shallow panting phase was the average of the respiration rate at warm conditions, rrwarm, and the maximum respiration rate, rrmax.

The maximum tissue insulation was calculated following Berman (2004) using equation 25:

(25)

where

Ib,max = Maximum tissue insulation, (m2 K W-1)

The minimum tissue insulation, Ib,min, was set at 0.016 m2 K W-1 (McGovern and Bruce, 2000).

The nonlinear tissue insulation function used in this model, based on respiration rate and the maximum and minimum tissue insulation values, is shown in figure 3 and found using equation 26:

(26)

where

Ib = tissue insulation at respiration rate rr (m2 K W-1)

a1 and b = constants found using equations 27 and 28, respectively:

(27)

where

a1 = Parameter in equation 26 to find tissue insulation as a function of respiration rate (m2 K W-1 bpm-1)

Ib,min = 0.016 m2 K W-1

rr50%resp = Respiration rate at 50% of first phase respiratory cooling (bpm).

(28)

(29)

where

b = Parameter in equation 26 to find tissue insulation as a function of respiration rate (m2 K W-1)

The minimum tissue insulation value was commonly used during heat stress conditions.

Latent Heat Loss from Skin

The latent heat loss flux due to sweat evaporation from the skin depends on sweat rate, skin temperature, and air temperature. A modified sweat rate relation was developed for this model.

Sweat rate is important because skin evaporative cooling can account for up to 85% of the total heat loss when air temperatures exceed 30C (Maia et al., 2005). Thompson et al. (2011) compared four empirical relationships relating sweat rate to skin temperature for Bos taurus (Gatenby, 1986; McArthur, 1987; Maia et al., 2005). Data scatter around the best-fit empirical relations was significant and showed large variation within breeds and between animals (Gatenby, 1986; McArthur, 1987; Maia et al., 2005). McArthur (1987) used linear relations, while Gatenby (1986), Maia et al. (2005), and Thompson et al. (2011) used exponential relations. For this model, the cows sweat rate was described with a new exponential function based on skin temperature, like Gatenby (1986), McArthur (1987), Maia et al. (2005), and Thompson et al. (2011).

Minimum Sweat Rate and Corresponding Lower Critical Skin Temperature

The first point needed to find the exponential relation between sweat rate and skin temperature was the minimum sweat rate due to diffusion, 14.4 g m-2 h-1 (McArthur, 1987), and a corresponding skin temperature. The skin temperature was called the lower critical skin temperature, Ts,LCT, and was found using equation 24. The Ts,LCT is the skin temperature when the metabolic heat production, M, was dissipated to the surroundings with maximum tissue insulation, Ib,max, minimum respiratory cooling, Qr,min, body storage set to zero (Gb = 0) and the body temperature was normal (Tb = TN). The metabolic heat production value used for the Ts,LCT calculation was the maintenance and milk production, Mmaint+lact, heat production found using equation 10. Maximum tissue insulation was calculated using equation 25. The minimum respiratory cooling value was found using equation 14 for a lower critical dry-bulb temperature of 0C, a lower critical dew-point temperature of -8C, an exhaled air temperature found using equation 15, and the warm respiration rate, rrwarm. The warm respiration rate was found using equation17, and the tidal volume at warm conditions was found using equation 22.

Sweat Rate at 50% Total Respiratory Cooling Capacity and Corresponding Skin Temperature

The second sweat rate and corresponding temperature needed to find coefficients for the exponential sweat rate relation were at 50% total respiratory cooling capacity and were a function of a maximum physiological sweat rate, , and a maximum potential environmental sweat rate, (eq. 31). The maximum physiological sweat rate was assumed to be = 660 g m-2 h-1 based on Gebremedhin et al. (2010). The sweat rate at 50% of total respiratory cooling capacity was found using equation 30:

(30)

where

= Sweat rate at 50% total respiratory cooling capacity (g m-2 h-1)

= Maximum physiological sweat rate, (660gm-2 h-1)

= Maximum potential environmental sweat rate, (g m-2 h-1).

The maximum potential environmental sweat rate relation, equation 31, was from Turnpenny et al. (2000b):

(31)

where

Cva = Average air volumetric specific heat using ambient air and skin temperature, (J m-3 K-1)

Pws,Ts = Saturation vapor press at skin surface temperature (Pa)

Pw = Partial vapor pressure at ambient temperature (Pa)

?ave = Average psychrometric constant value using ambient and skin temperatures (Pa K-1)

?ave = Average latent heat of vaporization using ambient air and skin temperatures (J kg-1)

rv = Resistance to vapor transfer in the coat (s m-1)

3.6 = unit conversion (3.6 s kg h-1 g-1) (Turnpenny et al., 2000b).

The psychrometric constant, ?, was calculated using equation 32:

(32)

where

P = Atmospheric pressure, (Pa)

? = Latent heat of vaporization (J kg-1)

0.622 = Ratio of the molecular weight of water vapor to dry air.

The resistance to vapor transfer, rv, relation needed in equation 31, was taken from Turnpenny et al. (2000b), who had obtained it from Cena and Monteith (1975) (eq. 33):

(33)

where

rv = Resistance to vapor transfer in the coat (s m-1)

lc = Coat thickness (mm)

?lw = Wind penetration depth (m)

D = Averaged diffusivity of water vapor in air at air and skin conditions (m2 s-1)

= Skin virtual temperature (K)

= Ambient air virtual temperature (K).

Virtual temperatures, . and for equation 33, were calculated using equation 16 for the corresponding skin temperature and saturated vapor pressure and the air temperature and partial vapor pressure, respectively. Water vapor mass diffusivity values were adjusted for air and skin temperatures and pressure using equation 34 based on a water diffusivity value of 0.26 x 10-4 m2 s-1 at 101,325Pa pressure and 298K temperature (Welty et al., 2008):

(34)

where TK was the Absolute air or skin temperature (K).

Coat depth reduction due to wind was found using equation 35 from McGovern and Bruce (2000) as specified by McArthur (1987):

(35)

where

?lw = Coat depth decrease due to wind (m)

u = Wind speed at cow height (m s-1)

ld = lc / (1000 mm/m) = Coat depth (m)

19.0 10-6 m = Given constant

11,000 s m-2 = Given constant.

The corresponding skin temperature at 50% total respiratory cooling capacity, Ts,50%Qr,n, was calculated using equation 24, which requires a body temperature. For this model, body temperature was a function of respiration rate and is discussed later. The respiration rate at 50% total respiratory cooling capacity was found using equation 36, the tidal volume at warm conditions, Vt,warm, (eq. 22), and the maximal tidal volume, Vt,severe, (eq. 20):

(36)

where Vt,50%resp = Tidal volume at 50% total respiratory cooling capacity (L breath-1).

The corresponding respiration rate at 50% total respiratory cooling capacity, rr50%totresp, was found using equation37 which uses the corresponding respiration at Vt,50%resp, and the calculated rrwarm and rrmax:

(37)

where

rr50%totresp = Respiration rate at 50% total respiratory cooling capacity (bpm)

With rr50%totresp the corresponding body temperature at 50% total respiratory cooling capacity needed in equation 24 was calculated using equation 58, which is discussed later. The minimum tissue insulation, Ib,min, value used was 0.016 m2 KW-1 (McGovern and Bruce, 2000) and body heat storage was set to zero.

Sweat Rate Relation as a Function of Skin Temperature

The exponential sweat rate relation as a function of skin temperature used was similar to those used by Maia et al. (2005) and Gatenby (1986). The relation used uses the minimum sweat rate and its corresponding skin temperature and the sweat rate at 50% total respiratory cooling capacity and its corresponding skin temperature, Ts,LCT, found in the previous sections. Coefficient a2 in equation 38 satisfies the boundary conditions given by the two points. The sweat rate relation is:

(38)

where

= Sweat rate as a function of skin temperature, Ts. (gm-2 h-1)

= Minimum sweat rate (14.4 g m-2 h-1) (McArthur, 1987)

Ts,LCT = Lower critical skin temperature (C)

The parameter a2 was computed using equation 39:

(39)

where

= Skin temperature at 50% total respiratory cooling (C)

= Sweat rate at 50% total respiratory cooling (gm-2 h-1)

= Minimum sweat rate = 14.4 g m-2 h-1.

The sweat rate relation given by equation 38 includes the impact of the surrounding environment on sweat rate, including higher air humidity, and reduced wind speed that leads to reduced sweat rates, as observed by Gebremedhin et al. (2010).

Equation 40 was used in this model to calculate the latent heat loss flux from the skin given the sweat rate, , and exposed cow surface area, Ae:

(40)

where

Es = latent heat loss flux from skin (W m-2)

?s = Latent heat of vaporization at the skin temperature (Jkg-1)

Ae = Percent surface area exposed (%).

Convective Heat Flux from Coat to Surface to Surrounding Air

Sensible Heat Flux from Skin Surface to Coat Surface

For this model the sensible heat flow through the coat, Qs,c, was written in terms of the skin to coat temperature differential (Ts - Tc):

(41)

where

hc = Coat conductance (W m-2 C)

Tc = Coat surface temperature (C).

Coat conductance was found using the recommended hair coat insulation, CI, value per mm of coat depth for Holsteins of 0.30 C m2 d Mcal-1 mm-1 (Berman, 2004) and adjusting for coat depth reduction due to wind using equation 35.

Coat conductance was found using equation 42:

(42)

where

CI = Coat insulation value for Holsteins (0.30 C m2 d Mcal-1 mm-1) (Berman, 2004)

Unit conversion

Convective Heat Flux from Coat Surface to Surrounding Air

The boundary layer convective heat flow from the coat surface to the surrounding air was written in terms of the coat surface temperature to air temperature differential (Tc - Ta). Equation 43 was consistent with Turnpenny et al. (2000a) and Thompson et al. (2014):

(43)

where

C = Convective sensible heat flux from coat surface to surrounding air (W m-2)

k = Thermal conductivity of air using average of air and coat surface temperatures (W m-1 K-1)

Nu = Nusselt number (dimensionless).

The Nusselt number relation used depended on whether the air velocity past the cow led to forced convection, natural convection, or the transition between them based on the Grashof number, Gr, and the Reynolds number, Re. For natural convection, when Gr > 16 Re2, the Nusselt number relation used was:

(44)

For forced convection conditions (when Gr < 0.1 Re2) the Nusselt number relation used was:

(45)

For intermediate conditions (0.1 Re2 < Gr < 16 Re2) the larger Nu value was used.

The Reynolds number at the coat and air interface was calculated using:

(46)

where

Recoat = Reynolds number at the coat and air interface (dimensionless)

vcoat air = Kinematic viscosity at the coat and air interface (m2 s-1).

The airs kinematic viscosity at the coat and air interface is the dynamic viscosity divided by the air density. The dynamic viscosity was adjusted for temperature using Sutherlands formula (Sutherland, 1893) and divided by air density using equation 47. The reference temperature used was 30C and the corresponding reference dynamic viscosity was 18.6x 10-6 N s m-2:

(47)

where

Tcoatave = Average of coat surface and air temperature (C)

?a = Average air density of coat surface and air temperature (kg m-3)

The Grashof number needed to find the Nu (using either equation 44 or 45) was found using equation 48 (Thompson et al., 2014) and an average coat partial pressure between the skin and the coat (eq. 49):

(48)

where

Grcoat = Grashof number at coat air interface.

(49)

where

Pwcoat = Average coat partial pressure using partial air pressure and skin saturation pressure (Pa).

Longwave Radiation Heat Flux Exchange with Surrounding Surfaces

The longwave radiant heat exchange for cows located outdoors without shade is between the cows coat and the sky and ground. For cows in barns or chambers, the longwave radiation is between the cows coat and surfaces such as the underside of the barn roof, walls, or adjacent cows.

Longwave Radiation in Barn or Chamber

The two temperatures used to describe the longwave radiation for cows in a barn or chamber were the coat temperature, Tc, and a mean radiant temperature, Tmrt, which represented the temperature of surfaces surrounding the cow (i.e., floor, walls, ceiling, or the underside of the roof, or adjacent animals). The cows coat temperature was calculated by the modified model. The mean radiant temperature was an assumed input and constant. The longwave radiation flux was adjusted based on the amount of exposed cow surface, Ae, in equation 50:

(50)

where

LWn = longwave radiation heat flux between cows coat and surrounding surfaces when inside a barn or chamber (Wm-2)

f = Longwave exposure factor, f =1 (Dimensionless)

ecoat = Emissivity of cows coat, 0.98 (Dimensionless)

s = Stefan-Boltzmann constant 5.67 x 10-8 W m-2 K-4

Tmrt = Surrounding surfaces mean radiant temperature (C).

Longwave Radiation for Cows Outside

The longwave radiation exchange between the cows coat and the surrounding ground and sky when the cow is outside was described using equation 51 (Thompson et al., 2014):

(51)

where

qlwrad = Longwave radiation flux between cows coat and surroundings when outside (W m-2)

eground = Ground emissivity (0.9)

esky = Sky emissivity (Dimensionless)

Tc = Coat surface temperature (K)

Tsurrounding = Temperature of the surroundings, which was an average of the ground temperature and the sky temperature(K).

Thompson et al. (2014) used equations from Brutsaert (1975) and Crawford and Duchon (1999) to calculate the sky emissivity. The equations were combined in equation 52:

(52)

where

esky = Cloudy sky emissivity, (Dimensionless)

CCf = Cloud cover factor (Dimensionless)

Pws = Saturation vapor pressure at air temperature (Pa)

Ta,K = Air dry-bulb temperature (K).

The CCf value used depended on whether the horizontal solar load, Snh, was an input or calculated value. When Snh was calculated using equation 59 based on location latitude, day of the year, and hour of the day, the CCf value was calculated using equation 53 as described by McGovern and Bruce (2000) using an input okta value, CSf (e.g., CSf = 0 for a clear sky to 8 for a completely cloudy sky):

(53)

When Snh was used as an input, the CCf was calculated using equation 54 from Crawford and Duchon (1999), with the theoretical solar load calculated using equation 59 and the input solar load. If the input solar load was greater than the theoretical solar load, CCf was set to 0, which corresponds to clear skies. The theoretical solar load was calculated using equation 59 and is discussed later:

(54)

The surrounding temperature for a cow out of doors, Tsurround, was the average of the sky temperature and the ground temperature (Thompson et al., 2014). The sky temperature was calculated using equation 55 from Gwarda et al. (2017), which uses a total clear sky emissivity relation from Berdahl and Martin (1984) based on the air dew-point temperature and the solar time:

(55)

where

Tsky,K = Surrounding sky temperature (K)

Ta,K = Air dry-bulb temperature (K)

Tdp = Air dew-point temperature (C)

t2 = Solar hour (h).

The ground temperature was calculated using equation 56 from Gwadera et al. (2017) which included terms for solar irradiation and water evaporation from the ground or plants:

(56)

where

Tground = Ground surface temperature (K)

Snh = Solar irradiance (shortwave radiation) on a horizontal surface (W m-2)

hg = Ground convective heat transfer coefficient (33 W m-2 K-1)

Evap = Evaporation flux (W m-2)

eground = Ground emissivity (0.9 dimensionless)

LWg = Ground longwave radiation flux (63 W m-2)

The constants for hg, e, and LWg were taken from Gwarda et al. (2017) for the simplest model of ground temperature. The evaporation flux was an input. The coat surface temperature was an unknown found by the modified model. The shortwave irradiation, Snh, is discussed in the next section.

Shortwave Radiation Heat Flux from the Sun

Shortwave radiation in this modified model describes radiant energy from the sun or solar lamps used in chamber studies. The modified model shortwave solar load, Rn (eq.60), is the product of the solar load applied times the amount of the cows surface area absorbing shortwave radiation per cow total exposed surface area, Ah/A, one minus the coat reflectivity, ecoat, and one minus the fraction of cow shaded, ss. For cows outside, the Ah/A value was calculated using the area of the shadow cast by the cow (eq. 57) as a function of the solar altitude, hour angle, solar declination angle, solar azimuth, and the azimuth angle of the animal (McGovern and Bruce, 2000):

(57)

where

Ah/A = Horizontal area absorbing shortwave radiation per unit area of cow (m2 m-2)

= Solar altitude which is a function of date, declination angle, latitude, and solar hour (degrees) calculated using equations in the ASHRAE (2021)

?T = Angular difference between the solar azimuth, ?s, and the azimuth of the cow, ?cow (degrees)

x =(2lt)/dt where lt is the cow trunk length and dt is the trunk diameter (dimensionless).

When the shortwave radiant source was directly overhead, as in chamber studies, Ah/A simplified to:

(58)

The horizontal shortwave load, Snh, needed to calculate the shortwave radiation heat flux, Rn, can be either measured or estimated by calculating the total irradiance as a function of date, latitude, solar hour, cloud cover, albedo, and an atmospheric dispersion factor using equation 59 from McGovern and Bruce (2000) and extraterrestrial normal irradiance values for the 21st of each month and solar angles needed from ASHRAE (2021):

(59)

where

Snh = Direct and diffuse horizontal shortwave radiation heat flux (W m-2)

S0 = Extraterrestrial normal irradiance (W m-2)

f = Albedo (Dimensionless)

CCf = Cloud cover factor (Dimensionless)

= Solar altitude (degrees)

aa = Atmospheric dispersion factor (Dimensionless).

This model was also set up to use input horizontal solar radiation or sun lamp values (W m-2). Whether using a calculated or given measured Snh, the Rn was found using equation 60:

(60)

where

Rn = Shortwave radiation heat flux (W m-2)

rcoat = Animal coat reflectivity (0.3 dimensionless)

ss = Sun shading factor (Dimensionless decimal value).

Shading

The amount of shading was adjustable using a decimal value, ss, to describe the fraction of solar shortwave radiation that reached the cow as one minus the decimal fraction blocked. In full sun, the fraction blocked was zero in equation 60. When the cow was assumed to be completely in a barn, the fraction of shading was set to 1, and no shortwave radiation reached the cow. When shading was set to 1, longwave radiant energy exchange was calculated using equation 50 and the surrounding surfaces mean radiant temperature, Tmrt.

Body Temperature and Respiration Rate

The body temperature relation used in the modified model was a linear function of respiration rate. After nearly complete vasodilation, the body temperature begins to rise (McArthur, 1987). Thompson et al. (2011) reported that elevated respiration rates and body temperatures were well correlated for Bos taurus, crossbreds, and Bos indicus in three studies (Allen, 1962; Thomas and Pearson, 1986; and Brown-Brandl et al., 2003). Some studies reported correlations between body temperatures and respiration rates of r = 0.046 (Kabuga, 1992) and r = 0.55 (Martello et al., 2010).

Recently, Li et al. (2020) collected data from 45 high-producing Holstein cows on 26 days between 1000-1100 h to develop a linear relation between mean rectal temperature, MRT, and mean respiration rate, MRR. The relation Li et al. (2020) reported was MRT = 0.021 MRR + 37.6. They also determined that heat stress was triggered when MRT = 38.6C, which corresponded to a respiration rate of 48 bpm (Li et al., 2020).

For this model, the regression relation relating respiration rate to body temperature for Bos taurus from Thompson et al. (2011) was rearranged to find body temperature given the respiration rate, rr, in bpm (eq. 61):

(61)

The reference body temperature was set at 38.0C, which corresponded to a respiration rate of 37 bpm. When respiration rates were calculated to be below 37 bpm, the body temperature was set at the reference value of 38.0C because the modified model was designed for warm and hot weather conditions, not cool conditions, and 38.0C was at the lower range of rectal temperatures for dairy cows (Andersson and Jonasson, 1993; Collier and Gebremedhin, 2015). At higher respiration rates above 37 bpm, the body temperature was found using equation 61.

Model Implementation

A spreadsheet was developed to solve the system of steady-state equations describing heat exchange between a cow and the surrounding environment by convection, shortwave radiation, longwave radiation, evaporation from the skin, and sensible and latent respiratory heat loss. The model inputs and coefficients listed in table 1 were adjustable. For this article, the steady state analysis was done over either 1 or 3 h time steps in which the weather, solar, cow, and other characteristics were assumed to be constant.

In the spreadsheet, the inputs and equation coefficients were used to calculate intermediate constants not dependent on respiration rate or body, skin, or coat temperatures (fig 2). The calculated constants include cow surface area (eq. 5), trunk diameter and length (eqs. 6 and 7), heat production for maintenance and milk production (eqs. 8, 9, and 10), respiration rate and tidal volume indices (eqs. 17-22), maximum tissue insulation (eq. 25), and exhaled air temperature (eq.15). A spreadsheet add-in program, Solver, was used to adjust the inspired volumetric respiration rate, Vresp, to find the modified model cow body, skin, and coat temperatures, respiration rate, and numerous intermediate values that balanced the overall heat balance given by equation 4.

Model Relations Assessed

The sweat rate and tissue insulation relations used were compared to published relations. The sweat rate relation developed for this model was compared to empirical relations relating sweat rate to skin temperature (Gatenby, 1986; McArthur, 1987; Maia et al., 2005; Thompson et al., 2011). Then the model tissue resistance, body temperature, and skin temperature results were compared to equations by McArthur (1987) relating tissue resistance to skin temperature and body temperature to skin temperature.

For each analysis, the modified model was run for a 560kg cow on pasture at the WCROC with a 32.9kgd-1 milk yield and 3.7% milk fat. Weather data was recorded automatically every 3 to 5 min at a University of Minnesota Morris weather station (Vantage Pro, Davis Instruments Corp., Hayward, CA) located at 4535'20? N and 95 54'7? W, less than 1 km from the WCROC pasture. The weather data was averaged over 3-h time blocks (i.e., 04:45-07:30, 07:45-10:30, 10:45-13:30, 13:45-16:30, 16:45-19:30) over six days (i.e., 10 to 15 July 2015). The recorded weather data used included the air dry-bulb temperature, dew-point temperature, air velocity at 10 m, horizontal solar insolation, and atmospheric pressure. Changes in milk production were not recorded or modeled. Finally, an example case is presented with inputs and results to illustrate the modified models use.

A subsequent article assesses thermal balance model results against published body temperatures, respiration rates, skin temperatures, and unpublished body temperatures of grazing cows for heat-stressed lactating cows.

Results and Discussion

Table 2 summarizes the weather conditions based on 360readings recorded every 15 min from 04:45 to 19:30 each day from 10 to 15 July 2015. Averaged values from 3-h periods were used as model inputs for the sweat rate and tissue insulation analyses.

Table 2. Summary of weather conditions at the WCROC between 04:45 and 19:30 from July 10-15, 2015.
Weather variableAverageStandard DeviationRange
Air dry-bulb temperature (C)26.44.214.6-33.4
Air dew-point temperature (C)20.72.612.2-27.8
Temperature Humidity Index754.960 - 85
Wind speed at 10 m height (m s-1)3.51.60.0-7.1
Wind speed at mid cow height hm = 1 m (m s-1)[a]2.81.20.0-5.8
Solar insolation (W m-2)4363120-1,074
Atmospheric pressure (Pa)100,76638599,976-101,331

    [a]Adjusted for mid cow height using Power law, ucow = uw (hm/10)p, for wind speed, uw, measured at 10 m height, input mid cow height, hm=1m, and exponent, p = 0.1, for stability class C and smooth surface (Cooper and Alley, 2011).

Based on calculated Temperature Humidity Index (THI) values, the weather conditions observed were not considered too heat stressful, but solar loads added to the heat stress experienced by cows on pasture in the sun.

Sweat Rate Analysis

Figure 4 is a plot of sweat rate versus skin temperature for four models (Gatenby, 1986; McArthur, 1987; Maia et al., 2005; Thompson et al., 2011) and results from the modified model for average weather conditions for the 30 3-h periods in July 2015. McArthur (1987) used a linear relation, while Gatenby (1986), Maia et al. (2005), and Thompson et al. (2011) used exponential relations. Thompson et al. (2011) compared the four empirical relationships against three data sets (Allen, 1962; Thomas and Pearson, 1986; Brown-Brandl et al., 2003) and concluded that their exponential model was an improvement over the other sweat rate to skin temperature relations for Bos taurus.

The sweating rates obtained using the modified model were consistently above those predicted by the relation developed by Thompson et al. (2011) but followed a very similar trend. The modified model sweat rates were well above those predicted by McArthur (1987) and Gatenby (1986) at skin temperatures below 34C and below those at skin temperatures above 36.5C. The modified model sweat rate calculations include the impact of several other factors in addition to skin temperature and show some scatter. Additional data and analyses will be needed to assess the impact of air humidity and air speed past the cow on sweat rates.

Figure 4. Empirical models relating sweat rate to skin temperature (Gatenby, 1986; McArthur, 1987; Maia et al., 2005; Thompson et al., 2011).

Tissue Resistance Analysis

Model tissue resistance and skin temperature results (open circles) for 30 3-hour time blocks over six days (i.e., 10 to 15 July) in 2015 were compared to McArthurs equation (solid line) in figure 5. The modified model tissue resistance values found using the non-linear relation describing tissue insulation as a function of respiration rate were consistently lower than tissue resistance values found using the empirical equation McArthur (1987) developed. The skin temperature where tissue resistance reaches a minimum was found to be around 36.5C, which was roughly 0.5C higher than the relation found by McArthur (1987).

McArthur (1987) also computed body temperature based on skin temperature above a minimum value using data for Holsteins. In this spreadsheet model, body temperature was related to respiration rate, following Thompson et al. (2011). Figure 6 shows the modified model body temperature versus skin temperature results (open circles) and the relation from McArthur (1987). The body temperatures from the modified model increased smoothly as the skin temperature increased. The modified model body temperatures were consistently above the body temperatures predicted by McArthur (1987) for skin temperatures below 36.5C but were very close above 36.5C skin temperatures. This threshold skin temperature was below the skin temperature at which tissue insulation stabilizes at a minimum at approximately 36.5C (fig.5). McArthur (1987) used a minimum body temperature of 38.2C, while the modified model used 38C. The overlap in physiological responses (i.e., tissue resistance decline [fig.5] and body temperature rise [fig. 6]) with rising heat load agrees with the discussion by McArthur (1987).

Example Model Inputs and Results

Spreadsheet model inputs in tables 1 and 3 were used to illustrate model use for a hypothetical cow standing in a pasture in the sunshine during the 3-h period between 13:45 and 16:30 on 12 July 2015. This 3-h period had the highest dry-bulb temperature and THI between 10-15 July 2015 at the WCROC. Stability class C was selected because the incoming solar radiation was strong and the wind speed at 10 m was greater than 6 m s-1 (Cooper and Alley, 2011).

Table 4 lists spreadsheet model outputs and select calculated intermediate values. Based on the THI value the cow was experiencing moderate-severe heat stress and would be expected to have a respiration rate greater than 85 bpm and rectal temperature greater than 40C (Renaudeau et al., 2012). The modified model results were a body temperature of 39.8C and a respiration rate of 86 bpm.

The heat production values in table 4 demonstrate the important role that milk yield plays in total metabolic heat production flux. For this case, milk yield accounted for 141 W m-2 (58%) of the total metabolic heat production. A cow with a lower milk yield would have a lower amount of metabolic heat to dissipate in warm and hot weather to maintain homeothermic conditions. Milk yield is not considered when assessing heat stress using the temperature-humidity-index (THI) or the equivalent temperature index for dairy cattle (ETIC) (Wang et al., 2018a, 2019).

The heat flux results indicate the importance of the different methods of heat exchange. For this case, the total flux the cow must dissipate was 386 W/m2, the sum of the shortwave radiant flux absorbed (144 W/m2) and the metabolic heat flux generated (242 W/m2). The shortwave radiation added 37% to the total flux to be dissipated. The evaporative heat loss flux dissipated over 43% of the total, while the convective heat loss flux dissipated over 21% of the total heat dissipated. The remaining heat flux was either dissipated by convection, respiration, and longwave radiation or stored in the body.

Figure 5. Modified model tissue resistance (dots) versus skin temperature and equation from McArthur (1987) (line).
Figure 6. Model body temperatures (dots) versus skin temperature and body temperatures (line) computed using empirical equation from McArthur (1987).

Conclusion

A modified steady-state process-based heat transfer model was developed, building on work by McArthur (1987), McGovern and Bruce (2000), Turnpenny et al. (2000a,b), Berman (2005), Thompson et al. (2014), and Gwadera et al. (2017) to describe the heat exchange between a lactating cow and the surrounding environment needed to maintain homeothermic conditions under warm or hot conditions. New and modified relations were used to describe tissue insulation as a function of respiration rate, sweat rate as a function of skin temperature, body temperature as a function of respiration rate, heat storage as a function of body temperature, longwave radiation exchange, and coat and convective boundary heat fluxes as functions of skin temperature, coat temperature, and air temperature. The modified model sweat rate, tissue resistance, and body temperature results versus skin temperatures were similar to those published by McArthur (1987) and Thompson et al. (2011).

Table 3. Model inputs for the 15:00 solar hour on July 12, 2015, not listed in table 1.
Variable descriptionValueVariable descriptionValue
Cow mass560 kgAtmospheric pressure100,230 Pa
Milk production32.9 kg d-1Air dry-bulb temperature32.1C
Milk fat3.7%Air dew-point temperature26.9C
Mid cow height1.0 mSolar load (measured)535 W m-2
Percent cow surface area exposed100%Date and Solar timeJuly 12 @ 15:00
Cow azimuth relative to south0 degSolar shading factor0
Wind speed at 10 m height6.5 m s-1Latitude for the WCROC45.35N
Wind speed at mid cow height of 1.0 m[a]5.1 m s-1Evaporation flux for finding ground temperature52 W m-2
[a] Adjusted for mid cow height using Power law, ucow = uw (hm/10)p , for wind speed, uw , measured at 10 m height, input mid cow height, hm = 1 m, and exponent, p = 0.1, for an assumed stability class C and smooth rural surface (Cooper and Alley, 2011).
Table 4. Model outputs.
Variable descriptionValueVariable descriptionValue
Cow characteristicsAmbient air conditions
Cow surface area5.2 m2Relative humidity74%
Cow trunk diameter0.71 mTemperature Humidity Index83
Cow trunk length1.98 mAir density1.13 kg m-3
Coat depth decrease1.56 mmHumidity ratio
Maximum tissue insulation0.112 m2 K W-1Specific heat of air1006 J kg-1 K-1
Tissue insulation20.5 s m-1Volumetric specific heat of air1136 J m-3 K-1
Coat conductance112 W m-2 C-1
Sweat rate252 g m-2 h-1Temperatures
Body temperature39.8C
Respiration valuesSkin temperature37.6C
Warm respiration rate limit12 bpmCoat temperature37.9C
Transition respiration rate limit126 bpmExhaled air temperature38.1C
Severe respiration rate limit88 bpmLower critical skin temperature14C
Tidal volume at warm conditions3.79 L breath-1Skin temperature at 50% respiration rate36C
Tidal volume at maximum respiration 2.12 L breath-1Ground temperature315 K
Maximal tidal volume limit based on cow mass3.98 L breath-1Sky temperature298 K
Inspiration volumetric rate193 L min-1Ambient air virtual temperature309 K
Respiration rate88 bpmExhaled air virtual temperature319 K
Skin virtual temperature318 K
Heat production
Maintenance9.2 Mcal d-1Heat flux values
Lactation15.1 Mcal d-1Heat flux from body core to skin surface137 W m-2
Maintenance plus lactation production flux227 W m-2Respiratory heat flux43 W m-2
Flux increase due to body temperature rise12 W m-2Convective heat flux from coat to air84 W m-2
Flux increase due to increased respiration3 W m-2Evaporative (sweat) heat flux169 W m-2
Total metabolic heat production flux242 W m-2Body storage61 W m-2
Shortwave radiation absorbed144 W m-2
Solar characteristicsLongwave radiation outbound29 W m-2
Hour angle45 degOutbound skin to coat sensible heat flux-31 W m-2
Solar declination angle22 degOutbound coat surface to surroundings -31 W m-2
Solar altitude46.6 deg
Solar azimuth72.8 degSweat rate values
Extraterrestrial irradiance for July1324 W m-2Parameter a2 in sweat rate function (eq. 38)8.37
Cows azimuth angle0 degSweat rate at 50% respiratory cooling218 g m-2 h-1
Azimuth angle of cow to sun72.8 degMaximum environmental sweat rate1280 g m-2 h-1
Horizontal area absorbing shortwave radiation0.385 m2 m-2Sweat rate at skin temperature252 g m-2 h-1
Solar irradiance on a horizontal surface206 W m-2Tidal volume at 50% total respiration cooling198 L min-1
Net shortwave radiation absorbed144 W m-2Resistance to vapor transfer in coat51 s m-1
Sun shading factor0.85
Cloud cover factor0.15Additional psychrometric values
Parameter a1 used in equation 261.405Exhaled air saturation humidity ratio
Parameter b used in equation 26-0.00441Density exhaled air1.11 kg m-3
Cloudy sky emissivity0.96Latent heat of vaporization ambient air2.42106 J kg-1
Latent heat of vaporization exhaled air2.41106 J kg-1
Convection related valuesLatent heat of vaporization skin temperature2.41106 J kg-1
Kinematic viscosity at coat and air interface1.6810-5 m2 s-1Ambient air partial vapor pressure3,546 Pa
Nusselt Number384Saturation vapor pressure at skin temperature6,490 Pa
Reynolds Number2.19105Psychrometric constant at air temperature66.9
Grashof Number2.64108Water vapor diffusivity in air at air temperature2.7210-5 m2 s-1
Water vapor diffusivity in air at skin temperature2.8010-5 m2 s-1
Air thermal conductivity at air temperature0.0265 W m-2 K-1
Air thermal conductivity at coat surface temperature0.0269 W m-2 K-1

The modified model, which can be solved with a spreadsheet, provides insight into factors that affect lactating cow heat exchange in addition to dry-bulb and dew-point temperature. A companion paper compares model results with published measured body temperatures, respiration rates, and skin temperatures. The companion paper also compares model results with unpublished reticular temperatures of cows on pasture. It is expected that the modified model can help assess the relative importance of factors affecting heat exchange (i.e., body mass, daily milk yield, coat thickness, coat reflectivity, tissue insulation, solar load, exposed surface area, and air velocity) and their impact on heat exchange, cow respiration rate, and body temperature.

Acknowledgments

The authors thank Kirsten Sharpe and Brad Heins for sharing weather and cow data from the University of Minnesota West Central Research and Outreach Center.

References

Abuajamieh, M. K., Stoakes, S. K., Kohls, D. S., Groen, B. J., Lahr, D. A., Kuehn, C. R., & Baumgard, L. H. (2013). Heat stress: Practical lessons learned in 2011. Proc. 4-State Dairy Nut. and Mgmt. Conf., (pp. 61-73).

Allen, T. E. (1962). Responses of Zebu, Jersey, and Zebu X Jersey crossbred heifers to rising temperature, with particular reference to sweating. Austrailian J. Agric. Res., 13(1), 165-179. https://doi.org/10.1071/AR9620165

Andersson, B. E., & Jonasson, H. (1993). Temperature regulation and environmental physiology. In M. J. Swenson, & W. O. Reece (Eds.), Dukes physiology of domestic animals. Ithaca, NY: Cornell University.

ASHRAE. (2021). 2021 ASHRAE Handbook Fundamentals. Atlanta, GA: American Society of Heating, Refrigeration and Air-Conditioning Engineers.

Berdahl, P., & Martin, M. (1984). Emissivity of clear skies. Sol. Energy, 32(5), 663-664. https://doi.org/10.1016/0038-092X(84)90144-0

Berman, A. (2003). Effects of body surface area estimates on predicted energy requirement and heat stress. J. Dairy Sci., 86(11), 3605-3610. https://doi.org/10.3168/jds.S0022-0302(03)73966-6

Berman, A. (2004). Tissue and external insulation estimates and their effects on prediction of energy requirements and of heat stress. J. Dairy Sci., 87(5), 1400-1412. https://doi.org/10.3168/jds.S0022-0302(04)73289-0

Berman, A. (2005). Estimates of heat stress relief needs for Holstein dairy cows. J. Anim. Sci., 83(6), 1377-1384. https://doi.org/10.2527/2005.8361377x

Blaxter, K. L. (1967). The energy metabolism of ruminants. Hutchinson & Co.

Brody, S. (1945). Bioenergetics and growth; with special reference to the efficiency complex in domestic animals. New York: Reinhold Publ.

Brown-Brandl, T. M., Nienaber, J. A., Eigenberg, R. A., Hahn, G. L., & Freetly, H. (2003). Thermoregulatory responses of feeder cattle. J. Therm. Biol., 28(2), 149-157. https://doi.org/10.1016/S0306-4565(02)00052-9

Brutsaert, W. (1975). On a derivable formula for long-wave radiation from clear skies. Water Resour. Res., 11(5), 742-744. https://doi.org/10.1029/WR011i005p00742

Cena, K., & Monteith, J. L. (1975). Transfer processes in animal coats. III. Water vapour diffusion. Proc. Royal Soc. London. Series B, Biol. Sci., 188(1093), 413-423. https://doi.org/10.1098/rspb.1975.0028

Chen, J. M., Schtz, K. E., & Tucker, C. B. (2015). Cooling cows efficiently with sprinklers: Physiological responses to water spray. J. Dairy Sci., 98(10), 6925-6938. https://doi.org/10.3168/jds.2015-9434

Collier, R. J., & Gebremedhin, K. G. (2015). Thermal biology of domestic animals. Annu. Rev. Anim. Biosci., 3, 513-532. https://doi.org/10.1146/annurev-animal-022114-110659

Cooper, C. D., & Alley, F. C. (2011). Air pollution control: A design approach (4th ed.). Long Grove, IL: Waveland Press.

Crawford, T. M., & Duchon, C. E. (1999). An improved parameterization for estimating effective atmospheric emissivity for use in calculating daytime downwelling longwave radiation. J. Appl. Meteorol. Climatol., 38(4), 474-480. https://doi.org/10.1175/1520-0450(1999)038<0474:AIPFEE>2.0.CO;2

Gatenby, R. M. (1986). Exponential relation between sweat rate and skin temperature in hot climates. J. Agric. Sci., 106(1), 175-183. https://doi.org/10.1017/S0021859600061888

Gebremedhin, K. G., Lee, C. N., Hillman, P. E., & Collier, R. J. (2010). Physiological responses of dairy cows during extended solar exposure. Trans. ASABE, 53(1), 239-247. https://doi.org/10.13031/2013.29499

Gorczyca, M. T., & Gebremedhin, K. G. (2020). Ranking of environmental heat stressors for dairy cows using machine learning algorithms. Comput. Electron. Agric., 168, 105124. https://doi.org/10.1016/j.compag.2019.105124

Gwadera, M., Larwa, B., & Kupiec, K. (2017). Undisturbed ground temperature Different methods of determination. Sustainability, 9(11), 2055. https://doi.org/10.3390/su9112055

Kabuga, J. D. (1992). The influence of thermal conditions on rectal temperature, respiration rate and pulse rate of lactating Holstein-Friesian cows in the humid tropics. Int. J. Biometeorol., 36, 146-150. https://doi.org/10.1007/BF01224817

Key, N., Sneeringer, S., & Marquardt, D. (2014). Climate change, heat stress, and U.S. dairy production: USDA-ERS Economic Research Report Number 175. USDA Economic Research Service. http://doi.org/10.2139/ssrn.2506668

Li, G., Chen, S., Chen, J., Peng, D., & Gu, P. (2020). Predicting rectal temperature and respiration rate responses in lactating dairy cows exposed to heat stress. J. Dairy Sci., 103(6), 5466-5484. https://doi.org/10.3168/jds.2019-16411

Maia, A. S., daSilva, R. G., & Battiston Loureiro, C. M. (2005). Sensible and latent heat loss from the body surface of Holstein cows in a tropical environment. Int. J. Biometeorol., 50(1), 17-22. https://doi.org/10.1007/s00484-005-0267-1

Martello, L. S., Savastano Junior, H., Silva, S. L., & Balieiro, J. C. (2010). Alternative body sites for heat stress measurement in milking cows under tropical conditions and their relationship to the thermal discomfort of the animals. Int. J. Biometeorol., 54, 647-652. https://doi.org/10.1007/s00484-009-0268-6

Mauger, G., Bauman, Y., Nennech, T., & Salath, E. (2015). Impacts of climate change on milk production in the United States. Prof. Geogr., 67(1), 121-131. https://doi.org/10.1080/00330124.2014.921017

McArthur, A. J. (1987). Thermal interaction between animal and microclimate: A comprehensive model. J. Theor. Biol., 126(2), 203-238. https://doi.org/10.1016/S0022-5193(87)80229-1

McGovern, R. E., & Bruce, J. M. (2000). AP Animal production technology: A model of the thermal balance for cattle in hot conditions. J. Agric. Eng. Res., 77(1), 81-22. https://doi.org/10.1006/jaer.2000.0560

Mondaca, M. R., Rojano, F., Choi, C. Y., & Gebremedin, K. G. (2013). A conjugate heat and mass transfer model to evaluate the efficiency of conductive cooling for dairy cattle. Trans. ASABE, 56(6), 1471-1482. https://doi.org/10.13031/trans.56.10178

Renaudeau, D., Collin, A., Yahav, S., de Basilio, V., Gourdine, J. L., & Collier, R. J. (2012). Adaption to hot climate and strategies to alleviate heat stress in livestock production. Animal, 6(5), 707-728. https://doi.org/10.1017/S1751731111002448

Sharpe, K. T., Heins, B. J., Buchanan, E. S., & Reese, M. H. (2021). Evaluation of solar photovoltaic systems to shade cows in a pasture-based dairy herd. J. Dairy Sci., 104, 2794-2806. https://doi.org/10.3168/jds.2020-18821

Shoshani, E., & Hetzroni, A. (2013). Optimal barn characteristics for high-yielding Holstein cows as derived by a new heat-stress model. Animal, 7(1), 176-182. https://doi.org/10.1017/S1751731112001085

St-Pierre, N. R., Cobanov, B., & Schnitkey, G. (2003). Economic losses from heat stress by US livestock industries. J. Dairy Sci., 86(Sup), E52-E77. https://doi.org/10.3168/jds.S0022-0302(03)74040-5

Sutherland, W. (1893). LII. The viscosity of gases and molecular force. London, Edinburgh, Dublin Philos. Mag. J. Sci., 36(223), 507-531. https://doi.org/10.1080/14786449308620508

Thomas, C. K., & Pearson, R. A. (1986). Effects of ambient temperature and head cooling on energy expenditure, food intake and heat tolerance of Brahman and Brahman Friesian cattle working on treadmills. Anim. Sci., 43(1), 83-90. https://doi.org/10.1017/S0003356100018353

Thompson, V. A., Barioni, L. G., Rumsey, T. R., Fadel, J. G., & Sainz, R. D. (2014). The development of a dynamic, mechanistic, thermal balance model for Bos indicus and Bos taurus. J. Agric. Sci., 152(3), 464-482. https://doi.org/10.1017/S002185961300049X

Thompson, V. A., Fadel, J. G., & Sainz, R. D. (2011). Meta-analysis to predict sweating and respiration rates for Bos indicus, Bos taurus, and their crossbreds. J. Anim. Sci., 89(12), 3973-3982. https://doi.org/10.2527/jas.2011-3913

Turnpenny, J. R., McArthur, A. J., Clark, J. A., & Wathes, C. M. (2000a). Thermal balance of livestock: 1. A parsimonious model. Agric. For. Meteorol., 101(1), 15-27. https://doi.org/10.1016/S0168-1923(99)00159-8

Turnpenny, J. R., Wathes, C. M., Clark, J. A., & McArthur, A. J. (2000b). Thermal balance of livestock: 2. Applications of a parsimonious model. Agric. For. Meteorol., 101(1), 29-52. https://doi.org/10.1016/S0168-1923(99)00157-4

Wang, X., Gao, H., Gebremedhin, K. G., Schmidt Berg, B., Van Os, J., Tucker, C. B., & Zhang, G. (2018a). A predictive model of equivalent temperature index for dairy cattle (ETIC). J. Therm. Biol., 76, 165-170. https://doi.org/10.1016/j.jtherbio.2018.07.013

Wang, X., Gao, H., Gebremedhin, K. G., Schmidt Berg, B., Van Os, J., Tucker, C. B., & Zhang, G. (2019). Corrigendum to A predictive model of equivalent temperature index for dairy cattle (ETIC). J. Therm. Biol., 82, 252-253. https://doi.org/10.1016/j.jtherbio.2018.12.012

Wang, X., Zhang, G., & Choi, C. Y. (2018b). Effect of airflow speed and direction on convective heat transfer of standing and reclining cows. Biosyst. Eng., 167, 87-98. https://doi.org/10.1016/j.biosystemseng.2017.12.011

Welty, J. R., Wicks, C. E., Wilson, R. E., & Rorrer, G. L. (2008). Fundamentals of momentum, heat and mass transfer. Hoboken, NJ: John Wiley & Sons, Inc.

West, J. W. (2003). Effects of heat-stress on production in dairy cattle. J. Dairy Sci., 86(6), 2131-2144. https://doi.org/10.3168/jds.S0022-0302(03)73803-X

Zhou, B., Wang, X., Mondaca, M. R., Rong, L., & Choi, C. Y. (2019). Assessment of optimal airflow baffle locations and angles in mechanically-ventilated dairy houses using computational fluid dynamics. Comput. Electron. Agric., 165, 104930. https://doi.org/10.1016/j.compag.2019.104930