The mechanism of spontaneous corrugation on the snowy and icy roads produced by the moving vehicles in cold regions

Traffic safety in cold regions is seriously affected by the snow and ice brought by the extreme climate. The snowy and icy road cannot provide enough friction for the safe operation of vehicles due to its smooth and uneven surface. In this research, we are going to focus on the uneven corrugation occurred on snowy and icy roads and to investigate the formation mechanism of this spontaneous corrugation which can seriously threaten the traffic safety. According to field observations, we found that the corrugation phenomenon generated by moving vehicles is a complicated thermal–mechanical coupled process. In order to simplify this complicated process, we are going to focus on the mechanical process of the formation of spontaneous corrugation only at this stage. Field observation by time-lapse cameras has been conducted to disclose its forming process directly. Then, we adopted sand as the material to reproduce the spontaneous corrugation in the laboratory which can eliminate the influence of the thermal process. By considering the compressibility and mobility of the surface material comprehensively, a numerical model has been successfully constructed for imitating the forming process of corrugation. Then based on this proposed numerical model, a preliminary discussion on the influence of natural frequency on the number of the corrugation has been conducted. The relationship between the natural frequency which is decided by the vehicle itself and the corrugation is promising to be utilized in optimizing the vehicle design to improve the performance on the snowy and icy roads.


Introduction
Road corrugation or washboard road, which often spontaneously forms on unpaved roads under the action of plows or wheels, typically results in an uneven road surface.This unevenness can cause discomfort for drivers and passengers, lead to deterioration of road conditions [1], and increase the costs of maintenance obviously.This phenomenon has garnered extensive attention and been widely studied by scholars around the world.Previous studies have primarily focused on experimental and numerical approaches.For instance, Bitbol et al. [2] investigated the critical speed of formation of corrugation and proposed a Froude number to control the instability of ripples.Both the Finite Element Method (FEM) [3,4] and Discrete Element Method (DEM) [5,6] have been used in the study of washboard road mechanisms.In terms of mechanistic research, some scholars [1,7] postulate that the ripples on the road surface are instigated by vehicular vibrations.Ikeda et al. [8] designed an indoor experiment employing a seesaw-shaped oscillator to replicate ripples.The spring-mass model is often utilized for simulating the genesis of washboard road [1,9].In addition to unpaved roads, asphalt pavements can also exhibit similar phenomena [10,11].Audrius Vaitkus [12] believes that the strength of the bond in the asphalt layer is an important factor affecting pavement corrugations and has assessed the interlayer bond strength using a direct shear test.Hao Wang [13] established a three-dimensional finite element model to study the effects of vehicle acceleration and braking on pavement corrugations.
Similarly, this phenomenon has been observed on icy and snowy roads (see Fig. 1).As shown in Fig. 2, the formation of a washboard road on icy and snowy surfaces often involves complex tire-road coupling, accompanied by thermal transfer, weather conditions, and other factors.Compared to unpaved roads, the corrugations on icy and snowy roads pose a heightened threat to driving comfort and transportation safety.Each year in China and Japan, a significant number of accidents occur on snowy and icy roads, resulting in numerous fatalities and injuries [14,15].A primary cause of these accidents is the loss of friction between tires and road surfaces.Research by some scholars on tire-road coupling on icy and snowy roads has provided insights for safe driving in winter.Takahashi N et al. [16] conducted experiments using a wheel tracking test machine and proposed a method for predicting ice layer thickness based on pavement type, weather conditions, and traffic volume.Carlson et al. [17] has studied how rolling resistance is affected by water and snow on the road surface, providing better recommendations for winter road maintenance.
However, the underlying mechanisms and complexity of such a mysterious phenomenon as the formation of washboard patterns on icy or snowy roads remain not entirely comprehended.Consequently, the capability to predict its emergence and prevent its occurrence, crucial for ensuring traffic safety, continues to be an elusive target.
With the purpose of increasing our understanding for this corrugation on snowy and icy roads, this paper is structured as follows: In "Field observation" section, we examine the conditions that facilitate the spontaneous corrugation on snowy roads by continuous field observation."Reproduction experiment" section presents the experimental apparatus used in the laboratory which can reproduce the corrugation properly.Subsequently, in "Numerical simulation" section, we discuss the compressibility and mobility of the sand and then propose a numerical model to simulate the forming process of corrugation observed during the indoor experiment.Finally, in "Numerical analyses" section, we utilize this numerical model to conduct a numerical experiment to further study spontaneous corrugation and its influencing factors.

