Predicting Topographic Effect Multipliers in Complex Terrain With Shallow Neural Networks

This study applies computationally efficient shallow neural networks to predict topographic effect multipliers directly from digital elevation data obtained from complex terrain, such as mountainous areas. Data were obtained from boundary layer wind tunnel (BLWT) modeling of surface wind flow over six regions in mainland Puerto Rico and its municipal islands. The results demonstrate an improvement over linear regression models, even for computationally efficient low neuron count and single hidden layer models. The paper proposes the development of a global BLWT data atlas to inform development of methods to predict topographic wind speedup for a diverse range of topography and surface roughness conditions. It also identifies knowledge gaps that could prevent standardization of data collected from different BLWT experimental designs.


INTRODUCTION
Prediction of wind speed in complex topography is integral to the design and operation of civil infrastructure (e.g., Baker et al., 1985;Huang and Xu, 2013;Ngo and Letchford, 2009), siting and optimization of wind energy resources (Barthelmie et al., 2009;Badger et al., 2014;Wilczak et al., 2015;Brahimi, 2019;Haupt et al., 2019;Olson et al., 2019;Santoni et al., 2020), forest management (e.g., Quill et al., 2020), and meteorological forecasting (e.g., Finnigan et al., 2020;Donadio et al., 2021). Within the context of wind loads on structures, surface gusts occurring in a flat expanse can double in magnitude upslope (Kondo et al., 2002), quadrupling the loads occurring on low-rise buildings. Therefore, modern wind load provisions (e.g., ASCE 7-16) account for the presence of hills, escarpments, and other topographic features through provisioning of rational engineering analysis for simple landforms or special wind region maps for areas that have complex topography.
Additional methods to account for complex terrain in structural design include wind tunnel testing and computational fluid dynamics (CFD) simulation. While over time CFD will likely become the dominant method to predict speedup, these methods are still computationally expensive, time intensive, and occasionally prone to not resolving flows in complex terrain well. Thus, this study explores the use of shallow neural networks, which are excellent universal function approximators when trained properly (Hornik et al., 1989).
Flow over isolated topographic features has been studied for decades (e.g., Taylor and Teunissen, 1987;Salmon et al., 1998), later extending to the study of complex landforms (e.g., Chock and Cochran, 2005;Rasouli et al., 2009;Kikuchi and Ishihara, 2012;McAuliffe and Larose, 2012;Lange, 2016). The present research adds to the body of work that uses data-driven methods to predict speedup, using six mountainous regions in mainland Puerto Rico and the municipal islands of Vieques and Culebra as the study testbed ( Figure 1). Multiple shallow neural network configurations are applied to surface velocity measurements obtained from boundary layer wind tunnel (BLWT) modeling of complex terrain, including coastal-tomountain transitions and inner mountain ranges. The results support the efficacy of this approach for rapid and computationally inexpensive prediction of wind speedup in complex topography using digital elevation model (DEM) data. We conclude that BLWT topographic modeling data obtained from different facilities should be warehoused in a single, openly accessible data atlas to validate, improve, and cross-compare the results of predictive models for a wide range of terrain and land cover types.
Motivating the study is the impact of Hurricane Maria (2017), which struck the island as a Saffir Simpson Hurricane Wind Scale Category 4 storm that caused more than 4,000 deaths (Kishore et al., 2018), $90B in economic loss (Pasch et al., 2019), and widespread power outages that occurred for more than 6 months (Kwasinski et al., 2019). The National Hurricane Center (Pasch et al., 2019) speculated that "winds of category 5 intensity were almost certainly felt at some elevated locations on the island," which was further evidenced by post-storm damage assessments led by the FEMA Mitigation Assessment Team (Bass et al., 2018) and Prevatt et al. (2018). The testbed encompasses many of the most heavily impacted areas, as data were collected at the request of the FEMA Building Science Branch for the development of a special wind region map that was later adopted by the 2018 Puerto Rico Construction Code. The selection of these six regions during this study considered Hurricane Maria's path, location of urban areas, and terrain complexity parameters (i.e., elevation, slope, curvature, and rugosity) calculated using ArcGIS and DEM Surface Tools for ArcGIS (Jenness, 2013).

