Introduction

Dynamical downscaling is a useful approach in weather and climate researches in order to assess the impacts of weather disturbances and extreme weather phenomena at mesoscales and local scales. Because deep cumulus convection plays a critical role in spawning severe weather such as heavy rainfalls and strong winds, a horizontal grid spacing that can adequately resolve the activities and effects of deep convection is required to better represent severe weather phenomena in numerical simulations. As an example of operational numerical weather prediction models, the Japan Meteorological Agency is now conducting regional weather forecasting at the 5 km and 2 km grid spacings. A numerical model with the grid spacing on the order of 1 km was considered to be a cloud-resolving model (CRM) (Grabowski et al. 1996), and Weisman et al. (1997) showed that the non-hydrostatic effects of convection are properly represented at such resolutions. On the other hand, Bryan et al. (2003) demonstrated that the grid spacing on the order of 100 m is necessary to adequately resolve turbulent behaviors of convective motions, and Takemi and Rotunno (2003) indicated an improper spatial filtering in simulating convective clouds at the grid spacing on the order of 1 km from a subgrid-scale-modeling point of view. Thus, the grid spacing of 1 km may not be sufficient to represent convective motion.

From a climate research point of view, Kitoh et al. (2016) emphasized that the importance of high-resolution downscaling in representing extremes lies in reproducing well-resolved topography in numerical models. Takayabu et al. (2016) also listed high-resolution topography as one of the factors that add values to dynamical downscaling. However, most of the previous studies employed grid spacings on the order of 1 km, and the numerical representations at the grid spacing of 100 m have not been well explored. This study deals with the issue on the benefits or added values of dynamical downscaling by conducting a higher-resolution grid spacing on the order of 100 m.

The representation of terrains is critically important in numerically reproducing the spatial and temporal variability of surface winds (Jimenez et al. 2008; Bonnardot and Cautenet 2009; Jimenez et al. 2010; Jimenez and Dudhia 2012). By using the Weather Research and Forecasting (WRF) model, Jimenez and Dudhia (2012) found systematic overestimations of surface wind speeds in valleys and underestimations at mountain tops. This result suggests the importance of terrain representations in wind simulations. Oku et al. (2010) investigated the representations of surface winds during a typhoon passage by comparing two different mesoscale meteorological models with the grid spacing set to 1 km and found that more accurately reproduced terrain features even at the same grid spacing result in capturing more highly variable natures of surface winds. Takemi (2013) conducted regional-scale meteorological simulations of surface winds in complex terrain and demonstrated that the fine-scale spatial variability of winds is reproduced better in the 400 m grid simulation than in the 2 km grid simulation, notifying that the influences of the complex terrain can be found for winds up to the 200 m height above the ground level. In this way, the representation of terrains requires special care in simulating airflows especially in complex terrain.

Strong winds are one of the major meteorological hazards in places with complex terrain features. In 2018, a major typhoon, called Typhoon Jebi (2018), landed in the central part of Japan and spawned severe damages induced by strong winds not only in urban areas but also in mountainous areas (Cabinet Office 2018). Typhoon Jebi (2018) produced extreme winds in a business district of Osaka City which is characterized by complex geometrical features of building distribution and height (Takemi et al. 2019). Much severer winds are considered to have occurred in mountainous areas, because of their higher elevations and more complex terrain features; such severe winds led to serious damages to trees, forest ecosystems, power lines, etc. In the past, major typhoons caused severe damages to trees and forest ecosystems in mountainous areas. For example, due to Typhoon Songda (2004), widespread forest damages occurred in Hokkaido (Sano et al. 2010; Hayashi et al. 2015; Takano et al. 2016; Morimoto et al. 2019). Because major damages to trees by typhoon’s winds have serious influences on forest ecosystems for a long period of time, assessing the impacts of typhoon winds on forest areas and identifying their climate change factors are critically important (Ito et al. 2016; Takemi et al. 2016a; Toda et al. 2018). In this sense, simulating airflows in complex terrain is a scientific challenge.