Field observation
In order to investigate the formation of the spontaneous corrugation on icy and snowy road, we selected an observation point that produced corrugation frequently in winter, to monitor the road surface for 2 winter seasons (2016-2017).The observation point selected is at a toll gate on the campus of Hokkaido University where we used two time-lapse cameras on either side of the gate.With a serious and through analysis on the photographs taken by time-lapse cameras, traffic data, and the meteorological data collected from the meteorological station, it is founded that three specific conditions are necessary [18] for the occurrence of corrugation on snowy and icy roads.We selected a representative week, spanning January 19 to January 25, 2016, for data Fig. 1 The spontaneous corrugation on the snowy and icy roads Fig. 2 Sketch map of corrugation on snowy and icy roads leading to the loss of friction analysis.During this week, corrugations were observed on two days.The detailed conditions are explained as follows: First, a substantial amount of snowfall is required because snow is the basic material for the corrugation.As shown in Fig. 3, on the two days when corrugations occurred, snow depth exceeded 55 cm on both occasions, and there was an average snowfall of 4 cm before these days.Subsequently, there should be a declining trend in accumulated snowfall, with the temperature falling between 0 and -2°C following the snowfall.As clearly visible in Fig. 4, when the temperature exceeded 0 ℃ and dropped below -2 ℃, there were no corrugation phenomena observed, even when the snowfall and traffic conditions were similar to when corrugations occurred.For example, on January 25th, corrugations did not appear.This temperature range indicates a relatively warm condition which can promote the thawing process of snow, especially considering the heat brought by the moving vehicles.Lastly, a certain amount of traffic can facilitate the forming and growing of corrugation.As our observation, the moving vehicles not only provide vertical compaction and horizontal pushing effects for snow, but also brings a large amount of heat to promote the phase change process of ice-water on the snowy road surface.
Next, we will discuss the formation of corrugations on a particular day, 19th January 2016.
Figure 5 shows the formation of road corrugations observed using a time-lapse camera (time period: 11:59 AM to 6:08 PM, temperature: -2 to 0℃).From Fig. 5a to b, it can be observed that there was no obvious corrugation occurred between 12 AM to 3 PM.However, by 5 PM, roughly five hours later, the corrugations become quite significant.By 6:08 PM, the corrugations become stable and no further growth was observed.During the period between Fig. 5a and b, traffic volume was approximately 70-80 vehicles per hour, while from Fig. 5c to d, the traffic volume exceeded 140 vehicles per hour.From this observation example, it is demonstrated that the three main factors influencing the formation of corrugations: snowfall, ambient temperature, and traffic volume.An increase in traffic volume indeed promotes the growth of corrugations.Besides the loading efforts of moving vehicles, we also noticed a difference in the temperature distribution between the road surface and the tires (see Fig. 6).When vehicles are idling or moving slowly, the heat from the tires or exhaust would transfer to the road surface, leading to a snow-water phase transition.In this process, we posit a hypothesis that transient heat transfer from the tires to the snow surface, whether Fig. 3 The data of fresh snowfall between January 19 to January 25, 2016, red boxes are the corresponding data when the corrugation occurred due to frictional heat generation or tire heat transfer, induces the surface snow melting.Under the influence of air temperature, this meltwater re-solidifies into new ice, redistributing the height of the ice and snow.After a certain number of vertical compaction and horizontal pushing produced by the moving vehicles, a washboard pattern in the ice and snow surface is thereby formed.In other words, the formation of corrugation on snowy and icy roads is a complicated thermal-mechanical coupled process accompanied by moving vehicles.To fully understand this complex and interconnected process, more detailed and small-scale experimental studies are desired.
Based on the results of field observation, the causes of the formation of ice and snow ripples may involve two aspects: one is the mechanical behavior of snow and ice under moving vehicle load, and the other is the phase transition effect of snow and ice in contact with the tires.While it is challenging to alter the environmental factors such as snowfall and ambient temperature by human intervention, these conditions can be predicted using weather forecasting models.In light of this, preventative measures could be strategically implemented in areas where the conditions tending to the formation of corrugations.Such measures may include traffic regulation, pavement heating, or modifications of vehicles, aiming to mitigate or even eradicate the adverse effects of these corrugations on the snowy and icy road.

Reproduction experiment
In our field observation, we discovered that the formation of corrugation results from a complex thermomechanical coupling process induced by the mechanical action of vehicles on snow surfaces and heat transfer.Consequently, in this experiment, we focus primarily on the mechanical behavior of snow and ice under vehicle load, considering non-phase transition material as a substitute for snow to more clearly investigate the mechanical mechanisms involved in this process more clearly.It can provide promising insight into adjusting vehicles to alleviate corrugation.