BACKGROUND Definition of Wind Speedup and Topographic Effect Multipliers
For structural design, topographic wind speedup effects are typically quantified as a non-dimensional multiplier that amplifies the velocity pressure used in the calculation of localized wind loads on the building's surface. Ngo and Letchford (2009) provide an overview of four such codified procedures from four internationally recognized building codes. As shown in Figure 2, the structural load multipliers are derived from a ratio that relates the gust wind speed (û zt ) at elevation z above the site of interest to the gust wind speed (û z ) occurring in a flat expanse under the same wind conditions. This ratio is termed the topographic effect multiplier (M zt ), which is equivalent to 1 + the FIGURE 1 | Regions of Puerto Rico selected for testing at BLWT at UF.
1 These data will be made publicly available prior to publication of the manuscript (remove this footnote later and update URL in text with more specific link via DOI). Requests for the data can be made to the corresponding author prior to its release. fractional speedup ratio (Δû z /û z ) given elsewhere in the literature.

Methods to Predict Topographic Wind Speedup in Strong Winds
Approaches to predict topographic wind speedup include analytical methods (e.g., Jackson and Hunt, 1975;Deaves, 1980), physical simulation through boundary layer wind tunnel modeling (e.g., Glanville and Kwok, 1997;Chock and Cochran, 2005;Lubitz and White, 2007;Rasouli et al., 2009;Kikuchi and Ishihara, 2012;McAuliffe and Larose, 2012;Lange, 2016), and CFD (e.g., Abdi and Bitsuamlak, 2014). The two physical simulation efforts most relevant to this study are Chock and Cochran (2005) and McAuliffe and Larose (2012). Chock and Cochran (2005) conducted a "microzonation" study to predict windspeed up factors for use in wind load provisioning for the Hawaiian Islands by applying linear regression to hot film anemometer measurements taken over a 1:6000 scale model of Oahu and its corresponding DEM data. The current study adopts a similar approach; however, it applies a neural network (i.e., logistic regression) for prediction. McAuliffe and Larose (2012) modeled a wind farm at 1:1500 scale to study aerodynamic scaling effects; its experimental design heavily informed the current study. Finnigan et al. (2020) provides a historical narrative on the development of analytical and numerical methods and supporting experiments in BLWTs and in the field, e.g., Askervein Hill (Taylor and Teunissen, 1987) and Perdigão (Mann et al., 2017). Here we turn to describing recent applications of machine learning approaches to improve wind speed prediction.
In the last two decades, artificial neural networks have found a broad range of applications in the study of boundary layer flows. Selected examples include the design and operation of wind tunnels (e.g., Križan et al., 2015), prediction of design wind speeds for civil infrastructure (e.g., Huang and Xu, 2013), and forecasting for wind energy resources (e.g., Burlando and Meissner, 2017;Chen et al., 2018;Liu et al., 2018;Yu et al., 2018;Brahimi, 2019). Starting in the 2000s, machine learning began to be applied in the prediction of topographic wind speedup. One of the first studies was led by Bitsuamlak et al. (2002), who successfully trained a model to predict speedup over shallow and steep sinusoidal hills to match analytical modeling. Bitsuamlak et al. (2007) subsequently incorporated CFD data to quantify speedup accounting for surface roughness conditions. Robert et al. (2012) produced monthly wind speed maps from a data architecture that include multiscale topographic features extracted from a digital elevation model and weather station data. Burlando and Meissner (2017) and Mayo et al. (2018) developed hybrid approaches involving CFD, numerical weather prediction (NWP), and surface wind measurements to inform neural network modeling. More recently, Donadio et al. (2021) developed a highly automated "prediction pipeline" that integrated global and regional NWP to predict wind speed and wind energy production.
These studies foretell of a promising future for regional-scale modeling of wind speedup in complex terrain accounting for the influence of wind directionality or by extension, local climatology. The current work continues to assess the premise that wind speedup in adiabatic boundary layer flow conditions can be predicted directly from elevation data, within the experimental limitations imposed by conducting geometrically scaled physical simulations in a BLWT. The work presented uses data collected over scaled models of complex terrain regions, to develop a feed forward neural network architecture. The literature review provides a summary of the efforts to implement neural networks for wind speed prediction, but no attempt has been made to train a neural network with the goal of predicting topographic wind speedup multipliers over real complex terrain.  Figure 3 for details on the BLWT configuration used during testing. In this study, the roughness elements in the "Terraformer" development section were oriented in the wide configuration to a uniform height of h = 7 mm or 20 mm to produce marine and suburban aerodynamic roughness lengths of z 0 = 0.005 and 0.300 m fullscale (FS) at a geometric scale of 1:3100, respectively. At the chosen geometric scale of 1:3100, matching z 0 to achieve the desired mean velocity profiles effected a reduction in the longitudinal turbulence intensity for the rougher case from an expected value of 25% to nominally 15%. Models were installed in a 4 m turntable that was rotated in one of 16 wind directions (N, NNE, . . . , NNW). The exposure was selected based on the specific upwind terrain of the model section being tested.
Turbulent Flow Instrumentation (TFI) Cobra Probes measured three components of velocity with a sampling frequency of 1250 Hz. Figure 4A shows the dimensions of the probe, which controlled how close to the surface measurements could be made without interfering with the flow. The lowest elevation that could be achieved without introducing speedup between the probe and the model was determined to be approximately three probe head diameters (~8 mm above the tunnel floor or 25 m at full scale). An automated gantry system ( Figure 4B) translated the probes laterally, vertically, and longitudinally using a user-defined sample location input file. Data were then filtered using a 3 rd order Butterworth filter with a 150 Hz cutoff frequency.
Gust values were quantified to conservatively estimate the surface wind field intensity. Given the complex flow conditions, the gust value was obtained from a moving average of the instantaneous wind speed magnitude (s) of the three-dimensional wind speed magnitude (i.e., s where the velocity components are defined as u = longitudinal, v = lateral, and w = vertical) instead of the longitudinal component to ensure that the scalar wind speed maxima was enveloped in the result. Second, the gust was computed asŝ S + gσ s , where S and σ s are the mean  and standard deviation of the instantaneous wind speed magnitude, respectively. The peak factor (g) was set to 3.0, which Balderrama et al. (2012) has shown to be a conservative upper bound for strong winds generated by Atlantic tropical cyclones (TC) landfalling in a variety of surface roughness exposures. The selected value of g is based on a T = 600 s duration record because TCs, as a rule, exhibit non-stationary behavior for longer time periods. The value would need to be re-referenced to T = 3,600 s for straight line wind applications. The automated gantry system and Cobra Probes were also used to measure velocity profiles immediately upwind of the models. Figure 5 presents representative profiles of the longitudinal mean velocity and turbulence intensity for the u, v, and w components. Data measured in the inertial sublayer (ISL) portion of the neutral boundary layer are compared to the logarithmic law modified by Sutton (1949) to account for a zeroplane displacement. Figure 6 presents plots of the longitudinal turbulence spectra measured at the end of the Terraformer in the centerline of the tunnel at heights of z 8 and 988 mm (i.e., fullscale dimensions of z = 25 and 3,063 m, respectively) after applying Welch's method (Bendat et al., 2000) to the unfiltered velocity records. Data are compared to the power spectra model in ESDU 74030 and 74031 (1974), which was first derived by von Kármán (1948) for isotropic turbulence: where f nL x u /U, n is the frequency, L x u is the longitudinal integral length scale of the turbulence, U is the longitudinal mean velocity, S uu (z, n) is the power spectral density function, and σ u is the longitudinal standard deviation of velocity. The data exhibit the same phenomenon identified by McAuliffe and Larose (2012) during topographic modeling of a 1:1,500 scale wind farm in Quebec. Similitude in the freestream flow is achieved, however the observed energy at the high frequency range near the surface is larger than the von Kármán model. McAuliffe and Larose (2012) attributed this phenomenon to eddy surface layer  behavior, where anisotropic spectral decay proportional to the inverse of the wavenumber ( k −1 ) is expected to occur. Although not the subject of this paper, it does draw attention to the validity of the upper bound of geometric scaling requirements, which is generally accepted as 1:5,000 (ASCE 49, 2012).

