Numerical Simulation of Volume Change of the Backshore Induced by Cross-Shore Aeolian Sediment Transport

: Predicting the morphological changes of the backshore is vital for appropriate beach management because the backshore plays a signiﬁcant role in the ecosystem and disaster prevention. In this study, a one-dimensional model was developed and applied to the Hasaki Coast in Japan to predict changes in backshore volume. The volume change was estimated from the di ﬀ erence between the aeolian sediment transport rates at the seaward and landward boundaries of the investigation area, considering the wind velocity and direction, sediment size, precipitation, and vegetation in the process. The model was calibrated and validated using the ﬁrst and second halves of beach proﬁle data obtained weekly at the Hasaki Coast over a 28-year period from 1987 to 2014. The validation suggests that the model can reasonably reproduce the cumulative volume change, which is the amount of volume change from the initial value, but it underestimates the time-varying ﬂuctuations of the weekly averaged volume-change rate. This can be attributed to the presence of small-scale features, such as dense vegetation and wrack, which are not taken into account in the model. Although the model performance for the cumulative volume change was good, it overpredicted the values in the second half of the validation process. This can be attributed to the fact that the model is not able to predict reductions in the aeolian sediment transport rate caused by an increase in beach steepness.


Introduction
The backshore, which is the area from the berm crest to the foredune foot [1], and the foredune both play an important role in protecting infrastructure and properties against flooding and wave overtopping. The backshore is an ecologically valuable region where materials are exchanged between the beach/surf zone and sand-dune zone, such as sediment, salt spray, groundwater, and organic materials [2,3]. Hence, predicting the morphological changes in the backshore and foredune is crucial for proper beach management.
Backshores and foredunes are developed by landward aeolian sediment transport; conversely, they are eroded by seaward sediment transport due to run-up waves during storms. To predict such morphological changes, numerical models have been developed to investigate dune and backshore evolution. For instance, van Dijk et al. [4], Luna et al. [5], and Duran and Moore [6] numerically investigated the influence of vegetation, sediment supply, and adaptation length (i.e., the distance to reach the equilibrium state) on dune development due to landward winds. The dune evolution due to winds and waves was examined by Davidson-Arnott et al. [7] using a model that assumes temporally constant amounts of sediment supply and wave-induced erosion. Larson et al. [8] developed a model incorporating the sediment exchanges among the dune, backshore, foreshore, and longshore bar, which was then applied to field sites by Palalane et al. [9] and Hallin et al. [10]. Cohn et al. [11] developed a process-based model to predict dune evolution, Windsurf, by coupling three models for wave-induced sediment transport and morphological change (XBeach) [12], aeolian sediment transport (Aeolis) [13], and dune development [6]. Roelvink and Costas [14] developed another process-based model by combining XBeach and Duna, the latter of which is a model for predicting wind-induced sediment transport and topographical change. These process-based models have also been applied to field sites.
The quantitative validations of the models developed by Palalane et al. [9], Hallin et al. [10], and Cohn et al. [11] have been confirmed using field data. However, the duration periods of the field data used for the validations were relatively short: approximately four years for Palalane et al. [9], approximately six years for Hallin et al. [10], and one year for Cohn et al. [11]. Although Palalane et al. [9] calculated morphological changes at another beach for 18 years and compared the calculation results with the field measurements, field data were sparse and consisted of just four measurements. Furthermore, in validating the models, the influence of waves and wind were combined, which is to say that the two factors were not individually tested, which, in turn, resulted in ambiguity with respect to model validation.
In this study, we focused on morphological change on the backshore induced by aeolian sediment transport. This study aimed to develop a model to simulate the backshore volume change induced by aeolian sediment transport and to examine its accuracy using field data obtained over a 28-year period.
Dune growth and the aeolian sediment transport rate are influenced by a number of factors, such as wind velocity and direction, sediment size, fetch and moisture [15][16][17][18][19][20][21][22][23], sediment size distribution [24], beach slope [25][26][27], vegetation [21,28,29], and wrack [30,31]. Accordingly, those factors were included in the abovementioned models for morphological change induced by aeolian sediment transport. In this study, we developed a one-dimensional model for volume changes of the backshore considering the influences of wind velocity and direction, sediment size, precipitation, and vegetation. The details of the model are explained in Section 3. The field data used for model calibration and validation consisted of beach profile data obtained at the Hasaki Coast in Japan over 28 years. The beach has a relatively wide backshore; only the area affected by wind was used for calibration and validation. The data used in this study are described in Section 2.

