^{1,a)}, David M. Walker

^{2}, Antoinette Tordesillas

^{2}and Chi K. Tse

^{3}

### Abstract

For a given observed time series, it is still a rather difficult problem to provide a useful and compelling description of the underlying dynamics. The approach we take here, and the general philosophy adopted elsewhere, is to reconstruct the (assumed) attractor from the observed time series. From this attractor, we then use a black-box modelling algorithm to estimate the underlying evolution operator. We assume that what cannot be modeled by this algorithm is best treated as a combination of dynamic and observational noise. As a final step, we apply an ensemble of techniques to quantify the dynamics described in each model and show that certain types of dynamics provide a better match to the original data. Using this approach, we not only build a model but also verify the performance of that model. The methodology is applied to simulations of a granular assembly under compression. In particular, we choose a single time series recording of bulk measurements of the stress ratio in a biaxial compression test of a densely packed granular assembly—observed during the large strain or so-called critical state regime in the presence of a fully developed shear band. We show that the observed behavior may best be modeled by structures capable of exhibiting (hyper-) chaotic dynamics.

Inferring dynamical structure of a deterministic nonlinear system from noisy time series data is a long standing problem. Here, we demonstrate the application of a range of methods to a problem of particular interest in the study of granular materials. An initially dense assembly of grains under biaxial compression with constant confining pressure reaches a near steady-state behavior (modulo fluctuations) in the large strain regime. The dynamics of this system at this limit are the focus of this paper—and of significant interest within the soil mechanics literature. Understanding the dynamical response of the granular structure in this regime is crucial for prediction and control in many industrial and geophysical systems. Our analysis indicates that the dynamics of this system in the presence of a fully developed shear band are consistent with a form of transient chaos—characterized by a bistable system which switches between two chaotic regimes, one with large slow oscillations and one with smaller and faster dynamics.

This work was funded by a Hong Kong University Grants Council grant under the General Research Fund: grant number PolyU 5262/11E. This work was partially supported by US Army Research Office (W911NF-11-1-0175), the Australian Research Council (DP0986876 and DP120104759) and the Melbourne Energy Institute (AT, DMW). MS is supported by an Australian Research Council Future Fellowship (FT110100896) and would like to thank the Melbourne Energy Institute for travel support. We thank Dr John Peters for useful discussions and Mr Tuan Tran for the software implementation to calculate local Lyapunov exponents.

I. INTRODUCTION

A. Modeling data and inferring dynamics

B. Dynamical behavior in granular assembles

II. DEM BIAXIAL COMPRESSION TEST

III. COMPUTATIONS

A. Embedding

B. Modeling

C. Model testing

IV. RESULTS

V. CONCLUSIONS

### Key Topics

- Time series analysis
- 44.0
- Chaotic dynamics
- 15.0
- Chaos
- 12.0
- Phase space methods
- 12.0
- Granular systems
- 11.0

## Figures

**The Canadian Lynx time series: chaotic or periodic?** The upper panel shows the original data (open circles) from which various global nonlinear models have been built. The lower panel depicts the dynamical behavior of three such models. The horizontal axis in the upper panel is the observation year: one data point per year. In the lower three panels the horizontal axis are in units of model time steps: 10 time steps per year. Hence, the upper panel covers a period of 115 years, the lower panels cover a period of 300 years. The vertical scales are arbitrary (derived from the total quantity of lynx pelts harvested in a given year). The dynamics observed in the lower three panels are: chaotic, almost periodic and exactly periodic.

**The Canadian Lynx time series: chaotic or periodic?** The upper panel shows the original data (open circles) from which various global nonlinear models have been built. The lower panel depicts the dynamical behavior of three such models. The horizontal axis in the upper panel is the observation year: one data point per year. In the lower three panels the horizontal axis are in units of model time steps: 10 time steps per year. Hence, the upper panel covers a period of 115 years, the lower panels cover a period of 300 years. The vertical scales are arbitrary (derived from the total quantity of lynx pelts harvested in a given year). The dynamics observed in the lower three panels are: chaotic, almost periodic and exactly periodic.