Topographic Models
Six 1:3,100 scale topographic models of 12.4 km diameter circular regions encompassing the municipalities of Aibonito, Culebra Humacao, Utuado, Vieques, and Yabucoa ( Figure 1) were constructed from DEM data obtained from the U.S. Geological Survey (USGS). Figure 7A shows the Yabucoa model installed in the BLWT turntable. Data represent the bare-earth surface at a resolution of 1/3 arc-second (~10 m) referenced to the North American Vertical Datum of 1983 (NAVD83) High Accuracy Reference Network (HARN) State Plane Puerto Rico Virgin Islands Federal Information Processing Standards (FIPS) 5200. A MultiCam APEX304 3-axis computer numerical control (CNC) routing system fabricated the topographic models by cutting contours into Type IV 25 psi XPS foam panels, which were layered to build the final model. A 30 m spatial filter was applied to the DEM data to ensure the CNC could follow the tool paths. The contour interval height was selected from results of trial-and-error testing of a surface flow over model escarpment routed at different cut heights (smooth, 1, 3, and 5 mm) across the width of the tunnel. Dimensions of the model escarpment were 3.0 m W × 0.3 m H × 1.7 m L. The front section of the model was routed into a curved ramp (rotation angle from 0-45°along a circular arc), and the back section was a tabletop configuration. Ultimately the 3 mm height was chosen because coarser routing produced an unrealistic rise in turbulence intensity over the escarpment and using a finer contour would have been time prohibitive given the project schedule.
The outer regions of the models were modified to minimize the effects of the abrupt transitions in flow caused by the truncation of mountainous areas at the edge of the model. As Figure 7B shows, each model was divided into three regions. Measurements were confined to the inner R = 1 m radius (5 km FS) circular region unless the outer region was flat. No modification was made to the annular "fetch" region (nominal width of R/3) surrounding it, ensuring that surface flows would travel at least 2 km full-scale over unmodified terrain before reaching the measurement location. Finally, any terrain projecting above a line angled 7°upward from the outermost model edge at the tunnel floor was excluded, with the angle of the cutting plane having been chosen to prevent boundary layer