Study Site
The study site was the Hasaki Coast located in eastern Japan facing the Pacific Ocean (Figures 1 and 2). Measurements were taken at the Hazaki Oceanographical Research Station (HORS), which has a 427-m long field-observation pier. Offshore waves are large during January-April and September-October due to tropical cyclones [32]. The bathymetry around the HORS is almost uniform alongshore, and the influence of pier pilings on the bathymetry is relatively minimal [33]. The sediments in the study site are mostly sands deposited during the Holocene [34] with a thickness of approximately 20 m [35]. The medium sediment size is approximately 0.2 mm [36]. The beach slope was approximately 1/25 near the shoreline and 1/120 at the tip of the HORS pier, where the water depth was approximately 5 m [37]. The high, mean, and low water levels were 1.25, 0.65, and −0.20 m, respectively.

Beach Profile
Along the pier, the profiles were measured at 5-m intervals using a rope with graduated depth markings and a 5-kg lead, from the pier. On the land side of the pier, the bed profile was directly measured by a level and a staff. All these bed profiles were obtained every workday from 12 March 1986 to 31 March 2011. Beginning on 5 April 2011, measurements were taken once a week, mostly on Mondays.
The beach profile data used in this study consisted of beach profiles measured once a week from 5 January 1987 to 3 November 2014. When readings were not taken, data were inferred by linearly interpolating the adjacent measured values.
The area investigated in this study extended from x = −115 to −65 m (Figures 2 and 3), as in Kuriyama et al. [29]; x is the seaward distance relative to the reference point near the pier entrance. Kuriyama and Yanagishima [37] applied rotated empirical orthogonal function analysis to the beach profile dataset and claimed that the backshore was located within x = −115 to −50 m. However, to eliminate the influence of waves, we set the seaward limit of the investigation area to x = −65 m, which is the seaward limit of the vegetation area and was reached by run-up waves twice during a one-year period from July 1995 to June 1996 [29].  Figure 4 shows the temporal variation of cumulative volume change from the initial profile recorded on 5 January 1987. As seen in Figure 4, the volume of the backshore increases almost monotonously with time, with some fluctuations occurring in relatively short time periods.

Winds
Wind velocities and directions were measured at the tip of the pier at an elevation of 10 m relative to the DL using a propeller anemometer for 10 min every 3 hr. The wind direction was defined as being relative to the seaward shore-normal direction, which is positive in the counterclockwise direction ( Figure 5). Accordingly, wind directions of 0 and 90 degrees indicate wind blowing landward in the shore-normal direction and blowing southward in the alongshore direction, respectively. Missing data were corrected as described in Appendix A.

Rainfall
Rainfall was measured at the HORS, but the measurement period was limited from 1995 to 1999. Thus, for the simulation, we used rainfall data obtained every three hours at the Choshi Weather Observation Station of the Japan Meteorological Agency, which is located 14 km south of the HORS. The root mean square (RMS) error between the rainfall values at the HORS and Choshi during the four-year period was 0.06 mm/hr. The coincidence ratio of rainfall at the HORS and Choshi was defined as the ratio of times in which measured precipitations at both stations were equal to, or lower than, the threshold value, or both measured values were greater than the same threshold value. When the threshold was set at 0 mm/hr, the coincidence ratio was 98.1%. When the threshold was set at 1 mm/hr, which was the critical value used in the simulation, the coincidence ratio was 98.5%. There is another weather observation station at Kashima, which is located 18 km north of the HORS. The RMS error here was 0.11 mm/hr; the coincidence ratios of precipitation with threshold values of 0 and 1 mm/hr were 95.3% and 96.4%, respectively.
The monthly averaged precipitation values show that the amount of rainfall was relatively large in June due to seasonal rain fronts as well as from September to October due to typhoons ( Figure 8). The yearly averaged values varied from 1000 to 2000 mm/year; the mean annual precipitation was 1473 mm/year ( Figure 9).

