A NUMERICAL MODEL OF CROSS-SHORE BEACH PROFILE EVOLUTION: THEORY, MODEL DEVELOPMENT AND APPLICABILITY

A practical numerical model was developed to simulate cross-shore profile evolution at two coastal sites in Iran. The model consists of three sub-models for calculating wave and current, sediment transport, and bed level changes. Validation and calibration of the model was carried out using the measured field data on the north and south coasts of Iran, where historic measurements of cross-shore beach profiles and wave conditions have been recorded. The model is formulated for calculating cross-shore sediment transports in and outside the surf zone by the product of timeaveraged suspended sediment concentration under three different mechanisms and undertow velocity. The comparisons between the model results and field data show reasonable agreement for both coastal sites and will be capable of applying it to other coastal sites with modifications to the free parameters.


INTRODUCTION
The knowledge on nearshore morphological processes such as prediction of sediment transport rates and beach profile evolution is still limited in spite of past research studies carried out over the last few decades.Even though recent developments led to 3D process-based beach profile evolution models, they are still at early stage of development.Therefore, simplified models of beach profile changes provide coastal engineers helpful guidance and the model developers still apply best-fit technique with measured field and laboratory data in order to calibrate free parameters and improve their predictive models.The process-based models based on continuity and momentum equations combined with models of sediment transport and bed level changes are an extremely valuable tool for evaluating bottom changes in a beach profile.
In the present study, a process-based cross-shore beach profile evolution model consisting of three sub-models is applied to validate the model with field data at Nowshahr and Zarabad beaches in Iran.In the first section of this paper the field measurements and beach locations are summarized.In the second section, the cross-shore beach profile evolution model is presented that calculates wave energy decay, surface roller and bottom friction, wave breaking heights, radiation stresses, wave set-up and set-down, velocity profiles, suspended sediment concentration profiles, avalanching mechanism and bed level changes.Finally, the model setup, calibration and comparison results are presented.

BEACH LOCATIONS AND FIELD DATA
The study areas and field data (Fig. 1) of the present study are as follows: 1. Nowshahr coastal site: This site is located on the south coast of the Caspian Sea (Point 1).The field site in Nowshahr is about 13 km alongshore and about 1 km offshore coastal area.The data set consists of the measured cross-shore profiles and wave conditions for the storm events of 30 th October 2013 to 18 th January 2014.2. Zarabad coastal site: This site is situated south of Iran at Sistan and Balouchestan province on the Gulf of Oman (Point 2).A field monitoring program of periodic hydrography surveys was carried out from 2006 to 2008.The hourly time series of offshore spectral waves were adopted from 22years' hindcast data of the Gulf of Oman.

NUMERICAL MODEL DEVELOPMENT
The proposed cross-shore beach deformation model consists of a combination of three sub-models.The following sub sections (a-c) explain each sub-model briefly.