Data Collection and Processing
Surface wind speed measurements were obtained for each region (model) from gridded locations for each of 16 wind directions in 30 s duration records. Table 2 provides a measurement (grid location) count for each region, the corresponding number of records that passed quality control, and the minimum and maximum elevations of the measurement locations. Records were rejected if the Cobra Probe failed to collect at least 80% of the record. Sample quality drops if the resultant wind velocity is outside of the 45°cone of acceptance around the probe head. This generally occurs in separated flow regions where topographic speedup effects are not a concern.
Measurement locations, shown in Figure 7B, were selected to coincide with lines parallel to the 16 wind compass rose aligned with true north (i.e., N, NNE, . . . , NNW), with two lines generally intersecting near a local maxima at the terrain and at least one of those lines following the highest ridgeline. Measurement locations extended between 300 and 4,400 m from the local maxima, generally terminating in a valley or a flat expanse. Spacing between locations varied from 100-500 m.
The testing sequence was as follows: 1) The Terraformer height was adjusted to create a compatible upwind fetch condition by setting the roughness element height (h) to a) h = 20 mm (z 0 = 0.300 m) for regions located fully inside a mountain range (e.g, Aibonito, Utuado) and for the case where terrain that extended beyond the model perimeter. This category of fetch is referred to as "rough" because of both the mountainous terrain and treed conditions that were prevalent in these locations b) h = 7 mm (z 0 = 0.005 m) where the model terminated at the floor before the perimeter (i.e., where upwind marine or flat, open conditions were present). This category of fetch is referred to as the "smooth-to-rough transition" 2) The BLWT fans were brought up to maximum operating speed (20 m/s freestream) 3) Reference velocity measurements were obtained at a height of 8 mm (25 m FS) at five locations spaced 750 mm apart at the trailing edge of the Terraformer development section 4) Data were sequentially collected at the model measurement locations at a height of 8 mm (25 m FS) for a 30 s duration 5) Fan speed was reduced to 0 RPM 6) The turntable was rotated 22.5°to the next wind direction 7) Steps 1-6 were repeated until all 16 wind directions were tested 8) Steps 1-7 were repeated until all six models were tested.