Vegetation
The vegetation data were obtained near the investigation area from 1995 to 1998 [29] and in the investigation area from 1995 to 2010. Approximately every month from 1995 to 1998, the species, height, and vegetation-coverage rate in square meshes (1 × 1 m) were recorded every 10 m from x = −110 to −70 m along three transects located 100, 120, and 140 m north of the study survey line. Figure 10 shows examples of vegetation in a square mesh, which was observed on 24 June 1996 and 14 January 1997 at x = −100 m along the transect 120 m north of the survey line. The obtained data indicate that the study site is occupied by perennial creeping beach grasses, including Carex kobomugi (Japanese sedge, Asiatic sand sedge), Calystegia soldanella (sea bindweed), and Elymus mollis (American dunegrass). The monthly averaged vegetation-coverage rate, C veg1 , suggests that beach grass grows from May to September, as shown in Figure 10a, but disappears from December to March, as shown in Figures 10b and 11. The vegetation-coverage rate represents seasonal variations in the ratio of the projected area of stems and leaves within a vegetation patch.  From 1995 to 2010, the vegetation along and around the survey line was visually observed using photos taken from a nearby building from June to August in each year. On the north side of the survey line, beach grass was removed for another field experiment in 1994. The time series of the ratios of the area covered by vegetation patches on the north and south sides of the survey line, which are denoted by C veg2,N and C veg2,S , respectively, show that it took approximately 13 years for the beach grass to fully recover ( Figure 12). Figure 13 shows the photos taken on 2 August 1995, when vegetation on the north side of the survey line was stripped, and on 1 June 2007, when vegetation recovered.

Simulation Model of Backshore Volume Change
The temporal change rate of the sediment volume in the investigation area per unit width, V, was determined as a balance of volume flux at the landside and seaside boundaries, which can be expressed as: where F S and F L denote the volume flux of aeolian sediment transport per unit width at the seaside and landside boundaries of the investigation area, respectively. These volume fluxes are positive in the landward direction. In a discretized form, the volume change in the investigation area per unit width at step time of i-th, V i , can be expressed as: where ∆t denotes the time interval, which, in this study, was set to 3 hr. The volume flux at both boundaries can be expressed as: where ρ s denotes the sediment density (2.65 × 10 3 kg/m 3 ); λ is the sediment porosity (0.3); q is the aeolian sediment transport rate on a bare sandy bed, which is positive in the landward direction; R b is the ratio of bare area and subscripts 1 and 2 denote the values in the investigated area and in the landside area, respectively; and c 1 and c 2 are the calibration coefficients. As seen in Equations (3) and (4), the aeolian sediment transport rate is multiplied by the ratio of bare area on the upwind side of the sediment transport to represent the volume fluxes, F L and F S . Note that the ratio of bare area on the seaside of the investigation area is assumed to be in unity when the aeolian sediment transport rate is positive, as seen in Equation (4). These bare-area ratios can be expressed as: where α w is the wind direction determined in Figure 3. The parameters C veg1 and C veg2 represent the monthly averaged ratio of the projected area of stems and leaves within a vegetation patch and the yearly averaged ratio of the area covered by vegetation patches, respectively. The coefficients, c 1 and c 2 , were only introduced for volume flux at the landside boundary. This is because the local wind conditions at the landside boundary (i.e., at the dune foot) should differ from those measured at the pier tip. Moreover, this study applied different c 1 and c 2 values because the dune-foot wind characteristics can differ depending on wind direction.
Following Kuriyama and Mochizuki [38], Kuriyama and Kamidozono [39], and Kuriyama et al. [40], the aeolian sediment transport rate was estimated using Kawamura's formula [41], which is described as: where K is a nondimensional coefficient, ρ a is the air density, g is gravity, u * is the wind-shear velocity, u * c is the threshold-shear velocity, u is the wind velocity at a height of z from the surface, d is the sand particle size (0.2 mm), and A is a nondimensional coefficient (0.1). The u * value was estimated considering the focal point (u', z'); it was assumed that the vertical distribution of wind velocity was logarithmic. In this study, since u z was determined from the measured velocity 10 m above the Hasaki DL, z was assumed to be 10 m. The aeolian sediment transport rate was assumed to be proportional to the square of R b because, under the equilibrium condition, the aeolian sediment transport rate can be expressed by the product of the sediment pickup rate and the saltation length, both of which are proportional to R b [39,40].
In 1994, the values of C veg2,N and C veg2,S were assumed to be the same as the values in 1995, which were 0.0 and 1.0, respectively. Both were set as 1.0 from 1987 to 1993 and from 2001 to 2014.
Twenty-four hours after precipitation ≥1 mm/hr, the aeolian sediment transport rate was assumed to be zero, following Kuriyama et al. [40], who employed the said assumption based on the field measurements of the water content of sediment on the backshore surface after precipitation.