Experiment apparatus
As shown in Fig. 7, the experimental apparatus consists of three parts: a self-rotating roadbed track, a seesaw-shaped oscillator, and a laser displacement sensor.Figure 7a shows a self-rotating track with a diameter of 50 cm.It can achieve a rotation speed of up to 20 revolutions per minute (rpm), corresponding to a velocity of 0.52 m/s.It has a width of 9 cm and can accommodate a 5-cm thick roadbed layer.The purpose of designing the rotating track is to create a relative motion between a position fixed oscillator and the track, simulating the forward movement of a vehicle.Different rotation speeds can simulate the effects of different vehicle speeds on the formation of ripples.The track is equipped with a dry sand bed, with the sand material selected being Toyoura sand.Toyoura sand has an average particle size of 0.2mm and is the standard test sand in the Japanese Industrial Standards (JIS).Table 1 shows the physical parameters of Toyoura sand.With the help of sand, we can remove the influence of phase change factor from materials of snow or ice and then focus on the mechanical process.Figure 7b shows a sectional view of a seesaw-shaped oscillator fixed above the track, capable of rotating around its fixed end.The vehicle can be simplified to the oscillator like a seesaw mainly contains three parts: counterweights, a plastic cylinder and a spring.One end of the seesaw arm is connected to a cylinder, which serves to simulate a vehicle's wheel without rolling.The other end of arm loads a changeable pendulum to adjust the natural frequency of the oscillator system.The middle part of the seesaw is mechanically suspended by a spring.Therefore, through the relative motion between the oscillator and the track, the system can imitate the interaction between the simplified running vehicle and the sand bed.A laser displacement sensor with a precision of 0.01 mm is fixed in front of the oscillator (see Fig. 7b) to measure the height of the sand surface as the track rotates.The measurement frequency is 100 Hz.

Experimental procedure
Before conducting the experiment, preparatory work is necessary.The sand surface is leveled using a scraper, and Fig. 6 The temperature distribution at the observation point with a moving vehicle using a thermography instrument the initial height of the sand surface is measured by the laser sensor.Subsequently, the oscillator is gently placed on the sand surface.The track rotates slowly at a speed of 3.0 rpm while measuring the height of the sand surface, continuing for two minutes.Without stopping the track, the speed is adjusted to the designated value.The track runs at a constant speed for 7 min, during which the height of the sand surface is continuously recorded to observe the development of the washboard road.After stopping the experiment, a top-view photograph of the track is taken, and the number of ripples present on the track is calculated.

Numerical simulation
This section aims at developing a numerical model to simulate the formation of spontaneous corrugations.Through the model, we can understand the process of corrugation formation better and carry out sensitivity analyses by adjusting different parameters conveniently.A better understanding of the formation mechanism of the corrugation in detail can facilitate to more effective measurements and corresponding strategies to eliminate the corrugations on the snowy and icy roads, to promote the safety and comfort of transportation.