Artificial Neural Networks
Single and double hidden layer shallow feedforward neural networks were constructed using the Deep Learning Toolbox in Matlab 2019a. The model architectures consisted of neurons initialized with random weights and biases. Figure 8 shows the generic single layer neural network architecture configuration used. Hidden layers were fully connected using rectifying linear unit activation functions, and the output layer used a linear activation function. The Levenberg-Marquardt training method (Hagan and Menhaj, 1994) was applied following Liu  (2010) to optimize the weights and biases based on the mean square error. During training, the network was limited to only 1,000 epochs and six validation checks to avoid overfitting the data. Data were separated 70% for training, 15% for testing, and 15% for validation at the beginning of each training process via random assignment. The 4,956 quality-controlled records ( Table 2) measured from 373 unique geographic locations in 16 wind directions were used to train the neural network.
The model input consisted of 2 km of upwind terrain elevation data and 1 km of downwind terrain elevation data at each site in each wind direction. The input was chosen after comparing results of training for different combinations of inputs that included elevation and slope. Inputting the slope as the single variable (as was done in Chock and Cochran, 2005) was also attempted. None of these inputs produced a better trained model than inputting the elevation data alone.
The model output was the gust topographic effect multiplier at each site in each wind direction. As discussed in the Methodology, the choice of calibrating the approach flow to match the mean velocity profile led to a reduction in turbulence intensity for the rough terrain case, leading to the turbulence intensity characteristics representative of open exposure fetch. Given that the M zt values are based on gusts (usually referenced to open exposure), the reference values were not adjusted for this comparison.