Calibration and Validation
The simulation model of volume change was calibrated using the field data obtained from 1987 to 2000, which is approximately half of the investigation period; the simulations were validated by comparing the results with the remaining data from the other half of the investigation period, that is, from 2001 to 2014.
In the calibration process, the coefficients, K, c 1 , and c 2 , were determined using the least-squares method, thereby minimizing the error between the measured and simulated weekly averaged volume-change rates. In the calibration and validation processes, the RMS errors and the correlation coefficients between the measured and simulated values of the weekly averaged volume-change rate and the cumulative volume change were estimated, as well as their Brier skill score (BSS) [42], which can be described as: where B s and B m denote the simulated and measured values, respectively. The angle brackets indicate the temporary averaged values. According to van Rijn et al. [43], when the BSS value of a model for morphological change is <0, 0-0.3, 0.3-0.6, 0.6-0.8, and 0.8-1.0, the model performance is categorized as bad, poor, reasonable, good, and excellent, respectively.

Calibration
The determined K, c 1 , and c 2 values were 1.00, 0.58, and 0.53, respectively. Although the model underestimated variation in the weekly averaged volume-change rate (Table 1 and Figure 14), it reasonably reproduced the cumulative volume change (Table 1 and Figure 15). The BSS indicates that model performance is poor for estimating relatively short-term volume-change rates according to the criterion of van Rijn et al. (2003); however, it is excellent for the cumulative volume change.

Validation
Similar to the results of the calibration process, in the validation process, the model had relatively poor prediction skills with respect to the weekly averaged volume-change rate (Table 2 and Figure 16), but it could reasonably reproduce the cumulative volume change (Table 2 and Figure 17). The model performance is bad for the volume-change rate, but good for the cumulative volume change. Table 2. The root-mean-square error, the correlation coefficient, and BSS for the weekly averaged volume change rate and the cumulative volume change in the validation process.

