Sinusoidal shaking in event-driven simulations

https://doi.org/10.1016/j.cpc.2011.09.005Get rights and content

Abstract

Event-driven algorithms are a powerful and efficient simulation method that can be used to numerically study the behavior of particulate systems, such as granular materials or powders. Its superior efficiency originates from only considering analytically calculable points in time, where physically relevant events occur. This, however, is only possible for pairs of trajectories whose relative distance can be expressed as a polynomial of order less than or equal to four. For more complex trajectories, like a sinusoidal motion of a wall, the collision times have to be calculated numerically, which dramatically reduces the efficiency. In this paper we present a tabulation method implemented in an event-driven algorithm, which is able to efficiently handle the one-dimensional motion of a periodically moving wall. By demonstrating this for a sinusoidally moving wall, we enhance the comparability of simulations with (usually sinusoidally driven) experiments. Because granular materials or powders are systems which are easily accessible through both experiments and simulations, our work enhances the growing field of statistical physics far from equilibrium.

Introduction

Whether in nature, as sand or coarse gravel, or in industrial processes and products like mining, food processing, or the pharmaceutical industry, granular materials or powders are omnipresent. Because these systems of particles are out of equilibrium, the concepts and tools known from classical thermodynamics, like phase diagrams or thermodynamic potentials, unfortunately do not universally apply in these cases. Due to the multitude of individual particles, it is a challenge to predict the behavior of such collective systems, and, for example, phase diagrams have to be extracted laboriously for any non-equilibrium system. Thus, it is advantageous to develop simulation techniques which allow a fast and efficient exploration of the desired bulk behavior.

A suitable and efficient1 simulation method is the event-driven molecular dynamics (EDMD) [1] simulation. However, as pointed out in [2], it has never become as popular as the less efficient time-driven molecular dynamics (TDMD) [1] simulation method. This seems to be due to an often higher implementation effort combined with less flexibility in EDMD. However, interesting implementations for EDMD were published in recent years, which extend the range of application of this method. For example, there are algorithms for parallelization [2], the use of non-spherical particles [3], [4], [5], a priority queue with complexity O(1) [6], the collision detection of moving spheres with unknown trajectories [7], and the use in a ball mill simulation [8]. EDMD simulations are not exclusively used in the field of physics but are also investigated in theoretical computer science in the context of Kinetic Data Structures [9].

In this paper we discuss a tabulation scheme to efficiently simulate moving walls of infinite mass, which are assumed to be flat. Its trajectories have to be mathematically smooth, periodic and finite. The motion of the wall itself is one-dimensional, but can be coupled to a system of any dimension. We demonstrate our method with a sinusoidally moving wall confining a system with circular (two-dimensional) or spherical (three-dimensional) particles. Sinusoidally moving walls in EDMD simulations were used already more than a decade ago [10], [11]. However, in a one-dimensional system the number of particles which can possibly collide with the boundaries is constant and less or equal to two while in higher dimensions all particles can possibly collide with the wall. Thus, in one-dimensional systems the sinusoidal motion does not create a performance issue, but in system with higher dimensions the computational complexity of the sinusoidal wall motion increases with the number of particles. Therefore, the method presented in this paper provides a significant speed-up for systems in dimensions larger than one. Due to the high implementation effort of a sinusoidal moving wall and the loss of computational efficiency the motion of a wall was often approximated by polynomials of the order two, or even less (see e.g. [12]). Sometimes (when sinusoidal motion is requested) even the much slower TDMD simulation method is used [13].

We increase the computational efficiency of simulating e.g. a sinusoidally moving wall in EDMD simulations by one or two orders of magnitude in order to strengthen the usage of sinusoidally moving walls compared to approximations with polynomials. Yet, the shape of the driving function becomes important at the point where one wants to compare simulations with experimental results: experiments on wet granulates [14], [15] have shown that the second derivative of the driving function is significant for the nonlinear response of the sample. Even though the method is fast its implementation is rather involved and only possible for periodically moving, straight walls.

The paper is structured as follows. In Section 2 we define the mathematical expression of which we need to find the roots. A first implementation, which partially reviews textbook knowledge, is shown in Section 3. The main novel work of this paper is contained in Section 4 where the tabulation method is presented in great detail. This is followed by Section 5 where necessary numerical improvements are discussed which assure the physical correctness of the results. The superior performance of our tabulation method is shown in Section 6 where numerical simulations are performed. Finally we summarize in Section 7.

Section snippets

The mathematical task

The basic idea of EDMD is first to determine the time where the next event occurs, where events are particle–particle collisions or collisions with a wall. The positions of the particles involved in these events are then easily computed in a second step of the evolution along their trajectories for the previously calculated time period. The velocities before the collision are calculated the same way and after applying the collision laws in the final step we obtain the final velocities after the

A first algorithm

After finding the roots we will introduce the minimal ingredients to successfully perform the EDMD simulation, and we will recognize at the end of this section that we lose a lot of the EDMD efficiency, which we will subsequently enhance in Section 4.

Improving the efficiency with look up tables

The main purpose of this task is to speed up the simple algorithm that was discussed in the previous section. Because the search of unique bracketing bounds is very expensive, this is the critical point. We have already seen that finding one bracketing bound cannot be further accelerated with the linear search method, so we reduce the number of calculations by storing each result in a look up table. This, of course, only saves computation time if we are able to reuse the stored values as often

Essential numerical improvements

In the previous section we have seen the essential concepts and methods that are necessary for simulating components with trajectories of arbitrary functional form, shown for a sinusoidally moving wall, with the help of look-up tables. In this section we will additionally show some essential numerical enhancements in order to couple this approximated dynamic of the wall collisions to the mathematically exact calculation of particle–particle interactions in the bulk.

Validity and comparison of performance

We have made a significant effort to create a robust algorithm. The algorithm described in this paper was already successfully used in several publications like in [15], [18], [19], [20] or in the limit of zero amplitude where it recovers exact results in [21], [22], [23]. Especially we stress [15], where we compared the results also with a TDMD simulation and found superior agreement even though the physical model also has to be changed slightly for the EDMD model. This shows that the

Summary

An algorithm is presented to simulate flat walls with periodic trajectories in event-driven molecular dynamics simulations. The details are shown for the example of a sinusoidally moving wall. It is based on a root finding method, which is briefly reviewed in Section 3.1. It is a combination of the Newton–Raphson method and Bisection method as described in [17]. The computationally most expensive task in this algorithm is to find a unique bracketing bound around the root, which leads to a

References (23)

  • S. Miller et al.

    J. Comput. Phys.

    (2004)
  • A. Donev et al.

    J. Comput. Phys.

    (2005)
  • A. Donev et al.

    J. Comput. Phys.

    (2005)
  • C. De Michele

    J. Comput. Phys.

    (2010)
  • G. Paul

    J. Comput. Phys.

    (2007)
  • T. Poeschel et al.

    Computational Granular Dynamics

    (2004)
  • H. Kim et al.

    Algorithmica

    (2005)
  • R. Reichardt et al.

    Granular Matter

    (2007)
  • S. Luding et al.

    Phys. Rev. E

    (1994)
  • S. Luding et al.

    Phys. Rev. E

    (1994)
  • Cited by (5)

    View full text