**Biaxial compression test.** (Top) The observed time series of the bulk measurement of the stress ratio with respect to axial strain . The strain interval of interest is indicated by the bolder blue trace and covers the high strain post-peak regime (so-called critical state), where the material has failed and the response is in an approximately steady-state, exhibiting characteristic stick-slip (jamming-unjamming) dynamics. (Bottom) Although not used in the modeling it is informative to soil mechanicians to observe the volumetric strain response evolution with respect to axial strain. The sample initially contracts and then dilates before reaching a more steady state response in the high strain region from where models in a reconstructed phase space of time-delay stress ratio variables are built.

**Biaxial compression test.** (Top) The observed time series of the bulk measurement of the stress ratio with respect to axial strain . The strain interval of interest is indicated by the bolder blue trace and covers the high strain post-peak regime (so-called critical state), where the material has failed and the response is in an approximately steady-state, exhibiting characteristic stick-slip (jamming-unjamming) dynamics. (Bottom) Although not used in the modeling it is informative to soil mechanicians to observe the volumetric strain response evolution with respect to axial strain. The sample initially contracts and then dilates before reaching a more steady state response in the high strain region from where models in a reconstructed phase space of time-delay stress ratio variables are built.

**False nearest neighbors.** The upper two panels depict the proportion of false nearest neighbors (embedding lag 1) for the data depicted in Fig. 2 . For there are no false nearest neighbors, the proportion drops to almost 0 for . The lower panel is the proportion of local false nearest neighbors, a measure of the locally sufficient embedding dimension, for local dimension (again, embedding lag is 1 and the prediction horizon is 3). The calculation is repeated for various n = local neighborhood sizes and plateau onset at a local dimension of is evident. The local Lyapunov exponent spectrum indicated four genuine exponents with one or two positive.

**False nearest neighbors.** The upper two panels depict the proportion of false nearest neighbors (embedding lag 1) for the data depicted in Fig. 2 . For there are no false nearest neighbors, the proportion drops to almost 0 for . The lower panel is the proportion of local false nearest neighbors, a measure of the locally sufficient embedding dimension, for local dimension (again, embedding lag is 1 and the prediction horizon is 3). The calculation is repeated for various n = local neighborhood sizes and plateau onset at a local dimension of is evident. The local Lyapunov exponent spectrum indicated four genuine exponents with one or two positive.

**Linearly filtered noise.** Comparison of nonlinear statistics computed with the Gaussian Kernel algorithm (correlation dimension, entropy and noise level), and higher order linear distribution statistics (skewness and kurtosis), for the original data in Fig. 2 and linear surrogates (monotonic nonlinear transformations of linearly filtered noise). As expected the linear statistics (lower right plots) show no difference. However, the nonlinear measures indicate that the linear model does not adequately described the dynamics (upper panels and lower left). For the nonlinear measures, the solid blue lines (no error bars) indicate the statistic values computed for the data (as a function of the embedding dimension—a parameter of the statistic). The tight error bars (red) are the mean andstandard deviation from 100 simulations, the larger error bars (green) are the full range (minimum to maximum). For the higher order linear statistics, a distribution of values is plotted, the value for the data are indicated as an asterisk.

**Linearly filtered noise.** Comparison of nonlinear statistics computed with the Gaussian Kernel algorithm (correlation dimension, entropy and noise level), and higher order linear distribution statistics (skewness and kurtosis), for the original data in Fig. 2 and linear surrogates (monotonic nonlinear transformations of linearly filtered noise). As expected the linear statistics (lower right plots) show no difference. However, the nonlinear measures indicate that the linear model does not adequately described the dynamics (upper panels and lower left). For the nonlinear measures, the solid blue lines (no error bars) indicate the statistic values computed for the data (as a function of the embedding dimension—a parameter of the statistic). The tight error bars (red) are the mean andstandard deviation from 100 simulations, the larger error bars (green) are the full range (minimum to maximum). For the higher order linear statistics, a distribution of values is plotted, the value for the data are indicated as an asterisk.

**Stable node.** The data displayed here are in the same format as Fig. 4 , except here the 100 simulations come from a model of the data driven by noise. In the absence of noise the model exhibits a stable node.

**Stable node.** The data displayed here are in the same format as Fig. 4 , except here the 100 simulations come from a model of the data driven by noise. In the absence of noise the model exhibits a stable node.

**Stable focus.** The data displayed here are in the same format as Fig. 4 , except here the 100 simulations come from a model of the data driven by noise. In the absence of noise the model exhibits a stable focus.