Discussion
The coefficient K in Equation (4) was estimated to be 1.00. Although this value is smaller than the K values of 2.3-3.1 obtained in the field using trench traps [41], it is equal to the K value of 1.0 in Kuriyama et al. [40] and close to that of 0.8 in Kuriyama and Mochizuki [38] and Kuriyama and Kamidozono [39]; K values of 1.0 and 0.8 [38][39][40] were obtained approximal to the study site by comparing the measured and simulated aeolian sediment transport rates. Hence, a K value of 1.00 is reasonable.
The obtained c 1 and c 2 values were 0.58 and 0.53, respectively. Based on the calculation results of aeolian sediment transport rate using reduced wind velocities, the obtained c 1 and c 2 values were equivalent to a wind-velocity reduction of 89% and 88% at the landward limit of the investigation area (x = −115 m) for landward and seaward winds, respectively. According to field measurements and numerical simulations, at the downwind side of the foredune, the wind velocity decreases and reverses close to the dune crest [44][45][46][47]. Even at the upwind side of the foredune, the wind velocity decreases near the dune toe [48,49] with sediment deposited in front of the dune toe [50]. Accordingly, calibrated values of c 1 and c 2 are consistent with the measured values with respect to aeolian sediment transport in the vicinity of the dune toe.
Sensitivity tests for K, c 1 , and c 2 , which were conducted by increasing and decreasing each parameter by 10%, suggest that the model is most sensitive to K and c 1 , and least sensitive to c 2 ( Figure 18). By changing the parameters by 10%, the simulated cumulative volume change at the end of the investigation period varied by approximately 10% for K and c 1 , and 3% for c 2 . The difference between the ratios of c 1 and c 2 indicates that, in the investigation area, landward winds contributed more to volume change than seaward winds.
Although the model developed in this study is not able to reproduce the weekly averaged volume-change rate of the backshore, it can reproduce the cumulative volume change. This can be attributed to the fact that the model reproduces the low-frequency components of the volume-change rate, as shown in the spectral analysis results (Figures 19 and 20). The coherence between the measured and simulated volume-change rates is high at frequencies lower than 0.00112 cycle/day (periods longer than 893 days) in the calibration process and, in the validation process, at frequencies lower than 0.00209 cycle/day (periods longer than 478 days); the coherence fluctuates at frequencies higher than the abovementioned frequencies.   The reproductivity of the model is consistent with the results obtained by de Vries et al. [27] and Strysteen et al. [51], who both demonstrated that the yearly averaged volume-change rate on dunes has minimal correlation with the yearly averaged potential aeolian sediment transport rate, which was estimated from wind velocity and direction as well as sediment size; however, only Strysteen et al. [51] demonstrated that the 17-year averaged volume-change rate on dunes has correlations with the 17-year averaged potential aeolian sediment transport rate.
According to de Vries et al. [27] and Strysteen et al. [51], the causes behind the discrepancy between the yearly averaged values are limiting factors, which include sediment size distribution, moisture, and beach geometry, for the aeolian sediment transport rate. We took into account the influences of precipitation and vegetation in our model.
In terms of fetch, this study assumed that aeolian sediment transport reached the equilibrium condition. This assumption is partially supported by the results of Kuriyama and Kamidozono [39] and Kuriyama et al. [40], who numerically calculated the cross-shore variations of the aeolian sediment transport rate under equilibrium and nonequilibrium conditions at Hasaki by assuming that the shoreline was located at x = −25 m. In these papers, the aeolian sediment transport rate under the nonequilibrium condition was estimated based on the difference between the amount of sediment picked up and the amount of sediment deposited; moreover, the difference between the aeolian sediment transport rate under equilibrium and nonequilibrium conditions from x = −65 to x = −115 m, which is the same as the investigation area in this study, is quite small. Although the reported fetch distances have a wide range of 10-200 m [23], the obtained results are consistent with those of Hotta [41], who demonstrated that the minimum fetch distance is 10 m. Thus, the equilibrium condition assumption was thought to be appropriate at Hasaki when the fetch length was larger than 40 m. However, as mentioned in Section 2.2, run-up waves reached the seaward boundary of the investigation area once or twice a year, reducing the fetch length to values smaller than 40 m. This was ignored in developing the model, which may have resulted in errors in the simulation.
The discrepancy between the measured and simulated volume-change rates at high frequencies can also be attributed to relatively small-scale features, such as dense vegetation and wrack, which were not considered in the model. For instance, dense vegetation and wrack can result in changes in the alongshore continuity, which act as boundaries that reduce the aeolian sediment transport rate [30,52]. Furthermore, wrack induces scour at the downwind side, which, in turn, increases the aeolian sediment transport rate [31].
The difference between the measured and simulated volume-change rates at high frequencies can also be explained by the field measurements along one survey line. Even though the morphology on the Hasaki beach is relatively uniform alongshore, small-scale features may have influenced the aeolian sediment transport rate, as mentioned above. Thus, to eliminate the influences of small-scale features, field measurements on relatively large areas extended alongshore may be required. Kuriyama et al. [40] measured beach profiles on the backshore at Hasaki almost every month for approximately four years along three survey lines, which were located 100, 120, and 140 m north of the survey line of this study. They compared the simulated and measured monthly averaged aeolian sediment transport rates at x = −67.5 m, which were estimated based on the beach profiles averaged for three survey lines. The simulated aeolian sediment transport rates did not correlate with the measured ones. This result suggests that the average of three survey lines is not sufficient to remove the influence of small-scale features.
Although the model performance for the cumulative volume change was good in the validation process, as mentioned in Section 4.2, the model overpredicted the values during the second half of the validation period. This can be attributed to the fact that the beach became steeper during the second half of the validation period, which decreased the aeolian sediment transport rate. Indeed, our model did not include the influence of the said change. In the investigation area, the absolute value of the beach slope increased from 1/28 at the beginning of the investigation period to 1/16 at the end of the investigation period ( Figure 3). The aeolian sediment transport rate is smaller for a steep-slope beach than for a mild-slope beach [25,26]. Following Hardisty and Whitehouse [25] and Roelvink and Costas [14], including the beach-slope effect should improve the prediction accuracy with respect to cumulative volume change.

