This document provides an overview of an R implementation of the NicheMapR endotherm model, via the sub-functions of ‘endoR’ and ‘endoR_devel’. This implementation modularises the original version such that many of the internal subroutines and functions can be called separately. The idea is that the user has more flexibility in setting the model up for the behavioural and physiological responses of a particular animal. All the underlying theory, equations and subroutines are described in Endotherm Theory Tutorial. Also the Endotherm Model Tutorial shows how to use ‘endoR’ (and ‘endoR_devel’) as a whole.
library(NicheMapR)
library(knitr)
library(plotrix)
This subroutine computes the three parameters needed for computing conduction and infrared radiation through the fur.
# environmental input
<- 20 # air temperature, for calculation of conductivity of air (°C)
TA
# shape input
<- 0.3 # maximum fraction of surface area that is ventral (fractional, 0-1)
PVEN
# fur properties
<- 3e-05 # hair diameter, dorsal (m)
DHAIRD <- 3e-05 # hair diameter, ventral (m)
DHAIRV <- 0.0239 # hair length, dorsal (m)
LHAIRD <- 0.0239 # hair length, ventral (m)
LHAIRV <- 0.009 # fur depth, dorsal (m)
ZFURD <- 0.009 # fur depth, ventral (m)
ZFURV <- ZFURV # depth of compressed fur (for conduction) (m)
ZFURCOMP <- 39680000 # hair density, dorsal (1/m2)
RHOD <- 27810000 # hair density, ventral (1/m2)
RHOV <- 0.301 # fur reflectivity dorsal (fractional, 0-1)
REFLD <- 0.301 # fur reflectivity ventral (fractional, 0-1)
REFLV <- 0.209 # hair thermal conductivity (W/m°C)
KHAIR
# call the subroutine
<- IRPROP(TA, DHAIRD, DHAIRV, LHAIRD, LHAIRV, ZFURD,
IRPROP.out
ZFURV, RHOD, RHOV, REFLD, REFLV, ZFURCOMP, PVEN, KHAIR)
# output
<- IRPROP.out[1:3] # effective thermal conductivity of fur array, mean, dorsal, ventral (W/mK)
KEFARA <- IRPROP.out[4:6] # term involved in computing optical thickness (1/m)
BETARA <- IRPROP.out[7:9] # optical thickness array, mean, dorsal, ventral (-)
B1ARA <- IRPROP.out[10:12] # fur diameter array, mean, dorsal, ventral (m)
DHAR <- IRPROP.out[13:15] # fur length array, mean, dorsal, ventral (m)
LHAR <- IRPROP.out[16:18] # fur density array, mean, dorsal, ventral (fibers/m2)
RHOAR <- IRPROP.out[19:21] # fur depth array, mean, dorsal, ventral (m)
ZZFUR <- IRPROP.out[22:24] # fur reflectivity array, mean, dorsal, ventral (fractional, 0-1)
REFLFR <- IRPROP.out[25] # test of presence of fur (length x diameter x density x depth) (-)
FURTST <- IRPROP.out[26] # effective thermal conductivity of compressed ventral fur (W/mK)
KFURCMPRS
<- c("KEFARA mean", "KEFARA dorsal", "KEFARA ventral",
IRPROP.lab "BETARA mean", "BETARA dorsal", "BETARA ventral", "B1ARA mean",
"B1ARA dorsal", "B1ARA ventral", "DHAR mean", "DHAR dorsal",
"DHAR ventral", "LHAR mean", "LHAR dorsal", "LHAR ventral",
"RHOAR mean", "RHOAR dorsal", "RHOAR ventral", "ZZFUR mean",
"ZZFUR dorsal", "ZZFUR ventral", "REFLFR mean", "REFLFR dorsal",
"REFLFR ventral", "FURTST", "KFURCMPRS")
kable(cbind(IRPROP.lab, IRPROP.out[1:26]))
IRPROP.lab | |
---|---|
KEFARA mean | 0.0293068804352327 |
KEFARA dorsal | 0.0296911080500765 |
KEFARA ventral | 0.0289267341187979 |
BETARA mean | 613.67390356235 |
BETARA dorsal | 674.176485875967 |
BETARA ventral | 553.171321248733 |
B1ARA mean | 5.52306513206115 |
B1ARA dorsal | 6.0675883728837 |
B1ARA ventral | 4.9785418912386 |
DHAR mean | 3e-05 |
DHAR dorsal | 3e-05 |
DHAR ventral | 3e-05 |
LHAR mean | 0.0239 |
LHAR dorsal | 0.0239 |
LHAR ventral | 0.0239 |
RHOAR mean | 36119000 |
RHOAR dorsal | 39680000 |
RHOAR ventral | 32558000 |
ZZFUR mean | 0.009 |
ZZFUR dorsal | 0.009 |
ZZFUR ventral | 0.009 |
REFLFR mean | 0.301 |
REFLFR dorsal | 0.301 |
REFLFR ventral | 0.301 |
FURTST | 0.25605504 |
KFURCMPRS | 0.0289267341187979 |
Example analysis of fur diameter on effective thermal conductivity of fur. Note this does not include the radiative component of fur conductivity, which is calculated in SIMULSOL. Empirical measurements of fur thermal conductivity will include both the radiative and conductive aspects.
<- seq(0, 150, 2) # hair diameters (micrometers)
DHAIRs <- NULL
KEFARAs for (i in 1:length(DHAIRs)) {
<- IRPROP(TA, DHAIRs[i] * 1e-06, DHAIRs[i] * 1e-06,
KEFARAs[i]
LHAIRD, LHAIRV, ZFURD, ZFURV, RHOD, RHOV, REFLD, REFLV,1]
ZFURCOMP, PVEN, KHAIR)[
}plot(KEFARAs ~ DHAIRs, type = "p", pch = 16, ylab = "effective fur conductivity, W K-1 m-1",
xlab = "hair diameter, um")
Example analysis of fur density on effective thermal conductivity of fur.
<- seq(0, 50000, 500) # hair densities (1/cm2)
RHOs <- NULL
KEFARAs for (i in 1:length(RHOs)) {
<- IRPROP(TA, DHAIRD, DHAIRV, LHAIRD, LHAIRV,
KEFARAs[i] * 10000, RHOs[i] * 10000, REFLD,
ZFURD, ZFURV, RHOs[i] 1]
REFLV, ZFURCOMP, PVEN, KHAIR)[
}plot(KEFARAs ~ RHOs, type = "p", pch = 16, ylab = "effective fur conductivity, W K-1 m-1",
xlab = "hair density, 1/cm2")
Example analysis of fur depth on effective thermal conductivity of fur.
<- seq(0, 50, 1) # hair depths (mm)
ZFURs <- NULL
KEFARAs for (i in 1:length(ZFURs)) {
<- IRPROP(TA, DHAIRD, DHAIRV, LHAIRD, LHAIRV,
KEFARAs[i] * 0.001, ZFURs[i] * 0.001, RHOD, RHOV, REFLD,
ZFURs[i] 1]
REFLV, ZFURCOMP, PVEN, KHAIR)[
}plot(KEFARAs ~ ZFURs, type = "p", pch = 16, ylab = "effective fur conductivity, W K-1 m-1",
xlab = "fur depth, mm")
This subroutine computes aspects of the animal’s geometry. Note that GEOM_ENDO requires outputs from IRPROP on fur properties. Note also that this routine is called ‘allom’ in past versions of ‘Niche Mapper’.
# input
<- 10 # animal mass (kg)
AMASS <- 1000 # animal density (kg/m3)
ANDENS <- 1 # is subcutaneous fat present? (0 is no, 1 is yes)
SUBQFAT <- 20 # body fat (%)
FATPCT <- 4 # shape, 1 is cylinder, 2 is sphere, 3 is plate, 4 is ellipsoid
SHAPE <- 2.7 # current ratio between long and short axis (-)
SHAPE_B <- SHAPE_B # current ratio of length:height (for plate)
SHAPE_C <- DHAR[1] # fur diameter, mean (m) (from IRPROP)
DHARA <- RHOAR[1] # hair density, mean (1/m2) (from IRPROP)
RHOARA <- ZZFUR[1] # fur depth, mean (m) (from IRPROP)
ZFUR <- 0 # fraction of body in contact with substrate (fractional, 0-1)
PCOND <- 0 # if 1, uses bird skin surface area scaling from Walsberg
SAMODE # and King 1978. JEB Biology 76:185–189, if 2, uses mammal
# surface area scaling from Stahl (1967) J. of App.
# Physiology, 453–460.
<- 0 # if 1, largest surface area normal to sun's ray's, if 2, largest surface parallel to sun's rays, if 0, average of normal/parallel posture
ORIENT <- 20 # zenith angle of the sun
Z <- pi/180 * Z # convert degrees to radians
ZEN
# call the subroutine
<- GEOM_ENDO(AMASS, ANDENS, FATPCT, SHAPE, ZFUR, SUBQFAT,
GEOM.out
SHAPE_B, SHAPE_C, DHARA, RHOARA, PCOND, SAMODE, ORIENT, ZEN)
# output
<- GEOM.out[1] # volume (m3)
VOL <- GEOM.out[2] # characteristic dimension for convection (m)
D <- GEOM.out[3] # mass body fat (kg)
MASFAT <- GEOM.out[4] # volume body fat (m3)
VOLFAT <- GEOM.out[5] # length (m)
ALENTH <- GEOM.out[6] # width (m)
AWIDTH <- GEOM.out[7] # height (m)
AHEIT <- GEOM.out[8] # total area at fur/feathers-air interface (m2)
ATOT <- GEOM.out[9] # silhouette area to use in solar calcs (m2) may be normal, parallel
ASIL # or average set via ORIENT
<- GEOM.out[10] # silhouette area normal to sun (m2)
ASILN <- GEOM.out[11] # silhouette area parallel to sun (m2)
ASILP <- GEOM.out[12] # mass (g)
GMASS <- GEOM.out[13] # area of skin (m2)
AREASKIN <- GEOM.out[14] # flesh volume (m3)
FLSHVL <- GEOM.out[15] # fat layer thickness (m)
FATTHK <- GEOM.out[16] # semimajor axis length (m)
ASEMAJ <- GEOM.out[17] # b semiminor axis length (m)
BSEMIN <- GEOM.out[18] # c semiminor axis length (m)
CSEMIN # (currently only prolate spheroid)
<- GEOM.out[19] # area of skin for evaporation (total skin area - hair area) (m2)
CONVSK <- GEOM.out[20] # area for convection (total area minus ventral area in contact
CONVAR # with the substrate, as determined by PCOND), m2
<- GEOM.out[21] # shape-specific core-skin radius in shortest dimension (m)
R1 <- GEOM.out[22] # shape-specific core-fur radius in shortest dimension (m)
R2
<- c("VOL", "D", "MASFAT", "VOLFAT", "ALENTH", "AWIDTH",
GEOM.lab "AHEIT", "ATOT", "ASIL", "ASILN", "ASILP", "GMASS", "AREASKIN",
"FLSHVL", "FATTHK", "ASEMAJ", "BSEMIN", "CSEMIN", "CONVSK",
"CONVAR", "R1", "R2")
kable(cbind(GEOM.lab, t(GEOM.out)))
GEOM.lab | |
---|---|
VOL | 0.01 |
D | 0.215443459147026 |
MASFAT | 2 |
VOLFAT | 0.00221975582685905 |
ALENTH | 0.529649771785497 |
AWIDTH | 0.196166582142777 |
AHEIT | 0.196166582142777 |
ATOT | 0.306057341103935 |
ASIL | 0.0640710561559047 |
ASILN | 0.0921179995701185 |
ASILP | 0.0360241127416909 |
GMASS | 10000 |
AREASKIN | 0.269773486865498 |
FLSHVL | 0.00778024417314096 |
FATTHK | 0.00980642236948426 |
ASEMAJ | 0.264824885892749 |
BSEMIN | 0.0980832910713883 |
CSEMIN | 0.0980832910713883 |
CONVSK | 0.262885895484012 |
CONVAR | 0.306057341103935 |
R1 | 0.0980832910713883 |
R2 | 0.107083291071388 |
Plot of current animal configuration.
if (SHAPE == 4) {
# ellipsoid
par(mfrow = c(1, 2))
plot(c(0, ASEMAJ * 2 + ZFUR * 2), c(0, ASEMAJ * 2 + ZFUR *
2), type = "n", main = "ellipsoid, sagittal section",
ylab = "minor axis, m", xlab = "major axis, m", asp = 1)
draw.ellipse(ASEMAJ + ZFUR, ASEMAJ + ZFUR, col = "black",
border = "black", a = ASEMAJ + ZFUR, b = BSEMIN + ZFUR)
draw.ellipse(ASEMAJ + ZFUR, ASEMAJ + ZFUR, col = "pink",
border = "pink", a = ASEMAJ, b = BSEMIN)
draw.ellipse(ASEMAJ + ZFUR, ASEMAJ + ZFUR, col = "red", border = "red",
a = ASEMAJ - FATTHK, b = BSEMIN - FATTHK)
plot(c(0, ASEMAJ * 2 + ZFUR * 2), c(0, ASEMAJ * 2 + ZFUR *
2), type = "n", main = "ellipsoid, transverse section",
ylab = "minor axis, m", xlab = "minor axis, m", asp = 1)
draw.ellipse(ASEMAJ + ZFUR, ASEMAJ + ZFUR, col = "black",
border = "black", a = BSEMIN + ZFUR, b = CSEMIN + ZFUR)
draw.ellipse(ASEMAJ + ZFUR, ASEMAJ + ZFUR, col = "pink",
border = "pink", a = BSEMIN, b = CSEMIN)
draw.ellipse(ASEMAJ + ZFUR, ASEMAJ + ZFUR, col = "red", border = "red",
a = BSEMIN - FATTHK, b = CSEMIN - FATTHK)
par(mfrow = c(1, 1))
}
Compute the solar radiation absorbed by different parts of the organism and its environment.
# environmental variables
<- 1000 # solar radiation, horizontal plane (W/m2)
QSOLR <- 0 # shade (%, 0-100)
SHADE <- 20 # zenith angle of sun (degrees from overhead)
Z <- 0.8 # solar absorptivity of substrate (fractional, 0-1)
ABSSB
# traits
<- ATOT # surface area for solar exchange, m2 (from GEOM)
AREA <- 0.8 # solar absorptivity of dorsal fur (fractional, 0-1)
ABSAND <- 0.8 # solar absorptivity of ventral fur (fractional, 0-1)
ABSANV <- ASIL # silhouette area to sun, m2 (from GEOM)
ASIL <- 0.15 # proportion of solar radiation that is diffuse (fractional, 0-1)
PDIF <- 0.4 # configuration factor to sky (-)
FASKY <- 0 # configuration factor to vegetation (-)
FAVEG
# solar radiation normal to sun's rays
<- pi/180 * Z # convert degrees to radians
ZEN if (Z < 90) {
# compute solar radiation on a surface normal to the
# direct rays of the sun
<- cos(ZEN)
CZ <- QSOLR/CZ
QNORM else {
} # diffuse skylight only
<- QSOLR
QNORM
}
<- SOLAR_ENDO(AREA, ABSAND, ABSANV, ABSSB, ASIL, PDIF,
SOLAR.out
QNORM, SHADE, QSOLR, FASKY, FAVEG)
<- SOLAR.out[1] # total (global) solar radiation (W)
QSOLAR <- SOLAR.out[2] # direct solar radiation (W)
QSDIR <- SOLAR.out[3] # diffuse solar radiation from sky (W)
QSSKY <- SOLAR.out[4] # diffuse solar radiation reflected from substrate (W)
QSRSB <- SOLAR.out[5] # total diffuse solar radiation (W)
QSDIFF <- SOLAR.out[6] # dorsal direct solar radiation (W)
QDORSL <- SOLAR.out[7] # ventral diffuse solar radiation (W)
QVENTR
<- c("QSOLAR", "QSDIR", "QSSKY", "QSRSB", "QSDIFF",
SOLAR.lab "QDORSL", "QVENTR")
kable(cbind(SOLAR.lab, t(SOLAR.out)))
SOLAR.lab | |
---|---|
QSOLAR | 90.436692916682 |
QSDIR | 46.3644357977154 |
QSSKY | 14.6907523729889 |
QSRSB | 29.3815047459777 |
QSDIFF | 44.0722571189666 |
QDORSL | 61.0551881707043 |
QVENTR | 29.3815047459777 |
Compute convective heat exchange.
# input
<- 33 # skin temperature (°C)
TS <- 20 # air temperature (°C)
TENV <- 10 # fur/air interface temperature (°C)
TFA <- 4 # shape, 1 is cylinder, 2 is sphere, 3 is plate, 4 is ellipsoid
SHAPE <- CONVAR # surface area for convection, m2 (from GEOM)
SURFAR <- 0 # FLUID TYPE: 0 = AIR; 1 = FRESH WATER; 2 = SALT WATER
FLTYPE <- FURTST # test of fur presence (-) from IRPROP
FURTST <- 1 # wind speed (m/s)
VEL <- 0 # elevation (m)
ELEV <- -1 # barometric pressure (Pa), negative means altitude is used
BP
# run subroutine
<- CONV_ENDO(TS, TENV, SHAPE, SURFAR, FLTYPE, FURTST,
CONV.out
D, TFA, VEL, ZFUR, BP, ELEV)
<- CONV.out[1] # convective heat loss (W)
QCONV <- CONV.out[2] # combined convection coefficient (W/m2 C)
HC <- CONV.out[3] # free convection coefficient (W/m2 C)
HCFREE <- CONV.out[4] # forced convection coefficient (W/m2 C)
HCFOR <- CONV.out[5] # mass transfer coefficient (kg/m2 s)
HD <- CONV.out[6] # free mass transfer coefficient (kg/m2 s)
HDFREE <- CONV.out[7] # forced mass transfer coefficient (kg/m2 s)
HDFORC <- CONV.out[8] # Nusselt number (-)
ANU <- CONV.out[9] # Reynold's number (-)
RE <- CONV.out[10] # Grasshof number (-)
GR <- CONV.out[11] # Prandtl number (-)
PR <- CONV.out[12] # Rayleigh number (-)
RA <- CONV.out[13] # Schmidt number (-)
SC <- CONV.out[14] # barometric pressure (Pa)
BP
<- c("QCONV", "HC", "HCFREE", "HCFOR", "HD", "HDFREE",
CONV.lab "HDFORC", "ANU", "RE", "GR", "PR", "RA", "SC", "BP")
kable(cbind(CONV.lab, t(CONV.out)))
CONV.lab | |
---|---|
QCONV | -42.330368639197 |
HC | 13.8308620490897 |
HCFREE | 4.18905907070396 |
HCFOR | 0 |
HD | 0.0166500411300533 |
HDFREE | 0.00389857961309036 |
HDFORC | 0.0127514615169629 |
ANU | 115.050182725776 |
RE | 14271.7032454331 |
GR | 14679635.9692975 |
PR | 0.712515322862147 |
RA | 10459465.5621628 |
SC | 0.595551256140649 |
BP | 101325 |
Compute evaporative heat exchange.
<- BP # barometric pressure (Pa) (from CONV)
BP <- 20 # air temperature (°C)
TA <- 20 # relative humidity (%)
RELHUM <- 37 # core temperature (°C)
TC <- 33 # skin temperature (°C)
TSKIN <- 11 # part of the skin surface that is wet (%)
PCTWET <- 0 # is flight occurring this hour? (imposes forced evaporative loss)
FLYHR <- 0 # is evaporation partly from bare skin? (0 = no, 1 = yes, % defined with PCTBAREVAP)
BAREVAP <- 0 # surface area for evaporation that is skin, e.g. licking paws (%)
PCTBAREVAP <- 0.03 # surface area made up by the eye (%) - make zero if sleeping
PCTEYES <- ZFUR # fur depth (m)
ZFUR <- 0 # part of the fur surface that is wet (%)
FURWET
<- SEVAP_ENDO(BP, TA, RELHUM, VEL, TC, TSKIN, ELEV,
SEVAP.out
PCTWET, FLYHR, CONVSK, HD, HDFREE, PCTBAREVAP, PCTEYES, ZFUR,
FURWET, TFA, CONVAR)
<- SEVAP.out[1] # skin evaporative heat loss (W)
QSEVAP <- SEVAP.out[2] # ocular evaporation (kg/s)
WEYES <- SEVAP.out[3] # forced cutaneous evaporation (kg/s)
WCUTHF <- SEVAP.out[4] # free cutaneous evaporation (kg/s)
WCUTF <- SEVAP.out[5] # total cutaneous evaporation (kg/s)
WCUT <- SEVAP.out[6] # total fur evaporation (kg/s)
WTFUR <- SEVAP.out[7] # fur evaporative heat loss (W))
QFSEVAP
<- c("QSEVAP", "WEYES", "WCUTHF", "WCUTF", "WCUT",
SEVAP.lab "WTFUR", "QFSEVAP")
kable(cbind(SEVAP.lab, t(SEVAP.out)))
SEVAP.lab | |
---|---|
QSEVAP | 9.01497254207792 |
WEYES | 4.2302327009646e-08 |
WCUTHF | 0 |
WCUTF | 3.63184056945864e-06 |
WCUT | 3.63184056945864e-06 |
WTFUR | 0 |
QFSEVAP | 0 |
Setup and the SIMULSOL algorithm to SIMULtaneously SOLve for skin temperature \(T_s\) and the temperature of the fur/feather interface \(T_fa\) given an environment, shape/posture and skin wetness. SIMULSOL is the core algorithm for computing the heat balance.
# environment
<- 0 # fluid type: 0 = air; 1 = fresh water; 2 = salt water
FLTYPE <- TA # reference air temperature (e.g. at 1.2 or 2m) (°C)
TAREF <- TA # ground temperature (°C)
TLOWER <- TA # sky temperature (°C)
TSKY <- TA # surface temperature for conduction (°C)
TCONDSB <- TA # bush temperature (°C)
TBUSH <- TAREF # assume vegetation casting shade is at reference (e.g. 1.2m or 2m) air temperature (°C)
TVEG <- 20 # relative humidity (%)
RH <- 1 # wind speed (m/s)
VEL <- BP # barometric pressure (Pa), negative means altitude is used (from CONV)
BP <- 0 # elevation (m)
ELEV <- 2.79 # substrate thermal conductivity (W/m°C)
KSUB
# physiology and morphology
<- 0 # part of the skin surface that is wet (%)
PCTWET <- 0.9 # initial thermal conductivity of flesh (0.412 - 2.8 W/mK)
AK1 <- 0.23 # conductivity of fat (W/mK)
AK2 <- 1 # fractional depth of fur at which longwave radiation is exchanged (0-1)
XR <- 0.99 # animal emissivity (-)
EMISAN <- 0 # is evaporation partly from bare skin? (0 = no, 1 = yes, % defined with PCTBAREVAP)
BAREVAP <- 0 # surface area for evaporation that is skin, e.g. licking paws (%)
PCTBAREVAP <- 0.03 # surface area made up by the eye (%) - make zero if sleeping
PCTEYES # behaviour
<- 0 # is flight occurring this hour? (imposes forced evaporative loss)
FLYHR
# configuration factors
<- 0.5 # configuration factor to ground (-)
FAGRD <- 0.5 # configuration factor to sky (-)
FASKY <- 0 # this is for overhead veg (at TAREF) (-)
FAVEG <- 0 # this is for veg below/around animal (at TALOC) (-)
FABUSH
# reference configuration factors
<- FABUSH # nearby bush (-)
FABUSHREF <- FASKY # sky (-)
FASKYREF <- FAGRD # ground (-)
FAGRDREF <- FAVEG # vegetation (-)
FAVEGREF
<- 0 # user-specified fur thermal conductivity (W/m C), not used if 0
FURTHRMK <- 4 # shape, 1 is cylinder, 2 is sphere, 3 is plate, 4 is ellipsoid
SHAPE
# Initial values
<- TC # current guess of skin temperature (deg C)
TS <- TA # current guess of fur/air interface temperature (deg C)
TFA
<- 0.001 # tolerance for SIMULSOL
DIFTOL
# vector to hold the SIMULSOL results for dorsal and
# ventral side
<- matrix(data = 0, nrow = 2, ncol = 16)
SIMULSOL.out
# repeat for each side, dorsal and ventral, of the animal
for (S in 1:2) {
# Calculating solar intensity entering fur. This will
# depend on whether we are calculating the fur
# temperature for the dorsal side or the ventral side.
# The dorsal side will have solar inputs from the
# direct beam hitting the silhouette area as well as
# diffuse solar scattered from the sky. The ventral
# side will have diffuse solar scattered off the
# substrate.
# Setting config factors and solar depending on whether
# the dorsal side (S=1) or ventral side (S=2) is being
# estimated.
if (QSOLAR > 0) {
if (S == 1) {
# proportion of upward view that is sky
<- FASKYREF/(FASKYREF + FAVEGREF)
FASKY # proportion of upward view that is vegetation
# (shade)
<- FAVEGREF/(FASKYREF + FAVEGREF)
FAVEG <- 0
FAGRD <- 0
FABUSH <- 2 * QSDIR + ((QSSKY/FASKYREF) * FASKY) # direct x 2
QSLR # because assuming sun in both directions, and
# un-adjusting QSSKY for config factor imposed
# in SOLAR_ENDO and back to new larger one in
# both directions doing ventral side
else {
} <- 0
FASKY <- 0
FAVEG <- FAGRDREF/(FAGRDREF + FABUSHREF)
FAGRD <- FABUSHREF/(FAGRDREF + FABUSHREF)
FABUSH <- (QVENTR/(1 - FASKYREF - FAVEGREF)) * (1 -
QSLR 2 * PCOND)) # un-adjust by
(# config factor imposed in SOLAR_ENDO to have
# it coming in both directions, but also
# cutting down according to fractional area
# conducting to ground (in both directions)
}else {
} <- 0
QSLR if (S == 1) {
<- FASKYREF/(FASKYREF + FAVEGREF)
FASKY <- FAVEGREF/(FASKYREF + FAVEGREF)
FAVEG <- 0
FAGRD <- 0
FABUSH else {
} <- 0
FASKY <- 0
FAVEG <- FAGRDREF/(FAGRDREF + FABUSHREF)
FAGRD <- FABUSHREF/(FAGRDREF + FABUSHREF)
FABUSH
}
}
# set fur depth and conductivity index for KEFARA, the
# conductivity, is the average (1), front/dorsal (2),
# back/ventral(3) of the body part
if (QSOLR > 0 | ZZFUR[2] != ZZFUR[3]) {
if (S == 1) {
<- ZZFUR[2]
ZL <- KEFARA[2]
KEFF else {
} <- ZZFUR[3]
ZL <- KEFARA[3]
KEFF
}else {
} <- ZZFUR[1]
ZL <- KEFARA[1]
KEFF
}
<- R1 # body radius (including fat), m
RSKIN <- R1 - FATTHK # body radius flesh only (no fat), m
RFLESH <- R1 + ZL # body radius including fur, m
RFUR <- 2 * RFUR # diameter, m
D <- RSKIN + (XR * ZL) # effective radiation radius, m
RRAD <- ALENTH # length, m
LEN
if (SHAPE != 4) {
#! For cylinder and sphere geometries
<- RSKIN + ZFURCOMP
RFURCMP else {
} <- RFUR #! Note that this value is never used if conduction not being modeled,
RFURCMP # but need to have a value for the calculations
}
if (SHAPE == 4) {
#! For ellipsoid geometry
<- BSEMIN + ZFURCOMP
BLCMP else {
} <- RFUR #! Note that this value is never used if conduction not being modeled,
BLCMP # but need to have a value for the calculations
}
# Correcting volume to account for subcutaneous fat
if (SUBQFAT == 1 & FATTHK > 0) {
<- FLSHVL
VOL
}
# Calculating the 'Cd' variable: Qcond =
# Cd(Tskin-Tsub), where Cd = Conduction
# area*ksub/subdepth
if (S == 2) {
<- ATOT * (PCOND * 2)
AREACND <- (AREACND * KSUB)/0.025 # assume conduction happens from 2.5 cm depth
CD else {
} # doing dorsal side, no conduction. No need to
# adjust areas used for convection.
<- 0
AREACND <- 0
CD
}
# package up inputs
<- c(LEN, ZFUR, FURTHRMK, KEFF, BETARA, FURTST, ZL,
FURVARS + 1], DHAR[S + 1], RHOAR[S + 1], REFLFR[S + 1],
LHAR[S
KHAIR, S)<- c(SHAPE, SUBQFAT, CONVAR, VOL, D, CONVAR, CONVSK,
GEOMVARS
RFUR, RFLESH, RSKIN, XR, RRAD, ASEMAJ, BSEMIN, CSEMIN,
CD, PCOND, RFURCMP, BLCMP, KFURCMPRS)<- c(FLTYPE, TA, TS, TBUSH, TVEG, TLOWER, TSKY, TCONDSB,
ENVVARS
RH, VEL, BP, ELEV, FASKY, FABUSH, FAVEG, FAGRD, QSLR)<- c(TC, AK1, AK2, EMISAN, FATTHK, FLYHR, FURWET,
TRAITS
PCTBAREVAP, PCTEYES)
# set IPT, the geometry assumed in SIMULSOL: 1 =
# cylinder, 2 = sphere, 3 = ellipsoid
if (SHAPE %in% c(1, 3)) {
<- 1
IPT
}if (SHAPE == 2) {
<- 2
IPT
}if (SHAPE == 4) {
<- 3
IPT
}
# call SIMULSOL
<- SIMULSOL(DIFTOL, IPT, FURVARS, GEOMVARS,
SIMULSOL.out[S, ]
ENVVARS, TRAITS, TFA, PCTWET, TS)
}
<- cbind(c(1, 2), SIMULSOL.out)
SIMULSOL.out colnames(SIMULSOL.out) <- c("SIDE", "TFA", "TSKIN", "QCONV",
"QCOND", "QGENNET", "QSEVAP", "QRAD", "QSLR", "QRSKY", "QRBSH",
"QRVEG", "QRGRD", "QFSEVAP", "NTRY", "SUCCESS", "KFUR")
<- t(SIMULSOL.out)
tSIMULSOL.out colnames(tSIMULSOL.out) <- c("DORSAL", "VENTRAL")
kable(tSIMULSOL.out)
DORSAL | VENTRAL | |
---|---|---|
SIDE | 1.0000000 | 2.0000000 |
TFA | 39.4826994 | 30.4086840 |
TSKIN | 37.8336962 | 34.7073375 |
QCONV | 83.1226823 | 44.1759445 |
QCOND | 0.0000000 | 0.0000000 |
QGENNET | -1.6448854 | 4.5234310 |
QSEVAP | 0.1423540 | 0.1150878 |
QRAD | 37.2004547 | 18.9954083 |
QSLR | 122.1103763 | 58.7630095 |
QRSKY | 37.2004547 | 0.0000000 |
QRBSH | 0.0000000 | 0.0000000 |
QRVEG | 0.0000000 | 0.0000000 |
QRGRD | 0.0000000 | 18.9954083 |
QFSEVAP | 0.0000000 | 0.0000000 |
NTRY | 1.0000000 | 1.0000000 |
SUCCESS | 1.0000000 | 1.0000000 |
KFUR | 0.0461079 | 0.0436035 |
Finally, find a value of QRESP (heat lost by respiration) to balance the heat budget by solving RESPFUN with ZBRENT.
# define basal metabolic rate
<- (70 * AMASS^0.75) * (4.185/(24 * 3.6)) # heat generation (W), from Kleiber 1947
QBASAL <- 0 # offset between air temperature and breath (°C)
DELTAR <- 20.95 # oxygen concentration of air (%)
O2GAS <- 79.02 # nitrogen concentration of air (%)
N2GAS <- 0.0412 # carbon dioxide concentration of air (%)
CO2GAS <- CO2GAS/100 # reference atmospheric dioxide concentration of air (proportion),
R_PCO2 # to allow for anthropogenic change (%)
<- 0.8 # respiratory quotient (fractional, typically 0.7 (fats)
RQ # to 1 (carbs), with 0.8 typical for protein)
<- 20 # O2 extraction efficiency (%)
EXTREF <- 100 # relative humidity of exhaled air (%)
RELXIT <- 1 # multiplier on breathing rate (-)
PANT
# Now compute a weighted mean heat generation for all the
# parts/components = (dorsal value *(FASKY+FAVEG))+(ventral
# value*FAGRD)
<- SIMULSOL.out[1, 6]
GEND <- SIMULSOL.out[2, 6]
GENV <- FASKYREF + FAVEGREF
DMULT <- 1 - DMULT # Assume that reflectivity of veg below equals reflectivity of soil
VMULT # so VMULT left as 1 - DMULT
<- GEND * DMULT + GENV * VMULT # weighted estimate of metabolic heat generation
X
# lung temperature and temperature of exhaled air
<- (SIMULSOL.out[1, 2] + SIMULSOL.out[2, 2]) * 0.5
TS <- (TC + TS) * 0.5 # average of skin and core
TLUNG <- min(TA + DELTAR, TLUNG) # temperature of exhaled air, °C
TAEXIT <- QBASAL
QMIN <- X - (5 * QMIN)
QM1 <- X + (10 * QMIN)
QM2 <- X
QSUM <- AMASS * 0.01
TOL <- c(TA, O2GAS, N2GAS, CO2GAS, BP, QMIN, RQ, TLUNG,
ZBRENT.in 1, TAEXIT, QSUM, PANT, R_PCO2)
GMASS, EXTREF, RH, RELXIT,
# call ZBRENT subroutine which calls RESPFUN
<- ZBRENT_ENDO(QM1, QM2, TOL, ZBRENT.in)
ZBRENT.out colnames(ZBRENT.out) <- c("RESPFN", "QRESP", "GEVAP", "PCTO2",
"PCTN2", "PCTCO2", "RESPGEN", "O2STP", "O2MOL1", "N2MOL1",
"AIRML1", "O2MOL2", "N2MOL2", "AIRML2", "AIRVOL")
<- t(ZBRENT.out)
tZBRENT.out colnames(tZBRENT.out) <- c("OUTPUT")
kable(tZBRENT.out)
OUTPUT | |
---|---|
RESPFN | -0.0795505 |
QRESP | 1.3829509 |
GEVAP | 0.0003645 |
PCTO2 | 0.2093880 |
PCTN2 | 0.7902000 |
PCTCO2 | 0.0004120 |
RESPGEN | 2.7426732 |
O2STP | 0.0010124 |
O2MOL1 | 0.0002259 |
N2MOL1 | 0.0008523 |
AIRML1 | 0.0010783 |
O2MOL2 | 0.0001807 |
N2MOL2 | 0.0008523 |
AIRML2 | 0.0010692 |
AIRVOL | 0.0241683 |