RESULTS
Sensitivity analyses were conducted on multiple configurations of shallow neural networks. Figures 9, 10 present results for single and double hidden layer network configurations, respectively. The number of neurons was parametrically increased from 10 to 150 neurons for the single hidden layers and 10 to 80 neurons for the double hidden layers, in steps of 20 neurons. Fine (10 m) and coarse (100 m) horizontal DEM input resolutions were evaluated separately for each network. Results are stratified by the upwind terrain (smooth-to-rough transition, rough) and horizontal fetch resolution. The box and whisker markers characterize coefficients of determination (R 2 ) based on 10 training sets, comparing predicted to measured values of M zt . Complementary values for root-mean-square-error (RMSE) are also given.  In all model configurations, R 2 and RMSE vary little for the lower neuron counts, ranging from approximately 0.65-0.75. Larger neuron counts produce a slight improvement for the rough terrain cases but reduce the accuracy of the model for three of four smooth-to-rough transition cases. This degradation in the model is the worst for the double hidden layer case ( Figure 10A). The only model that shows improvement as the number of neurons increases is the rough terrain case using 100 m horizontal fetch resolution ( Figure 9D), however the marginal improvement is far outweighed by the increased computational expense. For example, training a 150 neuron model takes approximately 30 times longer than training the 10 neuron model. Similarly, the use of a fine (10 m) resolution DEM input data slowed down the training process and did not improve the performance when compared to the use of a coarse (100 m) DEM. Figures 11, 12 present the testing, training, and validation results for two cases of the smooth-to-rough transition and rough terrain models using 30 and 150 neuron single hidden layer architectures, respectively. Data are stratified by the magnitude of the speedup value. Testing the neural network on data yields R 2 values of 0.67 and 0.74 (and RMSE values of 0.11 for both cases), exceeding the observed R 2 value of 0.5 by Chock and Cochran (2005). The neural network performance was compared to the linear regression method, using the same data for the training, and testing of the neural network. The linear model performance on the stratified data resulted in R 2 of 0.65 and 0.61. Notably, the neural network models exhibit less variability between runs as the speedup value increases.
The encouraging results of these sensitivity tests support the use of simple configurations of neural networks for predicting topographic effect multipliers in complex terrain directly from the elevation data. To simplify this approach further, single hidden layer models were trained using data from all terrain types. This modification has the effect of reducing any benefit of stratifying the upwind terrain types but beneficially increases the size of the data available for training, testing, and validation. Figure 13 presents the R 2 and RMSE box and whisker plots based on neuron count for the "all upwind terrain" case. Figure 14 presents the results from training, testing, and validation for a single model. The results yield similar R 2 and RMSE values as the terrain stratified model results given in Figures 9-12, indicating that stratifying the data by upwind terrain type did not produce any appreciable improvement in the modeling.
Returning to the primary goal of producing maximum topographic effect multiplier maps, a 30-neuron model with coarse horizontal fetch resolution and no stratification of approach flow condition was trained (R 2 = 0.819, RMSE = 0.09) and then applied to the six study regions. To create the maps, data were extracted in a 50 m horizontal resolution grid from a gridded interpolant of the DEM data. The same DEM interpolant was then used to extract the fetch elevation data at 100 m intervals (coarse resolution) ranging from 2 km upwind to 1 km downwind for 16 wind angles. For each grid location and each wind angle, the elevation fetch data were input to the neural network to predict the gust topographic effect multiplier specific to that location and wind angle. Finally, the largest multiplier over all wind angles was determined and plotted in the Figure 15.
Ground truth over the study regions is not available to draw a quantitative comparison to the neural network outputs, but qualitatively it is observed that the gradation of the topographic effect multipliers matches expectations. Values are less than unity in valley regions, approach unity in flat open expanses, and reach 1.6-1.7 for the upper slopes of the mountains. The upper bound compares favorably to previous studies of topographic effect multipliers for surface gusts, e.g., Kondo et al. (2002), Glanville and Kwok (1997), and Bowen and Clucas (1992), which give values of 1.6-1.8. Further, the computational expense of producing these maps, which cover 121 km 2 at a 50 m resolution, is approximately 5 mins using a single compute core on a laptop. With parallelization, it should be possible to compute a map for the entire island of Puerto Rico within a day or less.