a) Wave and current model
The wave propagation sub-model computes the wave transformation including wave energy decay, surface rollers (  ) and bottom friction (  ).The cross-shore wave energy flux is assumed equal to the energy decay and bottom friction where  is the wave angle,  is the cross-shore distance.The wave energy is given by  = 1 8  2 , where  is the sea water density,  is the gravitational acceleration, and  is the wave height.The group velocity is given by   =  [ ], where , , and ℎ are the local wave celerity, wave number, and water depth respectively.The energy decay (  ) is given by the model of Whitford (1988) [Eqs.(2)-(3)]. (2) where  is the wave frequency,  is an empirical coefficient of order 1, and  is the breaker index.
The bottom friction is computed by: where   is the friction coefficient, and   is the near bed root mean square orbital velocity.The energy balance for the rollers can be formulated as follows: where the roller energy (  ) is given by: where   is the wave height at breaking,  is the slope of the wave front.
The roller dissipation (  ) is, where  is the roller coefficient which varies through the surf zone.
Radiation stresses produce a shoreward increase of sea level (wave set-up).Breaking waves and surface rollers are mainly two sources of radiation stresses.
The wave set-up (), is computed using the depth-integrated and depth-averaged cross-shore momentum balance equation.
The breaking height (  ) is an essential requirement for prediction of wave height transformation, as well as computation of radiation stresses which can be extended with additional momentum associated with surface rollers.The breaking height is found by Battjes and Janssen (1978) and given by Eq. ( 11).
Velocity profiles () throughout the water column are computed for both non-breaking and breaking waves.To obtain the distribution of the undertow induced by breaking waves precisely, the surf zone is divided into the inner and outer zones.The undertow distribution is derived using a formula for computing the shear stress distribution and eddy viscosity coefficient proposed by Okayasu (1989), where  is the vertical elevation,   is the depth at wave trough. 1 and  2 are constants depend on inner and outer zones and expressed as follows: where   and   are the locations of breaking point and the boundary of the outer zone measured from nearshore, respectively.In order to calculate the vertically averaged velocities a set of formulae proposed by Rattanapitikon and Shibayama (1996) is used.
where  0 is a parameter as shown below: where   is the bed slope. 3 is a constant depends on the different zone of the coastal environment.
Separate formulations for the suspended load and bed load are used in proposed study.Sediment concentration (  ) is predicted using a set of explicit empirical formulas developed by Jayaratne and Shibayama (2007) for the suspension on the bottom boundary layer over sand ripples, suspension from sheet flow layer and suspension under breaking waves.The horizontal flux of suspended sediment through a section is the product of suspended sediment concentration and velocity over the entire water depth.
1. Sediment suspension over rippled bed: Figure 2 illustrates two distinctive type of suspension layers over rippled bed.For the lower suspension layer ( ≤ 2), the bed reference concentration   is given by Eq. ( 15).
where,  1 is a numerical constant, Θ is the mobility number,  is the kinematic viscosity,  is the specific gravity of sand,  is the gravitational acceleration, d is the median grain size, and  is the rippled height.
The diffusion coefficient   is calculated from Eq. ( 16).
where  2 is a numerical constant,   is the orbital amplitude near the bottom,   is the settling velocity of sediment particles,  *  is the shear velocity,  is the ripple length, and  * = (   2 ) 1/3 is the dimensionless grain diameter of Van Rijn (1984).
() is the concentration profile over rippled bed,  =  2 , is given by For the upper suspension layer ( > 2), the bed reference concentration   is given by Eq. ( 15).
2. Sediment suspension over sheet flow: Figure 3 shows two types of suspension layers over sheet flow regime.For upper sheet flow layer, the following set of formulae is applied: where c(z) is the concentration profile over sheet flow layer.

Sediment suspension under breaking waves:
The following set of formulae is used to calculate the bed reference concentration, diffusion coefficient and concentration profile under breaking waves.
where  ̂ is the local wave orbital velocity and  ̂ is the wave orbital velocity at the breaking point. . 1 ,  2 , … ,  12 are numerical constants (Jayaratne et al (2014)).
The suspended sediment transport rate is the product of concentration and velocity distributions, Where  is the boundary layer thickness, () the time-averaged velocity profile which is assumed to be equal to the current velocity.
The bed load transport flux is computed on the basis Watanabe (1982) and is expressed as The net sediment flux is the sum of bed load transport flux and the current-related suspended load transport flux,  =   +   . (33)

c) Bed level changes
Beach profile evolution due to both bed and suspended loads is computed by the mass balance equation.To solve the differential equation, a finite difference method (FDM) is used.
where  is the time step,  is the grid size and  is the sediment porosity.The avalanching concept is considered when the critical wet ( , ) and dry ( , ) slope is exceeded.Critical wet and dry slopes that control avalanching, used as a free calibration parameter.

MODEL SETUP AND CALIBRATION
A computer program was written and performed in MATLAB.The initial bed profile, regular grid size with constant spacing (∆) were given to the computer program.Free calibration parameters such as breaker index (), roller coefficient index (), reference concentration constants (  ),  , and  , were changed for each test case.The Brier Skill Score (BSS) defined as in Eq. ( 35), which evaluates model skill for beach profile change as a good measure of performance of the coastal morphological models.The BSS can be used to evaluate whether the predicted profile is closer to the measured or the initial profile.
where   ,   and   are the simulated profile from the numerical model the measured profile and the initial bed profile respectively.A BSS ≤ 0 implies a very poor predictive skill whereas a BSS=1 indicates that a simulated profile perfectly matches the field results.The initial input parameters used for Zarabad and Nowshahr coastal sites are given in Table 1 and 2 respectively.