This study investigates the representations of surface winds in complex terrain in numerical simulations at the 200-m grid spacing. By further downscaling from 1 km to 200 m under the experimental setups of Ito et al. (2016) who examined strong winds due to Typhoon Songda (2004), we explore the benefits of high-resolution numerical simulations of airflows in complex terrain. Identifying the benefits or added values in high-resolution modeling is a key to quantitatively assess the impacts of typhoon winds at local scales.

Numerical model and experimental design

The model used in this study was the WRF model, the Advanced Research WRF (ARW) version 3.6.1 (Skamarock et al. 2008). The model was configured in two-way nested computational domains: the outermost domain (D01) covers the Japanese islands and the surrounding areas in the western North Pacific at the 9 km horizontal grid spacing, downscaled to the area (D02) covering the northeastern part of Japan at the 3 km grid, to the Hokkaido area (D03) at the 1 km grid, and finally to the areas (D04) having the 200 m grid (Fig. 1). For the D04 areas, eight domains were separately placed within D03 in order to save computational costs. The reason why we set these specific regions as D04 is that in these D04 areas, significant damages to forests were reported (Takano et al. 2016; Morimoto et al. 2019) and hence evaluating the numerical representations of the typhoon winds on these forest areas is critically important. The model settings here were exactly the same with those used in Ito et al. (2016) except employing the 200 m grid areas. The number of grid points in each domain is 451 by 451 for D01, 301 by 301 for D02, 532 by 421 for D03, and 166 by 161 for D04. The number of vertical levels was 57, with the top set at the 20 hPa level.

Fig. 1
figure 1

The areas covered in the present numerical simulation. The computational areas for a the 9 km grid domain (D01), the 3 km grid domain (D02), and 1 km grid domain (D03) and b the 200 m grid domains (D04) nested in D03

Full size image

It is important to use high-resolution digital elevation model (DEM) data to set the model terrain used in numerical simulations. The model terrains for D01, D02, and D03 in this study were created with the use of the global DEM data having a horizontal grid spacing of 30 arc-seconds (about 1 km) provided by the United States Geological Survey, while the terrains for D04 were generated with the 50 m horizontal resolution DEM data provided by the Geospatial Information Authority of Japan.

The physics parameterizations chosen here were the same with those used in Ito et al. (2016): the Kain-Fritsch cumulus parameterization (Kain 2004) (applied only for D01), the WRF single-moment 6-class scheme (Hong and Lim 2006) for cloud microphysics, the Mellor-Yamada-Janjic (MJY) scheme (Janjic 2002) for boundary-layer mixing, a revised Monin-Obukhov scheme for surface flux parameterizations (Jimenez et al. 2012), a revised version of the Rapid Radiative Transfer Model shortwave/longwave scheme (Mlawer et al. 1997) for atmospheric radiative transfer, and simple five-layer thermal diffusion scheme (Skamarock et al. 2008) for land-surface processes. The MYJ boundary-layer mixing scheme was also chosen in Oku et al. (2010), Takemi (2009), Takemi et al. (2010), and Takemi (2013) in their study on the numerical representation of airflows in complex terrain. Corresponding to this MYJ boundary-layer scheme, the Monin-Obukhov surface scheme provides better turbulent momentum and heat fluxes on various land-use surfaces under various stability conditions. The surface roughness features are represented with roughness lengths depending on the land-use types. Because the surface roughness is not explicitly represented and resolved in the WRF model, the surface winds especially in complex geometry such as densely built urban districts are not decelerated in the surface layer and are overestimated in the WRF model (Nakayama et al. 2012; Takemi et al. 2019).

One may argue that the present simulations can be conducted in a large-eddy simulation (LES) mode. The main reason why we did not use the LES framework is that the grid spacing of 200 m is considered not to be a sufficiently high-resolution to conduct LES. There are two issues here. One issue is that the 200 m grid spacing will not appropriately resolve turbulent motion within an inertial subrange. Takemi et al. (2010) showed that even with the 80 m grid spacing, the WRF model is not able to reproduce realistic wind fluctuations. The second issue is that proper generation of turbulent flows to be ingested in the main computational domain must be performed to simulate turbulent flows; without turbulent ingestions, the flows look laminar (Nakayama et al. 2012). Therefore, the present study is not intended to employ an LES model but to use a standard boundary-layer mixing scheme in order to focus on the influences of terrain representation on airflow simulations.