Spring-mass model
For the above indoor experiment, it has been observed that the frequency of the oscillator triggers the formation of corrugation on the sand surface [8,19].Hereafter, we are aiming to develop a numerical simulation model to further simulate the behavior between oscillator and the surface.An oscillator with a damped spring, damper as well as mass M is used in this dynamic model (see Fig. 8a).Therefore, we introduced a one degree of freedom dynamic equation model to simulate motion of the oscillator, as shown in Fig. 8b.
Dynamic oscillation equation at the equilibrium condition can be written as where F (t) is the reaction force of the sand surface when the attachment touches the sand, and F S , F M and F C are spring force, inertia force and damping force respectively.Introducing each parameter and the time varying force, dynamic oscillation equation at the equilibrium condition and can be written as where M is the mass of the oscillator, c is the damping coefficient, k is spring constant.The schematic drawing of the reaction force is shown in Fig. 9, where x(t) repre- sents the vertical position of the attachment over time t.In motion process of the oscillator, we divide the motion into two cases.In case 1, the attachment loses contact with the sand surface and moves freely in the air.Therefore, there is no reaction force acting on the attachment (F (t) = 0) .On the contrary, in case 2, the attachment penetrates into the sand, which leads to a deformation of sand at the touch point, receiving an upward reaction F (t) from the sand.Now, we need to calculate the (1)  where r is the radius of attachment, z(t) is the height of the sand surface.We define that the deformation of the sand as d .The deformation d and the force F (t) have a linear relationship which is verified in the penetration test [8].The linear relation is introduced in this step and the equation is as follows: (3) where A S (unit: N/m) is the deformation coefficient of the sand measured by the penetration test which depends on the type of sand, and parameter b is a con- stant.Therefore, the deformation d(t) can be calculated by the Eq. ( 4).Additionally, compressibility and mobility of sand is considered.At the contact point n, the deformation of the sand d(t) leads to a new distribution of the sand surface, as shown in Fig. 10.The coordinates of each point surrounding the contact point n can be calculated via geometric relationships.Initially, we assume that the penetration deformation d and the rising height d′ of the sand surface near the touch point n are equal.Considering that the height of the sand pile will not increase infinitely during the collision process between the attachment and the sand surface, we introduce an equation [19] controlling the maximum height and depth of the sand surface.The equation is as follows where v is the velocity of the track (unit: rpm, round per minute), h(v) is the maximum height (unit: m) when the velocity is v.At the onset of the simulation, an initial excitation must be provided to the attachment.This can either be an initial upward force or an initially undulating surface, either of which can promote the formation of ripples on the sand bed.It is noteworthy that corrugations never form on an absolutely flat sand surface over time.Therefore, we use the following equation to control the initial height of the sand surface.The height of the sand surface at a given time t is denoted by where B is the amplitude of the sine wave, n 0 is the initial wavenumber, R is the radius of the track, v 0 is the velocity of the track (unit: m/s).

Model validation
The proposed model is implemented using the Finite Difference Method (FDM).The time step is set as 0.01 s, and the total duration of the experiment is 720 s.To examine the efficacy of the numerical model, we set parameters akin to those in the indoor experiment and simulated the wave number of corrugations at different velocities using numerical models for oscillators of three distinct (5) natural frequencies (1.37, 1.46, and 1.66 Hz).The results, as depicted in Fig. 11, reveal only a slight discrepancy in the wave number between the indoor and numerical experiments.Both indoor and numerical experimental outcomes exhibit consistent trends with velocity changes, and the two experimental datasets align closely within a limited range.This corroborates the high accuracy of our numerical model.

Numerical parameters
For a vehicle, the suspension system can be adjusted according to different driving modes, offering more flexibility compared to the vehicle's weight to provide a more safe and comfortable driving experience.Similarly, in the numerical model, we can modify the spring to alter the natural frequency.As shown in Table 2, the following parameters were used in the numerical analysis which were decided to close to the indoor experiment: mass (M) = 100 g, rotating velocity = 12 rpm (0.314 m/s), and damping coefficient (c) = 0.10 N • s/m.The natural frequency can be adjusted by varying the spring coefficient, ranging from 2.52 Hz to 5.63 Hz.The total duration of the experiment is set to 720 s, A S = 800, b = 0, and A S = 200.The results will be presented in "Results and discussion" section.