a) Application of the model to Zarabad coastal site
The construction of Zarabad fishery port was completed in 2006.For the period of 2006-2007, several severe storms occurred in Zarabad.The cross-shore profiles were surveyed from the top of the dune to the depth of closure.The input offshore wave characteristics during storms which were used to examine the model validity are shown in Figure 4.A median grain size was determined by sieving and found to be 0.02 cm.The offshore boundary of the model is located at the depth of closure, that is x= 600 m at about -7 m water depth.
The field measurements during the construction of Zarabad Fishery port indicated a positive gradient in the longshore sediment transport from west to east.As a consequence, sand is trapped on the east side of the port.The measured and initial beach profiles are shown in Fig. 5.It appears that net volume change is quite different from for the reference values of the profiles.In order to minimize this effect, which is considered to be due to longshore sediment transport (LST) gradients, the measured profile is adjusted by shifting the profile horizontally a distance of ∆ to yield zero net volume change.The value of ∆ can be calculated by the Eq. ( 36).
where subscripts  and  denotes initial and measured profile elevation, respectively. 0 and  ∞ are the offshore distance coordinates at the baseline and offshore profile change limit, respectively, and ℎ  is the total elevation of the measured profile.Figure 6 shows comparison between the simulated and measured profile after adjustment of the profile towards the equilibrium.

b) Application of the model for Nowshahr coastal site
Bathymetric surveys for eight beach profiles were collected between 30 th October 2013 and 18 th January 2014.To Evaluate the model an 18-day period from 3 rd December 2013 with   = 3.32 ,   = 9.45  was used (Figure 8).The average grain diameter of the sand is 0.02 cm and it remains nearly uniform along the cross-shore profile.The location of the beach and cross-shore beach profiles used in this study are given in Fig. 7.As the longshore sediment transport is interrupted by the port breakwaters, an accumulation of sand occurs on the updrift (west) side of the port.The model was set up for profile 8 (Fig. 9), which is located on the east side of the port and longshore sediment transport gradient was found to be small.Therefore, the most significant changes in the beach profile are linked to the cross-shore processes.

RESULTS AND DISCUSSION
The model application to two different coastal sites in Iran shows an overall good performance.The proposed model is sensitive to hydrodynamic and avalanching calibration parameters.Free parameters in reference concentration formulae have shown that they had an effective influence on the offshore sand bar formation.Calibration of beach profile evolution model indicates its high sensitivity to profile shape.As it is shown by existing beach profile models that systematically need calibration from one site to another, in the present numerical model performs in the same manner, The key limitation of the implementation of model was the lack of accurate field data with small gradients in the longshore sediment transport.Of the two study sites, profile measurements with negligible longshore gradient, were only available in data at Nowshahr site.Among eight profile measurements, there was only one where the gradient in the longshore transport is small.Thus, model specific parameters were calibrated based on limited data available.For Zarabad coastal site, the crossshore profile was adjusted by assuming that the longshore variations of coastline are only gradual and the cross-shore profile remains the same when the coastline shifts.More measurements are required to calibrate the model parameters accurately and validate the model performance.However, based on the numerical simulation results of Zarabad and Nowshahr beaches carried out using the proposed model, the existing model proved to be a useful practical tool to predict beach profiles at week to months' time scale.

Figure 4 .
Figure 4. Offshore measured time series of wave heights, wave periods, wave directions (June-August 2007) and wave rose for Zarabad coastal site.

Figure 5 .Figure 6 .
Figure 5. Location of selected cross-shore profile and comparison between initial and measured profile at Zarabad.

Figure 7 .
Figure 7. Location of measured cross-shore profiles at Nowshahr site.

Figure 8 .
Figure 8. Offshore measured time series of wave heights, wave periods, wave directions and wave rose for Nowshahr coastal site.