The initial and boundary conditions for the numerical experiments were prescribed by the Japanese 55-year reanalysis data (JRA55, Kobayashi et al. 2015) having the spatial resolution of 1.25° and the time interval of 6 h. It is critically important to closely reproduce the actual track of a typhoon for impact assessment studies, because the typhoon hazards at local scales strongly depend on the track as well as the intensity of the typhoon. For this purpose, a spectral nudging technique (Skamarock et al. 2008) was applied to the horizontal winds above the 700 hPa level in D01. In this nudging technique, the maximum wavenumber and the nudging coefficient were chosen as 2 and 0.00028 s−1, respectively. Note that the spectral nudging plays a critical role in better reproducing the tracks of tropical cyclones (Takemi et al. 2016b; Takemi 2018a; Takemi 2019). Furthermore, a typhoon bogus scheme which prescribes a Rankine-type vortex was employed to initialize a typhoon vortex at the initial time of the time integration. The maximum wind speed of 70 m s−1 was set at the radius of 75 km from the typhoon center, based on the best-track data of the RMSC Tokyo-Typhoon Center. With these approaches, we are able to favorably reproduce the track and intensity of the typhoon as explained in the next section. Although the ratio of the resolution between the JRA55 data the D01 domain seems to be large (Denis et al., 2003), the successful reproduction of such typhoon properties is considered to warrant the present experimental setting.

The time integration was conducted starting at 0000 UTC 3 September 2004 for D01, 0000 UTC 5 September for D02, and 0000 UTC 7 September for D03 and D04 until 0000 UTC 9 September. The model outputs in D03 and D04 were obtained at the 10 min interval and were used for the present analyses.

Results

Comparison with the observations

Ito et al. (2016) examined the representation of the track and central surface pressure of the simulated typhoon by comparing with the Regional Specialized Meteorological Center (RSMC) best-track data and also by referring to the analysis by Kitabatake et al. (2007) and confirmed the validity of the typhoon simulation. This study further compares the simulated results with the observations of surface winds. The wind data measured at automated weather stations (AMeDAS) of Japan Meteorological Agency (JMA) were used. The wind speeds observed at the AMeDAS stations are temporal mean values averaged for 10 min and are available at the 10 min interval. In order to evaluate the performance of the simulation over Hokkaido, we validate the results simulated in D03 with the data obtained at the AMeDAS stations included in D03. The time series of the surface winds at each station during the period from 0000 UTC 7 September to 0000 UTC 9 September 2004 are used to compute some statistics.

Figure 2a shows the difference of the simulated mean wind speed against the observed one at each location. In general, the simulation overestimates the observation but mostly reproduces the observed value within a few m s−1. There are some places in which overestimation is large. As shown shortly, the larger overestimation appears at locations where mean wind speed is higher.

Fig. 2
figure 2

Comparison of the surface winds simulated in D03 with the observations at surface stations. a The difference of the temporally averaged wind speed from the simulation to the observation. b The ratio of the variance of the simulation wind speed against that of the observation. c The Spearman’s rank correlation coefficient of the temporal change of wind speed between the simulation and the observation. The color shading in all the panels indicates the ground elevation (m)

Full size image

The representation of the wind speed variability in the simulation is examined in terms of the ratio of the variance of the simulated wind speed to that of the observation. Note that the ratio of 1 means that the observed wind variability is successfully reproduced in the simulation. The spatial distribution of this variance ratio is demonstrated in Fig. 2b. At most of the stations, the ratios are seen to be within the range of 1 to 3, meaning that the observed variability is well captured in the simulation. There are some stations indicating larger values of this ratio, in areas to the east or to the south of mountains.