**Stable focus.** The data displayed here are in the same format as Fig. 4 , except here the 100 simulations come from a model of the data driven by noise. In the absence of noise the model exhibits a stable focus.

**Transient chaos.** The data displayed here are in the same format as Fig. 4 , except here the 100 simulations come from a model of the data driven by noise. In the absence of noise the model exhibits a chaotic transient (typically over a time scale longer than the observed data) and a stable fixed point.

**Transient chaos.** The data displayed here are in the same format as Fig. 4 , except here the 100 simulations come from a model of the data driven by noise. In the absence of noise the model exhibits a chaotic transient (typically over a time scale longer than the observed data) and a stable fixed point.

**Motif frequency figure.** For each of the three model dynamics we computed the motif super-family which occurred most frequently—shown here as a percentage. Motif family ADBCEF is the corresponding family for the data (the right-most batch of columns). In models exhibiting a stable focus, for example Motif ABDCEF occurred for all but one of the model simulations, that exception was ADBCEF. According to this measure, stable node and transient chaos models exhibited behavior most similar to the data. *Inset*: the structure of the six four node motifs. The motif-superfamily is a ranking by frequency of sub-graphs of order 4. The inset enumerates the six different motifs of order 4 which are possible. The horizontal axes of the upper panel is a ranking of these six classes (from most frequent to least).

**Motif frequency figure.** For each of the three model dynamics we computed the motif super-family which occurred most frequently—shown here as a percentage. Motif family ADBCEF is the corresponding family for the data (the right-most batch of columns). In models exhibiting a stable focus, for example Motif ABDCEF occurred for all but one of the model simulations, that exception was ADBCEF. According to this measure, stable node and transient chaos models exhibited behavior most similar to the data. *Inset*: the structure of the six four node motifs. The motif-superfamily is a ranking by frequency of sub-graphs of order 4. The inset enumerates the six different motifs of order 4 which are possible. The horizontal axes of the upper panel is a ranking of these six classes (from most frequent to least).

**Sample trajectories.** The top panel depicts the original data. The following four panels are iterated free run simulations, with dynamical noise for four different simulations from the model exhibiting transient chaos and the same motif super-family as the data.

**Sample trajectories.** The top panel depicts the original data. The following four panels are iterated free run simulations, with dynamical noise for four different simulations from the model exhibiting transient chaos and the same motif super-family as the data.

**Sample trajectory.** The top panel depicts a single noise free simulation exhibiting transient chaos dynamics. The lower panel is an embedding of the same data (color coded to depict different dynamical regimes). Note that the system appears to switch between two different dynamical behaviors before eventually collapsing to a fixed point. The time scale is the same as Fig. 9 .

**Sample trajectory.** The top panel depicts a single noise free simulation exhibiting transient chaos dynamics. The lower panel is an embedding of the same data (color coded to depict different dynamical regimes). Note that the system appears to switch between two different dynamical behaviors before eventually collapsing to a fixed point. The time scale is the same as Fig. 9 .

**Sample trajectory.** The left panel depicts a short section of the embedded data from Fig. 10 —the section is representative of the transient chaotic state and is the purely deterministic model output (no noise). On the right is the complex network constructed from this embedding.

**Sample trajectory.** The left panel depicts a short section of the embedded data from Fig. 10 —the section is representative of the transient chaotic state and is the purely deterministic model output (no noise). On the right is the complex network constructed from this embedding.

## Tables

**DEM** 2D biaxial compression test parameters and material properties.

**DEM** 2D biaxial compression test parameters and material properties.

**Complex network statistics.** For the original data and each of the three model classes, we computed observed values of the various complex network statistics: mean path-length (aka diameter), clustering, and assortativity. In each case, the mean value and standard deviation are reported. Also reported is the number of simulations (out of 100 attempts) which remained bounded. According to these measures both transient chaos and stable node dynamics are in agreement for two out of the three measures and differ marginally for the third.

**Complex network statistics.** For the original data and each of the three model classes, we computed observed values of the various complex network statistics: mean path-length (aka diameter), clustering, and assortativity. In each case, the mean value and standard deviation are reported. Also reported is the number of simulations (out of 100 attempts) which remained bounded. According to these measures both transient chaos and stable node dynamics are in agreement for two out of the three measures and differ marginally for the third.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content