DISCUSSION AND CONCLUSION
The study provides initial proof-of-concept that shallow neural networks may offer an alternative and computationally inexpensive approach to produce topographic wind speedup maps such as the "special wind regions" in ASCE 7. A shallow neural network was successfully trained using digital elevation data (inputs) and surface velocity measurements taken at 25 m FS over 1:3,100 scale models of regions of Puerto Rico and its municipal islands (outputs). Training a neural network requires voluminous amounts of data, thus data collection was automated to conduct high throughput experimentation, making it possible to collect hundreds of measurements per day. Six 4 m diameter models were tested in BLWT for approximately 2 weeks, producing 5,968 records measured from 373 unique geographic locations in 16 wind directions. Feed forward neural networks were trained varying the number of hidden layers, comparing different horizontal resolutions of the DEM data (input) and by stratifying the approach flow conditions. It was ultimately found that low neuron count, single hidden layer networks are sufficient to train models to predict topographic effect multipliers, which opens the door to applying the model to similar terrain in other parts of Puerto Rico. Sensitivity studies revealed that stratifying the approach flow conditions before training and using as input a finer (10 m) resolution of the DEM data did not cause a significant improvement in the predictive performance of the neural network. Lastly, six topographic wind speedup maps provide predictions over six complex terrain regions, demonstrating the computational advantage of using this methodology.
We conclude that this approach could be extended to other terrain and land-use/land-cover (LULC) types if the BLWT dataset was expanded to other geographical regions, particularly those that possess different gross physical features, taller landforms, and landform elements atypical to the terrain evaluated in this study. Synthetic topography could also be generated from target information about elevation, slope, curvature, and rugosity, with variations in surface roughness introduced to model LULC effects. The results from these studies could be combined to create a shared global topographic wind speedup database that warehouses BLWT experimental data to validate, generalize, and cross-compare predictive modeling and ultimately reduce the need for the extremely laborious physical testing required to collect these data. Alternatively, these data could augment existing field measurement atlases such as the New European Wind Atlas (Mann et al., 2017).
Future research will be required to standardize experimental methods and data processing to achieve this vision. Our experience with this study, a thorough review of the literature, and follow up discussions with other BLWT modelers led us to identify multiple key areas that warrant further investigation. Examples of topics include: • Revisiting the acceptable range of geometric scales. The accepted upper bounds for geometric scaling (~1:2,500-1: 6,000) appear to originate from early BLWT tests on low hills free of obstructions (e.g., trees). While the Reynolds number requirement u p k s /] > 70 can be relaxed without loss of simulation accuracy (Taylor and Teunissen, 1987), it is questionable if the relaxation will hold for more complex environments (Bowen, 1983)-particularly for mountainous landforms that are (partially) covered in tree canopies. Further, the effect of geometric scale on the observed power spectra has received scarce discussion in the literature. As McAuliffe and Larose (2012) and the current (1:3,100 scale) study observed, analysis of the power spectra clearly indicate that surface (10-25 m FS) measurements sample the eddy surface layer, not the inertial sublayer. • Defining a non-arbitrary reference wind speed in the physical simulation. The boundary layer depth varies over the model as a function of the elevation profile, and mountain heights can exceed the boundary layer depth by a factor of five or more. Adopting the surface wind speed upwind of the model also introduces variability between tests of different models. Internal boundary layers form between the end of the development section and the model (causing wind speed recovery), and the presence of the land mass retards the wind speed well before it reaches the landform. • Modeling topographic features that extend to or beyond the perimeter of the model, e.g., a region within a mountain chain. Most of the study subjects in the experimental research literature are isolated landforms that terminate upwind into a flat expanse. While we would assert that the taper-fetch-test partitioning (Figure 7) is a rational approach to prevent an abrupt transition, additional work is required to move this heuristic approach into an analytical framework. This guidance should also address the modelling of the upper part of a tall landform. • Modeling (heterogenous) surface roughness and, as a practical matter, doing so without prohibitively increasing fabrication time/cost. Today topographic models are usually fabricated through subtractive manufacturing, e.g., the 3-axis CNC used in this study. This tried-and-true approach should either be replaced with large-scale additive manufacturing (e.g., 3D printing) or a 4-or 5-axis CNC to simultaneously impart the desired shape and surface roughness.

DATA AVAILABILITY STATEMENT
Data are available through the National Science Foundation Engineering Research Infrastructure DesignSafe site at https:// doi.org/10.17603/ds2-2gab-dw61.

AUTHOR CONTRIBUTIONS
RC led the model fabrication, wind tunnel simulations, and data collection efforts. ARS performed the pilot study to evaluate the efficacy of the machine learning algorithm. JXS-H improved on this algorithm and co-authored the paper, under the supervision of LDA-B and FM. BP scaled up the algorithm for the wind map production and provided input to the manuscript. FM provided overall oversight of the project and co-wrote the manuscript.
FUNDING Support for this research was provided by the NSF EAGER award (No. CMMI-1841979) with additional support for experimentation through the NSF Natural Hazards Engineering Research Infrastructure (NHERI) program (CMMI-1520843). We would also like to acknowledge the NSF Research Experiences for Undergraduates (REU) program for supporting co-author ARS for performing the initial study that led to the work documented herein. Data collection was performed with the support of the Strategic Alliance for Risk Reduction II (STARR II) program, which was funded by the Federal Emergency Management Agency.

ACKNOWLEDGMENTS
The authors wish to recognize the Powell Family Structures and Materials Laboratory staff, with special thanks to Eric Agostinelli, Kevin Stultz, and Shelby Brothers for their contribution in wind tunnel testing. Special thanks are extended to Fernando Carmona-Esteva for his feedback and assistance during the neural network development. Any opinions, findings, and conclusion or recommendations expressed in this paper are those of the authors and do not necessarily reflect the views of the sponsors, partners, and contributors.