The correlation of the simulated wind speed with the observation is further investigated. For this purpose, the Spearman’s rank correlation coefficient between the temporal changes of the simulation and the observation was chosen. The idea of using this correlation is based on a thought that the temporal change of wind speed strongly depends on the track and translation speed of a typhoon in addition to local effects such as topography and roughness and hence the exact correspondence of the simulation with the observation would never be expected. In contrast, the rank of the observed wind speed has a chance to be reproduced in the simulation. The spatial distribution of the Spearman’s rank correlation coefficients over Hokkaido is exhibited in Fig. 2c, indicating that the coefficients are mostly in the range of 0.6 to 1.0. Thus, the temporal change of the simulated wind speeds correlates favorably with the observed change.

Lastly, the bias of the simulated wind speed against the observation is shown. The wind bias at each station is calculated by:

$$ \mathrm{Bias}=\frac{1}{N}{\sum}_{t=1}^N\left(\sqrt{U^2+{V}^2}-{U}_{obs}\right), $$

where U and V are the simulated wind components at the 10 m height in the eastward and northward directions, respectively, Uobs is the observational wind speed, t is a time step, and N is the number of the time steps. How the bias of the simulated winds changes with the mean wind speed at each station is indicated in Fig. 3a. The bias generally increases with the mean wind speed, and the linear regression line indicates this increasing trend. The larger overestimations shown in Fig. 2a are considered to be due to larger biases at higher wind speeds. The relative bias at each station compared to the overall mean bias in D03 is demonstrated in Fig. 3b. Here, the mean bias averaged over all the stations in D03 was computed, and the difference of the bias at each station from this area-mean bias was then derived. It does not appear that there is any clear feature in the geographical characteristics of the relative bias. However, it is seen that the relative biases at inland locations (especially included in the D04 areas) are generally small, suggesting that local topography has no preferred influences on the wind representation and bias. This is a positive aspect in investigating the wind representation in complex terrain.

Fig. 3
figure 3

The characteristics of the wind speed bias. a The relationship between the wind speed bias against the simulated mean wind speed. The linear regression line is indicated by red line. b The difference of the wind speed bias at each station against the spatially averaged bias over all the stations in Hokkaido. The color shading in (b) indicates the ground elevation (m)

Full size image

Based on the characteristics demonstrated in this subsection, we will investigate the representations of winds simulated in D04 and how the horizontal grid spacing affects the wind representations in the following subsections.

Sensitivity to horizontal resolution

In order to examine the sensitivity of the wind representations to the horizontal grid spacing, we choose one area among the three D04 domains, respectively, in northern, eastern, and southern Hokkaido. The analysis areas chosen here are denoted as A, B, and C in northern, eastern, and southern Hokkaido (Fig. 4a). These areas A, B, and C include forests that were damaged by the typhoon (Takano et al. 2016; Morimoto et al. 2019), and each area A, B, and C is referred to as Nakagawa, Abashiri, and Niseko (see Fig. 1 of Takano et al. 2016). The differences in the terrain representation between the 1 km grid and the 200 m grid case are compared in Fig. 4b, c, d for the area A, in Fig. 4e, f, g for the area B, and in Fig. 4h, i, j for the area C. As clearly seen, the detailed features of the terrain height are more realistically reproduced in the 200 m grid case than in the 1 km grid case. Peaks, ridges, and valleys are represented in a more enhanced way in the higher-resolution case. With these terrain characteristics in mind, the wind representations will be investigated in the following.

Fig. 4
figure 4

Terrain features of the analysis areas. a Terrain height and the D04 areas in D03. Red boxes denoted as A, B, and C are the analysis areas. The terrain height in the area A (b) at the 1-km grid spacing in D03 and (c) at the 200 m grid spacing in D04. d The frequency distribution of terrain height in the area A (left) at the 1 km grid and (right) at the 200 m grid. e, f, g The same as (b, c, d) except for the area B. h, i, j The same as (b, c, d) except for the area C

Full size image