Results and discussion
Based on the parameters in Table 2 as well as using the numerical model, we have carried out numerical experiments to study the relationship between natural frequency and corrugation.
In the current experimental configuration, the selfrotating track operates at a speed of 12 rpm, achieving a full rotation every 5 s.This numerical simulation was conducted for 720 s, during which the track made 144 rotations.The circumference of the track used aligns with that of the indoor experiment described in "Reproduction experiment" section, approximately 157 cm.The road profile data, obtained from the first and the last rotations, were graphed to generate Fig. 12. From Figs. 12a to 12d, as the natural frequency of the oscillator gradually increases, a corresponding increase in the number of waves can also be observed.Moreover, it was discovered that the final height and depth of the sand piles formed by oscillators of four different frequencies were nearly identical.This is due to the fact that in our model, the function limiting maximum height is a velocity-dependent function.
To observe the changes in the sand surface throughout the experimental process, time series data at a rotation speed of 12 rpm with natural frequencies of 2.52 Hz, 3.56 Hz, 4.36 Hz, and 5.03 Hz were converted into the frequency spectrum, as shown in Fig. 13.In this figure, the position of the self-rotating track is represented by the horizontal axis, which covers a single rotation period,  From Fig. 13, it is intuitively evident that both the transient state and forward-moving mode have been successfully simulated, further validating the accuracy of the model.As the frequency increases, the transient state (see Fig. 13 red dashed line box) shortens, suggesting that an increase in frequency accelerates the growth rate of the ripples.Simultaneously, it also increases the density of the ripples and decreases their wavelength.These phenomena can be attributed to an insufficient rise of the oscillator caused by an increase in its spring coefficient, which, in turn, raises the probability of collisions between the oscillator and the sand surface in each cycle.Ultimately, these factors contribute to the emergence of the observed interesting corrugation.
The natural frequency is related to the mass of the oscillator and its spring coefficient, as shown in Eq. ( 7), where f is the natural frequency of the oscillator, m is the mass of the oscillator, k is the spring constant of the oscillator.For a fixed-mass oscillator, its natural frequency is positively correlated only with the spring constant.Therefore, we altered the natural frequency solely by changing the spring coefficient of the oscillator.Figure 14 illustrates the relationship between the wave number and the oscillator frequency under three different rotation speeds.Between 1 and 5 Hz, the data points can be fitted as follows: where n is the wave number, α and β are constants, with α = 10.5170 and β = 0.8664 at 9rpm, α = 8.0985 and β = 0.8491 at 12rpm, as well as α = 6.2785 and β = 0.8748 (7) 14, under conditions of specific velocity (9rpm), a positive correlation exists between the wave number and the natural frequency; specifically, the wavenumber escalates concurrently as the natural frequency increases.Conversely, with a fixed natural frequency, the wave number and velocity display a negative relationship which suggests that an elevated rotational speed results in a reduced wavenumber.Through the aforementioned numerical experiments, it is inferred that the natural frequency of the oscillator significantly influences the formation of spontaneous corrugations.In this simulation, we manipulated the natural frequency by adjusting the spring coefficient.If the aim is to minimize the occurrence of corrugations on the sand surface, applying a spring with a relatively small constant is advisable.In real-world terms, the springs equate to the suspension system in a vehicle.Therefore, if feasible, softening the vehicle's suspension system may alleviate the formation of corrugations on the road if only the mechanical process was taken into account.

Conclusions
This study clarified the occurrence conditions of spontaneous corrugation on snowy and icy road.On the basis of the understanding of field observation, we developed an experimental apparatus to examine the effect of natural frequency on spontaneous corrugation.Sand as the test material effectively avoids the influence of phase change factors on the experiment.Then, by considering the compressibility and mobility of the surface material, we successfully developed a numerical simulation model to imitate the corrugation process.The validity of the model was examined through comparison between the calculation results and the indoor test results.
The findings in this study are enumerated below: 1) The occurrence of ripples on icy and snowy road is related to the amount of snowfall, ambient temperature, and traffic volume, and is accompanied by a complicated thermal-mechanical coupled process accompanied by moving vehicles.2) By considering the compressibility and mobility of the surface material, we successfully constructed a numerical model for imitating the corrugation process even if it may need some modification on the movement of surface material for more precise simulation.3) By numerical simulation, it is clarified that natural frequency can affect the number of corrugations which is a critical factor for traffic safety and comfortability, especially in cold regions.
The numerical simulation successfully reproduced the process of gradually generating ripples on the sand surface, with a high accuracy in a certain range.We found that the natural frequency of the oscillator is an important influencing factor.In the future, we will attempt to use snow as experimental material to study the similarities and differences in the formation process of sand ripples and snow ripples.When more causes of formation are discovered, it will help eliminate or avoid the corrugations on the snowy and icy roads during winter.

Fig. 4 Fig. 5
Fig.4 The data of traffic and ambient temperature between January 19 to January 25, 2016, green boxes are the corresponding data when the corrugation occurred

Fig. 7
Fig. 7 The experimental apparatus: a shows the self-rotating track; b shows the profile of oscillator system

Fig. 8 Fig. 9
Fig. 8 Sketch drawing of oscillator system in spring-mass model (a) and free body diagram of the top attachment of oscillator (b)

Fig. 10 1 Fig. 11
Fig. 10 Schematic drawing of the deformation surrounding the touching point n

Fig. 12
Fig.12 comparison graph of sand surface with four different natural frequencies, where the vertical axis represents the height of the sand surface, and the horizontal axis represents time, with a cycle of 5 s

Fig. 13
Fig. 13 Four samples of running spectrum of 2.52 Hz (a), 3.56 Hz (b), 4.36 Hz (c) and 5.63 Hz (d), respectively, at the track rotating velocity of 12.0 rpm after 720 s of simulation

Fig. 14
Fig. 14 Relationship between natural frequency and wave number

Table 1
Physical properties of Toyoura sand reaction force F (t) at time t .Here, introducing a theoreti- cal coefficient A S ( unit: N / m) to control the force F (t),

Table 2
Parameters of the numerical simulation