Conclusions
In this study, a one-dimensional model was developed to predict backshore volume changes. The model estimated the volume change by taking the difference between the aeolian sediment transport rates at the seaward and landward boundaries of the investigation area. The aeolian sediment transport rates were calculated based on Kawamura's formula, which takes into account wind velocity and direction as well as sediment size. The influences of precipitation and vegetation on the aeolian sediment transport rate were also considered in the model.
The three free parameters, K, c 1 , and c 2 , were determined using beach profile data obtained weekly at Hasaki, Japan, from 1987 to 2000, thereby minimizing the error between the measured and simulated weekly averaged volume-change rates. The model validation, which was implemented using the beach profile data obtained from 2001 to 2014, demonstrated that the model cannot reproduce weekly averaged volume-change rates; however, it can reproduce the cumulative volume change. Based on the BSS in the validation process, the model performance was categorized as bad for the volume-change rate but good for the cumulative volume change.
This discrepancy can be attributed to the fact that the model does not take into account small-scale features, such as dense vegetation and wrack, which may have affected the aeolian sediment transport rate at high frequencies. Although the model performance for the cumulative volume change was good in the validation process, the model overpredicted values in the second half. This can be attributed to the fact that the model did not take into account the influence of beach slope: aeolian sediment transport rate is smaller for a steep-slope beach than for a mild-slope beach. Acknowledgments: The authors would like to thank two anonymous reviewers for their useful and constructive comments to improve the original manuscript. We are grateful to Hitoshi Tamura, the Port and Airport Research Institute, for helping with the analysis of wind data, and all the staff members of the Hazaki Oceanographical Research Station for their contributions to the field measurements.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Correcting Missing Wind Data
First, the missing values of the cross-shore and alongshore wind velocities were obtained by linearly interpolating the measurements when values three hours before and after a missing value were available. Then, the other missing values were corrected using the linear correlations from January 1987 to November 2014 between the wind velocities measured at Hasaki and those reproduced from past weather at the point of N36.25, E141.25, which is closest to the Hasaki Coast. The reproduced values were part of the results of the JRA-55 project conducted by the Japan Meteorological Agency [53]. The correlations were expressed as: where u is the wind velocity; the subscripts c and l indicate the cross-shore and alongshore values, respectively; and the subscripts m and r indicate the measured and reproduced values, respectively. The coefficient determinations for cross-shore and alongshore velocities were 0.41 and 0.25, respectively, and the RMS errors for cross-shore and alongshore velocities were 4.10 and 4.66 m/s, respectively. The correlations between the reproduced and measured values for cross-shore and alongshore wind velocities are shown in Figures A1 and A2, respectively. Because the temporal interval of the reproduced values was six hours and not three, the missing values between the corrected values were replaced with the means of the corrected values. The numbers of all data, the values obtained using the interpolation, and those using the correlation between the measured and reproduced values were 81, 815, 62 (0.08% of all data), and 529 (0.65%), respectively.   Tables   Table 1. The RMS error, the correlation coefficient, and BSS for the weekly averaged volume-change rate and the cumulative volume change in the calibration process. Table 2. The RMS error, the correlation coefficient, and BSS for the weekly averaged volume-change rate and the cumulative volume change in the validation process.                 Figure 18. The cumulative volume changes at the end of the investigation period simulated using the values of K, c1, and c2 that were increased and decreased by 10%. The broken line shows the volume simulated using the best-fit parameter values.   Figure A1. Correlation between reproduced and measured cross-shore wind velocities. The black line represents the linear correlation estimated by the least-squares method. Figure A2. Correlation between reproduced and measured alongshore wind velocities. The black line represents the linear correlation estimated by the least-squares method.