Figure 5 compares the representation of the maximum winds between the 1 km grid and the 200 m grid case in terms of the frequency distribution of the maximum surface wind speed during the simulated time period (i.e., from 0000 UTC 7 September to 0000 UTC 9 September 2004) at each grid point in each analysis area. In the 200 m grid case, the higher (lower) extremes become stronger (weaker) and more frequent, and the distribution is more widely spread. Because the detailed terrain features are more enhanced in the higher-resolution case, the contrasts between higher and lower elevations and between sharper and flatter slopes affect the magnitude of the maximum winds. Takemi (2018b) demonstrated that the difference in terrain representation between higher-resolution and lower-resolution topography data is very significant even at the same grid spacing, leading to major differences in quantitative rainfall simulations. This study reveals such influences of the terrain representation on the maximum winds.

Fig. 5
figure 5

The frequency distribution of the maximum surface wind speeds in the analysis areas. The relative frequency distribution (in %) of the maximum wind speeds during the simulated time period at the grid points in the area (a) A, (b) B, and (c) C for (left) the 1 km grid case and (right) the 200 m grid case

Full size image

Representations of surface winds in the areas A, B, and C in the simulations with the 1 km grid and the 200 m grid are investigated by referring to the slope angle of the model terrain. Figure 6 shows how the pointwise temporal mean wind speed changes with the terrain slope angle. Because the number of grid points is different between the 1 km grid and the 200 m grid case, the number of points in Fig. 6 is also different between these different grid cases. In all the analysis areas, the range of the scatter of the mean wind speed is wider in the 200 m grid case than in the 1 km grid case. The slope angle is also spread more widely between 0 and 30° in the higher-resolution case than in the lower-resolution case, in agreement with the widely distributed terrain height in the higher-resolution case (see Fig. 5). The linear regression line has a positive slope, which means that in general the mean wind speed increases with the increase in the slope angle. However, the rate of the increase of the mean wind speed is largely different between the two resolution cases and is gentler in the higher-resolution case. In other words, at a certain slope angle of the terrain, the mean winds can be both stronger and weaker.

Fig. 6
figure 6

Relationships between terrain slope and mean wind speed in the analysis areas. Relationships of the temporal mean wind speeds against the terrain slope angle for the 1 km grid case (black dots) and the 200 m grid case (gray dots) in the area (a) A, (b) B, and (c) C. Red solid and dotted lines indicate the linear regression line for the 200 m grid and the 1 km grid case, respectively

Full size image

The relationship between the maximum wind speed during the simulated time period and the terrain slope angle is shown in Fig. 7. The features identified for the mean wind speed indicated in Fig. 6 are similarly seen for the maximum wind speed. The maximum winds as well as the terrain slope angles are widely spread in all the areas. The increasing tendency of the maximum wind with the increase in the slope angle is also seen, but the difference in the regression line’s slope between the 200 m grid and the 1 km grid case appears to be smaller for the maximum wind speed than for the mean wind speed.

Fig. 7
figure 7

Relationships between terrain slope and maximum wind speed in the analysis areas. The same as Fig. 6, except for the maximum wind speed during the simulated time period

Full size image

Table 1 summarizes spatially evaluated statistics computed for the pointwise mean and maximum wind speeds within each analysis area. The area-mean values for both the temporal mean and maximum winds in all the analysis areas are not so different between the 1 km grid and 200 m grid cases. In contrast, the area-maximum values are 1.6 to 2.1 times larger in the 200 m grid case than in the 1 km grid case. The spatial variances are also larger in the 200 m grid than in the 1 km grid case; for the temporal maximum wind speed, the variances are 2.6 to 5.7 times larger in the 200 m grid than in the 1 km grid case. These statistics are the reflection of the scatter characteristics indicated in Fig. 6.

Table 1 Statistics of the mean and maximum wind speeds in the analysis areas. The results of mean and maximum wind speeds during the simulated time period obtained in the 1 km grid and 200 m grid cases for the areas A, B, and C are summarized. The statistics are the area-mean, area-maximum, and spatial variance values computed from the temporal mean and maximum winds in each analysis area
Full size table

From the analysis shown in Figs. 5, 6, and 7, the characteristics of the temporal and spatial variability of winds appear differently in the 1 km and 200 m grid cases. In the next subsection, we will further show how the representation of terrain features in the higher and lower-resolution cases affect the variability of surface winds.

Relationship with the variability of terrain height

The representation of the model terrain at a higher resolution compared with a lower resolution is investigated in terms of the subgrid-scale orography used in Jimenez and Dudhia (2012). Subgrid-scale orography can be represented as how terrain heights vary spatially within a certain area and is here defined as the standard deviation of terrain heights (σSSO) at the 200 m grid spacing within a 1 km by 1 km area. In this study, we evaluate the subgrid-scale orography as the area-mean standard deviation of terrain heights over each D04 area as follows:

$$ \mathrm{Subgrid}-\mathrm{scale}\ \mathrm{orography}=\sum \limits_{i=1}^{N_{1\ \mathrm{km}}}{\sigma}_{SSO,i}=\sum \limits_{i=1}^{N_{1\ \mathrm{km}}}{\left\{\sqrt{\frac{1}{n_{200\ \mathrm{m}}}\sum \limits_{j=1}^{n_{200\ \mathrm{m}}}\left({h}_j-\overline{h}\right)}\right\}}_i, $$

where hj and \( \overline{h} \) are, respectively, the terrain height at the 200 m grid and its mean at a grid point j within a 1 km by 1 km area, n200 m is the number of 200 m resolution grid points within a 1 km by 1 km area (i.e., 25), N1 km is the number of 1 km resolution grid points within a D04 area, and i is a certain grid point within the D04 area.

As an example, Fig. 8 indicates the spatial distribution of the standard deviation of the 200 m grid terrain within the area B (see Fig. 4f). The area-mean standard deviation for this area becomes 60.7. In this way, the area-mean standard deviation of terrain heights is computed for the 8 D04 areas. In the followings, the relationships between the subgrid-scale orography and the wind representations are indicated.

Fig. 8.
figure 8

The standard deviation of subgrid-scale orography. The spatial distribution of the standard deviation of the subgrid-scale orography computed in a 1 km by 1 km area from the model terrain at the 200 m grid spacing, for the area B (see Fig. 4f) as an example

Full size image

Figure 9 demonstrates how the area-means, area-maxima, and spatial variances of the temporal mean wind speeds within the D04 areas vary against the subgrid-scale orography that is defined as the area-mean standard deviation of the terrain heights at the 200 m grid spacing. In Fig. 9a, the area-mean value of the temporal mean wind speeds at the grid points in each D04 area is shown. It is seen that the area-mean values from the 200 m grid and the 1 km grid case are very similar with each other but slightly higher in the 1 km grid case than in the 200 m grid case. The reason why the mean wind speeds are stronger in the 1 km grid case than in the 200 m grid case is due to the difference in roughness. That is, the more complex feature of the terrain represented in the 200 m grid case results in larger roughness on average in a certain area, compared to the terrain represented in the 1 km grid case. It is also found that the mean winds gradually increase with the increase in the subgrid-scale orography.

Fig. 9
figure 9

The relationships between the subgrid-scale orography and the wind representations. The relationships of the area-mean standard deviation of terrain heights within each D04 area with a the area-mean wind speed, b the area-maximum wind speed, and c the spatial variance of wind speed in the same D04 area. Black and gray points indicate the results of the 200 m grid case and the 1 km grid case, respectively

Full size image

Figure 9b examines the relationship of the area-maxima among the time-mean wind speeds in each D04 area against the subgrid-scale orography. Overall, the area-maxima are stronger in the 200 m grid case than in the 1 km grid case. Furthermore, in the 200 m grid case, the area-maximum wind becomes stronger with a larger subgrid-scale orography, while in the 1 km grid case, the area-maximum does not change much with the change in the subgrid-scale orography. These features can also be found for the spatial variances of the time-mean wind speeds in the D04 areas, as indicated in Fig. 9c.

In this way, it was found that the area-maxima and spatial variances of the time-mean surface winds appear more pronouncedly in the higher-resolution case than in the lower-resolution case. We have also examined the features of the area-means, area-maxima, and spatial variances for the temporal maximum wind speeds in the D04 area and found that the area-mean, area-maxima, and variances all increase with the increase in the subgrid-scale orography (not shown). In other words, the complexity of the terrain features enhances the mean and maximum wind speeds over the simulated time period as the forms of the statistics examined above.

Discussion

We have investigated the representations of typhoon-induced surface winds in complex terrain in high-resolution simulations of a specific event in Hokkaido, the northern part of Japan, as a case study for Typhoon Songda (2004). The comparison of the simulated winds at the 1 km grid and the 200 m grid indicates that the range of the simulated wind speeds becomes wider and that the higher and lower extrema become more extreme in the higher-resolution case. We have also identified that the mean and maximum wind speeds during the simulated time period increase with the increase in the terrain slope angle and that the scatter becomes wider with a higher slope angle in the higher-resolution case. In other words, the response of the simulated winds to the slope angle appears differently between the 1 km grid and 200 m grid cases. By using a quantity representing the subgrid-scale orography, we have shown that the area-maxima and spatial variances over the D04 areas become more enhanced with the increase in the subgrid-scale orography in the 200 m grid case than in the 1 km grid case.

From these analyses, it is suggest that the temporal and spatial variability of surface winds becomes larger in the higher-resolution case, which is closely related to the more complex features of the terrain reproduced in the model. An implication of this result is that high-resolution simulations at the grid spacing of the order of 100 m are able to highlight fluctuating natures of surface winds and can be regarded as one of the benefits or added values of the high-resolution simulations of airflows in complex terrain. Highly fluctuating winds lead to enhanced extremes both in the higher and lower ranges, and hence more hazardous winds will be represented in the higher-resolution case. Therefore, this point is critically important in assessing the strong wind hazards in complex terrain.

Because the current modeling framework does not take into account explicitly the effects of turbulence through an LES approach, the simulations are not expected to reproduce turbulent flows. In spite of this deficiency, it is emphasized that the present meteorological modeling with the use of high-resolution DEM data can reproduce enhanced fluctuating natures of surface winds in time and space, which are generated by fine-scale features of terrains. Namely, high-resolution simulations of airflows in complex terrain can benefit from detailed representations of topographic features of terrain in representing extreme winds.

Note here that the present simulations employ the MYJ boundary-layer scheme, i.e., a type of Reynolds averaging Navier-Stokes model, and hence the simulated winds by the WRF model should be regarded as means averaged for a certain period of time. Nakayama et al. (2012) compared the surface winds among the observations, the WRF simulations at the 60 m grid spacing, and the LES results and found that the WRF winds correspond closely to the 10 min means. Although the temporal averaging scale of the WRF winds may depend on the surface roughness features and the meteorological conditions, the WRF winds are regarded as the temporal means for around 10 min. Even if we make WRF outputs at a higher frequency (e.g., 1 min), the WRF winds does not fluctuate in time. Therefore, we do not evaluate turbulent fluctuations in this study. What we are evaluating here is regarded as the temporal variation of winds induced by undulations of topography, corresponding to the time-scale of 10 min or longer.

Previous studies on airflows in complex terrain (Jimenez et al. 2010; Jimenez and Dudhia 2012; El-Samra et al. 2018) conducted numerical simulations with the horizontal grid spacings of 1 km or 2 km. Although such grid spacings appear to be insufficient to resolve actual terrains, the geometrical features of the terrains in their studies may be characterized with larger-scale undulations, as compared with the terrain features of the present analysis areas in Hokkaido, Japan. Bonnardot and Cautenet (2009) demonstrated that the simulations at 1 km and 200 m resolutions better reproduced local circulations than those at 5 km resolution because of a better representation of the local terrain. By distinguishing meteorological situations, they further showed that under a weak synoptic flow condition, local circulations prevailed and the high-resolution simulation at the 200 m grid better reproduced meteorological characteristics. Thus, it is considered that depending on the background meteorological conditions, the response of the meteorological changes to terrain features appears differently.

In this study, we examined a strong wind case induced by a severe typhoon. It was found that gust factors are higher in tropical cyclone cases than in extratropical cyclone cases (Paulsen and Schroeder 2005). Therefore, wind fluctuations including gust factors occur more clearly in typhoon cases, and hence reproducing highly fluctuating natures of typhoon winds in high-resolution simulations is an advantage in investigating the physical processes of strong winds in complex terrain and moreover assessing the impact and hazard of strong winds at local scales.

Conclusions

This study investigated the representation of surface winds in complex terrain in numerical simulations with the horizontal grid spacing of 200 m. The strong wind event in Hokkaido during the passage of Typhoon Songda (2004) was chosen as a case study, because the strong winds by this typhoon caused widespread damages to forests in Hokkaido. We used the WRF model to conduct regional-scale simulations at the 200 m grid by employing the same computational settings with those in Ito et al. (2016) except setting the 200 m domains. The results of the 200 m grid simulations were examined and evaluated in comparison to those of the 1 km grid case. After validating the simulated results with the surface observations, we explored the benefits of the high-resolution simulations of surface winds in complex terrain.

The 200 m grid simulations were able to clearly represent more enhanced and more frequent extremes both in the stronger and weaker ranges of surface winds than the 1 km grid simulations. This result corresponded to the difference in the terrain representations at these two grid spacings. It was also shown that the mean and maximum winds during the simulated time period increase with the increase in the slope angle of the terrain both in the 200 m grid and 1 km grid cases. Because the scatter of the mean and maximum winds is much wider in the 200 m grid case than in the 1 km case, the slopes of the linear regression lines are smaller in the 200 m case. The characteristics of the temporal and spatial variability of surface winds appeared differently in the 200 m and 1 km grid cases.

To assess the influences of the terrain complexity on wind variabilities, we introduced subgrid-scale orography which was defined as the area-mean value of the standard deviations of the terrain heights at 200 m grid within a 1 km by 1 km area. With this quantity, it was demonstrated that the area-maxima and spatial variance of surface winds over the 200 m grid computational domains are more enhanced with the increase in the subgrid-scale orography in the higher-resolution case than in the lower-resolution case.

Based on the present results, it is suggested that the high-resolution simulations at the 200 m grid are able to highlight and enhance the fluctuations of surface winds in complex terrain, because the complex topographic features of the terrain is better represented in the 200 m grid than the 1 km grid case. This is considered to be the benefits or added values of the high-resolution simulations of airflows in complex terrain. It is also considered that the topographic features of Hokkaido examined in this study are highly complex and thus the grid spacing on the order of 100 m is desirable.

Because the mountainous regions, characterized by abundant forest resources as well as complex topographic features, prevail over the Japanese islands, the temporal and spatial variability of surface winds largely affect the forest ecosystem. A strong wind event can be a great threat to the regions. Especially in the northern part of Japan, the landfall or approach of typhoons is rare, and hence the impacts of typhoons once occurred are huge. We have investigated the case of Typhoon Songda (2004), but there were a number of typhoons that affected northern Japan. For example, Typhoon Mireille (1991) landed on Hokkaido and caused severe damages on agriculture and forestry. Actually, this specific typhoon spawned the most costly insurance loss among the tropical cyclones in the Pacific region during the period from 1970 to 2017 (Takemi et al. 2016a). From viewpoints of disaster prevention and mitigation, it is critically important to quantify the impacts of typhoons in complex terrain where trees and forest are vulnerable to strong winds and also assess the impacts of climate change on such typhoon hazards. Climate change impacts on strong winds and heavy rainfalls in Hokkaido during typhoon passages were investigated by Ito et al. (2016), Takemi et al. (2016a), Kanada et al. (2017), and Nayak and Takemi (2019a, 2019b). Such typhoon hazard information, once obtained with the resolution on the order of 100 m, needs to be fully used to impact assessment studies on forest ecosystems.