# Publications

We demonstrate the utility of the equation-free methodology developed by one of the authors (IGK) for the study of scalar conservation laws with disordered initial conditions. The numerical scheme is benchmarked on exact solutions in Burgers turbulence corresponding to Levy process initial data. For these initial data, the kinetics of shock clustering is described by Smoluchowski's coagulation equation with additive kernel. The equation-free methodology is used to develop a particle scheme that computes self-similar solutions to the coagulation equation, including those with fat tails.

The exploration of epidemic dynamics on dynamically evolving ("adaptive") networks poses nontrivial challenges to the modeler, such as the determination of a small number of informative statistics of the detailed network state (that is, a few "good observables") that usefully summarize the overall (macroscopic, systems-level) behavior. Obtaining reduced, small size accurate models in terms of these few statistical observables - that is, trying to coarse-grain the full network epidemic model to a small but useful macroscopic one - is even more daunting. Here we describe a data-based approach to solving the first challenge: the detection of a few informative collective observables of the detailed epidemic dynamics. This is accomplished through Diffusion Maps (DMAPS), a recently developed data-mining technique. We illustrate the approach through simulations of a simple mathematical model of epidemics on a network: a model known to exhibit complex temporal dynamics. We discuss potential extensions of the approach, as well as possible shortcomings.

Spike-timing-dependent plasticity is the process by which the strengths of connections between neurons are modified as a result of the precise timing of the action potentials fired by the neurons. We consider a model consisting of one integrate-and-fire neuron receiving excitatory inputs from a large number-here, 1000-of Poisson neurons whose synapses are plastic. When correlations are introduced between the firing times of these input neurons, the distribution of synaptic strengths shows interesting, and apparently low-dimensional, dynamical behaviour. This behaviour is analysed in two different parameter regimes using equation-free techniques, which bypass the explicit derivation of the relevant low-dimensional dynamical system. We demonstrate both coarse projective integration (which speeds up the time integration of a dynamical system) and the use of recently developed data mining techniques to identify the appropriate low-dimensional description of the complex dynamical systems in our model.

We consider a harmonically driven acoustic medium in the form of a (finite length) highly nonlinear granular crystal with an amplitude-and frequency-dependent boundary drive. Despite the absence of a linear spectrum in the system, we identify resonant periodic propagation whereby the crystal responds at integer multiples of the drive period and observe that this can lead to local maxima of transmitted force at its fixed boundary. In addition, we identify and discuss minima of the transmitted force ("antiresonances") between these resonances. Representative one-parameter complex bifurcation diagrams involve period doublings and Neimark-Sacker bifurcations as well as multiple isolas (e.g., of period-3, -4, or -5 solutions entrained by the forcing). We combine them in a more detailed, two-parameter bifurcation diagram describing the stability of such responses to both frequency and amplitude variations of the drive. This picture supports a notion of a (purely) "nonlinear spectrum" in a system which allows no sound wave propagation (due to zero sound speed: the so-called sonic vacuum). We rationalize this behavior in terms of purely nonlinear building blocks: apparent traveling and standing nonlinear waves.

It is anticipated that in future generations of massively parallel computer systems a significant portion of processors may suffer from hardware or software faults rendering large-scale computations useless. In this work we address this problem from the algorithmic side, proposing resilient algorithms that can recover from such faults irrespective of their fault origin. In particular, we set the foundations of a new class of algorithms that will combine numerical approximations with machine learning methods. To this end, we consider three types of fault scenarios: (1) a gappy region but with no previous gaps and no contamination of surrounding simulation data, (2) a space-time gappy region but with full spatiotemporal information and no contamination, and (3) previous gaps with contamination of surrounding data. To recover from such faults we employ different reconstruction and simulation methods, namely the projective integration, the co-Kriging interpolation, and the resimulation method. In order to compare the effectiveness of these methods for the different processor faults and to quantify the error propagation in each case, we perform simulations of two benchmark flows, flow in a cavity and flow past a circular cylinder. In general, the projective integration seems to be the most effective method when the time gaps are small, and the resimulation method is the best when the time gaps are big while the co-Kriging method is independent of time gaps. Furthermore, the projective integration method and the co-Kriging method are found to be good estimation methods for the initial and boundary conditions of the resimulation method in scenario (3).

*N*neurons where this heterogeneity is characterized by a prescribed random variable, the oscillatory single-cluster state can transition—through

*N -*1 (possibly perturbed) period-doubling and subsequent bifurcations—to a variety of multiple-cluster states. The clustering dynamic behavior is computationally studied both at the detailed and the coarse-grained levels, and a numerical approach that can enable studying the coarse-grained dynamics in a network of arbitrarily large size is suggested. Among a number of cluster states formed, double clusters, composed of nearly equal sub-network sizes are seen to be stable; interestingly, the heterogeneity parameter in each of the double-cluster components tends to be consistent with the random variable over the entire network: Given a double-cluster state, permuting the dynamical variables of the neurons can lead to a combinatorially large number of different, yet similar “fine” states that appear practically identical at the coarse-grained level. For weak heterogeneity we find that correlations rapidly develop, within each cluster, between the neuron’s “identity” (its own value of the heterogeneity parameter) and its dynamical state. For single- and double-cluster states we demonstrate an effective coarse-graining approach that uses the Polynomial Chaos expansion to succinctly describe the dynamics by these quickly established “identity-state” correlations. This coarse-graining approach is utilized, within the equation-free framework, to perform efficient computations of the neuron ensemble dynamics.

The Koopman operator is a *linear* but infinite-dimensional operator that governs the evolution of scalar observables defined on the state space of an autonomous dynamical system and is a powerful tool for the analysis and decomposition of nonlinear dynamical systems. In this manuscript, we present a data-driven method for approximating the leading *eigenvalues, eigenfunctions, and modes* of the Koopman operator. The method requires a data set of snapshot pairs and a dictionary of scalar observables, but does not require explicit governing equations or interaction with a “black box” integrator. We will show that this approach is, in effect, an extension of dynamic mode decomposition (DMD), which has been used to approximate the Koopman eigenvalues and modes. Furthermore, if the data provided to the method are generated by a Markov process instead of a deterministic dynamical system, the algorithm approximates the eigenfunctions of the Kolmogorov backward equation, which could be considered as the “stochastic Koopman operator” (Mezic in Nonlinear Dynamics 41(1–3): 309–325, 2005). Finally, four illustrative examples are presented: two that highlight the quantitative performance of the method when presented with either deterministic or stochastic data and two that show potential applications of the Koopman eigenfunctions.

In recent years, individual-based/agent-based modeling has been applied to study a wide range of applications, ranging from engineering problems to phenomena in sociology, economics and biology. Simulating such agent-based models over extended spatiotemporal domains can be prohibitively expensive due to stochasticity and the presence of multiple scales. Nevertheless, many agent-based problems exhibit smooth behavior in space and time on a macroscopic scale, suggesting that a useful coarse-grained continuum model could be obtained. For such problems, the equation-free framework [16-18] can significantly reduce the computational cost. Patch dynamics is an essential component of this framework. This scheme is designed to perform numerical simulations of an unavailable macroscopic equation on macroscopic time and length scales; it uses appropriately initialized simulations of the fine-scale agent-based model in a number of small "patches", which cover only a fraction of the spatiotemporal domain. In this work, we construct a finite-volume-inspired conservative patch dynamics scheme and apply it to a financial market agent-based model based on the work of Omurtag and Sirovich [22]. We first apply our patch dynamics scheme to a continuum approximation of the agent-based model, to study its performance and analyze its accuracy. We then apply the scheme to the agent-based model itself. Our computational experiments indicate that here, typically, the patch dynamics-based simulation needs to be performed in only 20% of the full agent simulation space, and in only 10% of the temporal domain. (C) 2015 IMACS. Published by Elsevier B.V. All rights reserved.

We demonstrate that the Koopman eigenfunctions and eigenvalues define a set of intrinsic coordinates, which serve as a natural framework for fusing measurements obtained from heterogeneous collections of sensors in systems governed by nonlinear evolution laws. These measurements can be nonlinear, but must, in principle, be rich enough to allow the state to be reconstructed. We approximate the associated Koopman operator using extended dynamic mode decomposition, so the method only requires time series of data for each set of measurements, and a single set of "joint" measurements, which are known to correspond to the same underlying state. We apply this procedure to the FitzHugh-Nagumo PDE, and fuse measurements taken at a single point with principal-component measurements. Copyright (C) EPLA, 2015

Transient activation of the highly conserved extracellular-signal-regulated kinase (ERK) establishes precise patterns of cell fates in developing tissues. Quantitative parameters of these transients are essentially unknown, but a growing number of studies suggest that changes in these parameters can lead to a broad spectrum of developmental abnormalities. We provide a detailed quantitative picture of an ERK-dependent inductive signaling event in the early Drosophila embryo, an experimental system that offers unique opportunities for high-throughput studies of developmental signaling. Our analysis reveals a spatiotemporal pulse of ERK activation that is consistent with a model in which transient production of a short-ranged ligand feeds into a simple signal interpretation system. The pulse of ERK signaling acts as a switch in controlling the expression of the ERK target gene. The quantitative approach that led to this model, based on the integration of data from fixed embryos and live imaging, can be extended to other developmental systems patterned by transient inductive signals.

We study the coarse-grained, reduced dynamics of an agent-based market model due to Omurtag and Sirovich [18]. We first describe the large agent number, deterministic limit of the system dynamics by performing numerical bifurcation calculations on a continuum approximation of their model. By exploring a broad parameter space, we observe several interesting phenomena including turning points leading to unstable stationary agent density distributions as well as a type of "termination point." Close to these deterministic turning points we expect the stochastic underlying model to exhibit rare event transitions. We then proceed to discuss a coarse-grained approach to the quantitative study of these rare events. The basic assumption is that the dynamics of the system can be decomposed into fast (noise) and slow (single reaction coordinate) dynamics, so that the system can be described by an effective, coarse-grained Fokker-Planck(FP) equation. An explicit form of this effective FP equation is not available; in our computations we bypass the lack of a closed form equation by numerically estimating its components - the drift and diffusion coefficients - from ensembles of short bursts of microscopic simulations with judiciously chosen initial conditions. The reaction coordinate is first constructed based on our understanding of the continuum model close to the turning points, and it gives results reasonably close to those from brute-force direct simulations. When no guidelines for the selection of a good reaction coordinate are available, data-mining tools, in particular Diffusion Maps, can be used to determine a suitable reaction coordinate. In the third part of this work we demonstrate this "variable-free" approach by constructing a reaction coordinate simply based on the data from the simulation itself. This Diffusion Map based, empirical coordinate gives results consistent with the direct simulation.

Epithelial tissue, in which cells adhere tightly to each other and to the underlying substrate, is one of the four major tissue types in adult organisms. In embryos, epithelial sheets serve as versatile substrates during the formation of developing organs. Some aspects of epithelial morphogenesis can be adequately described using vertex models, in which the two-dimensional arrangement of epithelial cells is approximated by a polygonal lattice with an energy that has contributions reflecting the properties of individual cells and their interactions. Previous studies with such models have largely focused on dynamics confined to two spatial dimensions and analyzed them numerically. We show how these models can be extended to account for three-dimensional deformations and studied analytically. Starting from the extended model, we derive a continuum plate description of cell sheets, in which the effective tissue properties, such as bending rigidity, are related explicitly to the parameters of the vertex model. To derive the continuum plate model, we duly take into account a microscopic shift between the two sub-lattices of the hexagonal network, which has been ignored in previous work. As an application of the continuum model, we analyze tissue buckling by a line tension applied along a circular contour, a simplified set-up relevant to several situations in the developmental contexts. The buckling thresholds predicted by the continuum description are in good agreement with the results of stability calculations based on the vertex model. Our results establish a direct connection between discrete and continuum descriptions of cell sheets and can be used to probe a wide range of morphogenetic processes in epithelial tissues.

Understanding the mechanisms by which proteins fold from disordered amino-acid chains to spatially ordered structures remains an area of active inquiry. Molecular simulations can provide atomistic details of the folding dynamics which complement experimental findings. Conventional order parameters, such as root-mean-square deviation and radius of gyration, provide structural information but fail to capture the underlying dynamics of the protein folding process. It is therefore advantageous to adopt a method that can systematically analyze simulation data to extract relevant structural as well as dynamical information. The nonlinear dimensionality reduction technique known as diffusion maps automatically embeds the high-dimensional folding trajectories in a lower-dimensional space from which one can more easily visualize folding pathways, assuming the data lie approximately on a lower-dimensional manifold. The eigenvectors that parametrize the low-dimensional space, furthermore, are determined systematically, rather than chosen heuristically, as is done with phenomenological order parameters. We demonstrate that diffusion maps can effectively characterize the folding process of a Trp-cage miniprotein. By embedding molecular dynamics simulation trajectories of Trp-cage folding in diffusion maps space, we identify two folding pathways and intermediate structures that are consistent with the previous studies, demonstrating that this technique can be employed as an effective way of analyzing and constructing protein folding pathways from molecular simulations. (C) 2015 AIP Publishing LLC.

Progress of development is commonly reconstructed from imaging snapshots of chemical or mechanical processes in fixed tissues. As a first step in these reconstructions, snapshots must be spatially registered and ordered in time. Currently, image registration and ordering are often done manually, requiring a significant amount of expertise with a specific system. However, as the sizes of imaging data sets grow, these tasks become increasingly difficult, especially when the images are noisy and the developmental changes being examined are subtle. To address these challenges, we present an automated approach to simultaneously register and temporally order imaging data sets. The approach is based on vector diffusion maps, a manifold learning technique that does not require a priori knowledge of image features or a parametric model of the developmental dynamics. We illustrate this approach by registering and ordering data from imaging studies of pattern formation and morphogenesis in three model systems. We also provide software to aid in the application of our methodology to other experimental data sets.

As microscopic (e.g. atomistic, stochastic, agent-based, particle-based) simulations become increasingly prevalent in the modeling of complex systems, so does the need to systematically coarse-grain the information they provide. Before even starting to formulate relevant coarse-grained equations, we need to determine the right macroscopic *observables*—the right variables in terms of which emergent behavior will be described. This paper illustrates the use of data mining (and, in particular, diffusion maps, a nonlinear manifold learning technique) in coarse-graining the dynamics of a particle-based model of animal swarming. Our computational data-driven coarse-graining approach extracts two coarse (collective) variables from the detailed particle-based simulations, and helps formulate a low-dimensional stochastic differential equation in terms of these two collective variables; this allows the efficient quantification of the interplay of “informed” and “naive” individuals in the collective swarm dynamics. We also present a brief exploration of swarm breakup and use data-mining in an attempt to identify useful predictors for it. In our discussion of the scope and limitations of the approach we focus on the key step of selecting an informative metric, allowing us to usefully compare different particle swarm configurations.

In this paper we present a novel approach to the computation of convex invariant sets of dynamical systems. Employing a Banach space formalism to describe differences of convex compact subsets of Rn by directed sets, we are able to formulate the property of a convex, compact set to be invariant as a zero-finding problem in this Banach space. We need either the additional restrictive assumption that the image of sets from a subclass of convex compact sets under the dynamics remains convex, or we have to convexify these images. In both cases we can apply Newton's method in Banach spaces to approximate such invariant sets if an appropriate smoothness of a set-valued map holds. The theoretical foundations for realizing this approach are analyzed, and it is illustrated first by analytical and then by numerical examples.

We propose and illustrate an approach to coarse-graining the dynamics of evolving networks, i.e., networks whose connectivity changes dynamically. The approach is based on the equation-free framework: short bursts of detailed network evolution simulations are coupled with lifting and restriction operators that translate between actual network realizations and their appropriately chosen coarse observables. This framework is used here to accelerate temporal simulations through coarse projective integration, and to implement coarse-grained fixed point algorithms through matrix-free Newton-Krylov. The approach is illustrated through a very simple network evolution example, for which analytical approximations to the coarse-grained dynamics can be independently obtained, so as to validate the computational results. The scope and applicability of the approach, as well as the issue of selection of good coarse observables are discussed.

The adoption of detailed mechanisms for chemical kinetics often poses two types of severe challenges: First, the number of degrees of freedom is large; and second, the dynamics is characterized by widely disparate time scales. As a result, reactive flow solvers with detailed chemistry often become intractable even for large clusters of CPUs, especially when dealing with direct numerical simulation (DNS) of turbulent combustion problems. This has motivated the development of several techniques for reducing the complexity of such kinetics models, where, eventually, only a few variables are considered in the development of the simplified model. Unfortunately, no generally applicable a priori recipe for selecting suitable parameterizations of the reduced model is available, and the choice of slow variables often relies upon intuition and experience. We present an automated approach to this task, consisting of three main steps. First, the low dimensional manifold of slow motions is (approximately) sampled by brief simulations of the detailed model, starting from a rich enough ensemble of admissible initial conditions. Second, a global parametrization of the manifold is obtained through the Diffusion Map (DMAP) approach, which has recently emerged as a powerful tool in data analysis/machine learning. Finally, a simplified model is constructed and solved on the fly in terms of the above reduced (slow) variables. Clearly, closing this latter model requires nontrivial interpolation calculations, enabling restriction

(mapping from the full ambient space to the reduced one) and lifting (mapping from the reduced space to the ambient one). This is a key step in our approach, and a variety of interpolation schemes are reported and compared. The scope of the proposed procedure is presented and discussed by means of an illustrative combustion example.

Interacting particle systems constitute the dynamic model of choice in a variety of application areas. A prominent example is pedestrian dynamics, where good design of escape routes for large buildings and public areas can improve evacuation in emergency situations, avoiding exit blocking and the ensuing panic. Here we employ diffusion maps to study the coarse-grained dynamics of two pedestrian crowds trying to pass through a door from opposite sides. These macroscopic variables and the associated smooth embeddings lead to a better description and a clearer understanding of the nature of the transition to oscillatory dynamics. We also compare the results to those obtained through intuitively chosen macroscopic variables.

By applying an out-of-phase actuation at the boundaries of a uniform chain of granular particles, we demonstrate experimentally that time-periodic and spatially localized structures with a nonzero background (so-called dark breathers) emerge for a wide range of parameter values and initial conditions. We demonstrate a remarkable control over the number of breathers within the multibreather pattern that can be "dialed in" by varying the frequency or amplitude of the actuation. The values of the frequency (or amplitude) where the transition between different multibreather states occurs are predicted accurately by the proposed theoretical model, which is numerically shown to support exact dark breather and multibreather solutions. Moreover, we visualize detailed temporal and spatial profiles of breathers and, especially, of multibreathers using a full-field probing technology and enable a systematic favorable comparison among theory, computation, and experiments. A detailed bifurcation analysis reveals that the dark and multibreather families are connected in a "snaking" pattern, providing a roadmap for the identification of such fundamental states and their bistability in the laboratory.

Using the helix-coil transitions of alanine pentapeptide as an illustrative example, we demonstrate the use of diffusion maps in the analysis of molecular dynamics simulation trajectories. Diffusion maps and other nonlinear data-mining techniques provide powerful tools to visualize the distribution of structures in conformation space. The resulting low-dimensional representations help in partitioning conformation space, and in constructing Markov state models that capture the conformational dynamics. In an initial step, we use diffusion maps to reduce the dimensionality of the conformational dynamics of Ala5. The resulting pretreated data are then used in a clustering step. The identified clusters show excellent overlap with clusters obtained previously by using the backbone dihedral angles as input, with small-but nontrivial-differences reflecting torsional degrees of freedom ignored in the earlier approach. We then construct a Markov state model describing the conformational dynamics in terms of a discrete-time random walk between the clusters. We show that by combining fuzzy C-means clustering with a transition-based assignment of states, we can construct robust Markov state models. This state-assignment procedure suppresses short-time memory effects that result from the non-Markovianity of the dynamics projected onto the space of clusters. In a comparison with previous work, we demonstrate how manifold learning techniques may complement and enhance informed intuition commonly used to construct reduced descriptions of the dynamics in molecular conformation space. (C) 2014 AIP Publishing LLC.

We present a hybrid numerical approach for modeling surface reactions in the framework of a lattice-gas model with lateral interactions between adsorbed particles. A hybrid multiscale algorithm, which we refer to as Quasi-Equilibrium Kinetic Monte Carlo (QE-KMC), comprises traditional Metropolis Monte Carlo (MMC) simulations of equilibrium systems and standard numerical methods for deterministic ordinary differential equations (ODEs). The functional dependence of these ODEs on the macroscopic state variables (adsorbate coverages) is not explicitly known, but their right-hand sides can be evaluated "on the fly" with prescribed accuracy by means of the MMC simulations. At the time scale of these ODEs it is assumed that an equilibrium statistical distribution of adsorbed particles on an infinite lattice is attained at every moment in time due to infinitely fast surface diffusion. QE-KMC and conventional KMC simulations are used to study the temperature-programmed desorption (TPD) spectra of adsorbed particles. We critically discuss results of previous studies that applied Monte Carlo simulations to describe the TPD spectra in the case of fast adsorbate diffusion and strong lateral interactions. We show that the quasi-equilibrium TPD spectra can be quickly and accurately estimated by the QE-KMC algorithm, while the KMC simulations require much more extensive computational resources to obtain the same results. (C) 2013 Elsevier Ltd. All rights reserved.

We introduce radial basis functions (RBFs) whose time-varying coefficients determine not only the amplitude and position of each RBF but also their shape. The intended use of these Time Varying-RBFs (TV-RBFs) is in the local-in-time representation of low-dimensional approximations of functions that arise in solving spatiotemporal evolution problems; in particular, for time-varying spatially localized solutions with a temporal translation component such as traveling waves, modulated pulses or soliton-like solutions of evolutionary differential equations. This paper is restricted to the one-dimensional spatial case. We also present an algorithm that places the Time Varying-RBFs (TV-RBFs) over spatiotemporal data that may come from experiments, from finely discretized PDE simulations, or even from multiscale, particle-based simulations. It first approximates the function at a single time instant (a temporal snapshot) as a sum of RBFs using a novel weighted minimization that causes each RBF to primarily approximate one of the localized parts of the function). It then extends that approximation to TV-RBFs over a sequence of snapshots of the function at different times. We conclude by discussing the potential uses of these TV-RBFs. (C) 2014 Elsevier B.V. All rights reserved.

A macroscopic first-principles mathematical model of a solid oxide fuel cell (SOFC) with a BaCe1-xSmxO3-a type electrolyte is developed. This class of electrolytes exhibits both proton and oxygen-anion conductivity. The existence of steady-state multiplicity with respect to bifurcation parameters such as cell outlet voltage, cell current, ohmic external load, and cell power density is investigated. The analyses with respect to the first two bifurcation parameters represent potentiostatic and galvanostatic modes of operation, respectively. The cell can have up to three steady states with respect to the external load resistance and cell outlet voltage, and a unique steady state with respect to the cell current. The cell can have four steady states with respect to cell power density (when power density is the bifurcation parameter). The same qualitative steady-state behavior had been observed in oxygen-ion conducting and proton conducting SOFCs. Furthermore, thermal and concentration multiplicities coexist in this class of SOFCs; ignition in the solid temperature is accompanied by (a) extinction in the fuel and oxygen concentrations, and (b) ignition and extinction in concentrations of water in the anode and cathode sides, respectively. The cell exhibits steady-state multiplicity with respect to the fuel and air flow rates under fixed ohmic load, potentiostatic, and fixed cell power density modes of operation, but does not show this behavior under a galvanostatic mode of operation. The dynamic response of the cell also shows that the cell is highly nonlinear and some of the cell variables can show inverse response.

Reconciling competing desires to build urban models that can be simple and complicated is something of a grand challenge for urban simulation. It also prompts difficulties in many urban policy situations, such as urban sprawl, where simple, actionable ideas may need to be considered in the context of the messily complex and complicated urban processes and phenomena that work within cities. In this paper, we present a novel architecture for achieving both simple and complicated realizations of urban sprawl in simulation. Fine-scale simulations of sprawl geography are run using geographic automata to represent the geographical drivers of sprawl in intricate detail and over fine resolutions of space and time. We use Equation-Free computing to deploy population as a coarse observable of sprawl, which can be leveraged to run automata-based models as short-burst experiments within a meta-simulation framework.

We consider a statically compressed diatomic granular crystal, consisting of alternating aluminum and steel spheres. The combination of dissipation, driving of the boundary, and intrinsic nonlinearity leads to complex dynamics. Through both numerical simulations and experiments, we find that the interplay of nonlinear surface modes with modes caused by the driver create the possibility, as the driving amplitude is increased, of limit cycle saddle-node bifurcations beyond which the dynamics of the system becomes chaotic. In this chaotic state, part of the applied energy can propagate through the chain. We also find that the chaotic branch depends weakly on the driving frequency, and speculate a connection between the chaotic dynamics with the gap openings between the spheres. Finally, we observe hysteretic dynamics and an interval of multi-stability involving stable periodic solutions and chaotic ones. Copyright (C) EPLA, 2013

We present a model based on the lattice Boltzmann equation that is suitable for the simulation of dynamic wetting. The model is capable of exhibiting fundamental interfacial phenomena such as weak adsorption of fluid on the solid substrate and the presence of a thin surface film within which a disjoining pressure acts. Dynamics in this surface film, tightly coupled with hydrodynamics in the fluid bulk, determine macroscopic properties of primary interest: the hydrodynamic slip; the equilibrium contact angle; and the static and dynamic hysteresis of the contact angles. The pseudo-potentials employed for fluid-solid interactions are composed of a repulsive core and an attractive tail that can be independently adjusted. This enables effective modification of the functional form of the disjoining pressure so that one can vary the static and dynamic hysteresis on surfaces that exhibit the same equilibrium contact angle. The modeled fluid-solid interface is diffuse, represented by a wall probability function that ultimately controls the momentum exchange between solid and fluid phases. This approach allows us to effectively vary the slip length for a given wettability (i.e., a given static contact angle) of the solid substrate. DOI: 10.1103/PhysRevE.87.013302

A macroscopic first-principles mathematical model of a solid oxide fuel cell (SOFC) with a BaCe1-xSmxO3-a type electrolyte is developed. This class of electrolytes exhibits both proton and oxygen-anion conductivity. The existence of steady-state multiplicity with respect to bifurcation parameters such as cell outlet voltage, cell current, ohmic external load, and cell power density is investigated. The analyses with respect to the first two bifurcation parameters represent potentiostatic and galvanostatic modes of operation, respectively. The cell can have up to three steady states with respect to the external load resistance and cell outlet voltage, and a unique steady state with respect to the cell current. The cell can have four steady states with respect to cell power density (when power density is the bifurcation parameter). The same qualitative steady-state behavior had been observed in oxygen-ion conducting and proton conducting SOFCs. Furthermore, thermal and concentration multiplicities coexist in this class of SOFCs; ignition in the solid temperature is accompanied by (a) extinction in the fuel and oxygen concentrations, and (b) ignition and extinction in concentrations of water in the anode and cathode sides, respectively. The cell exhibits steady-state multiplicity with respect to the fuel and air flow rates under fixed ohmic load, potentiostatic, and fixed cell power density modes of operation, but does not show this behavior under a galvanostatic mode of operation. The dynamic response of the cell also shows that the cell is highly nonlinear and some of the cell variables can show inverse response.

We process snapshots of trajectories of evolution equations with intrinsic symmetries, and demonstrate the use of recently developed eigenvector-based techniques to successfully quotient out the degrees of freedom associated with the symmetries in the presence of noise. Our illustrative examples include a one-dimensional evolutionary partial differential (the Kuramoto-Sivashinsky) equation with periodic boundary conditions, as well as a stochastic simulation of nematic liquid crystals which can be effectively modeled through a nonlinear Smoluchowski equation on the surface of a sphere. This is a useful first step towards data mining the symmetry-adjusted ensemble of snapshots in search of an accurate low-dimensional parametrization and the associated reduction of the original dynamical system. We also demonstrate a technique (Vector Diffusion Maps) that combines, in a single formulation, the symmetry removal step and the dimensionality reduction step. (C) 2013 Elsevier Ltd. All rights reserved.

Finding informative low-dimensional descriptions of high-dimensional simulation data (like the ones arising in molecular dynamics or kinetic Monte Carlo simulations of physical and chemical processes) is crucial to understanding physical phenomena, and can also dramatically assist in accelerating the simulations themselves. In this paper, we discuss and illustrate the use of nonlinear intrinsic variables (NIV) in the mining of high-dimensional multiscale simulation data. In particular, we focus on the way NIV allows us to functionally merge different simulation ensembles, and different partial observations of these ensembles, as well as to infer variables not explicitly measured. The approach relies on certain simple features of the underlying process variability to filter out measurement noise and systematically recover a unique reference coordinate frame. We illustrate the approach through two distinct sets of atomistic simulations: a stochastic simulation of an enzyme reaction network exhibiting both fast and slow time scales, and a molecular dynamics simulation of alanine dipeptide in explicit water. (C) 2013 AIP Publishing LLC.

A numerical approach for the approximation of inertial manifolds of stochastic evolutionary equations with multiplicative noise is presented and illustrated. After splitting the stochastic evolutionary equations into a backward and a forward part, a numerical scheme is devised for solving this backward-forward stochastic system, and an ensemble of graphs representing the inertial manifold is consequently obtained. This numerical approach is tested in two illustrative examples: one is for a system of stochastic ordinary differential equations and the other is for a stochastic partial differential equation.

Model reduction is an important systems task with a long history in traditional chemical engineering modeling. We discuss its interplay with modern data-mining tools (such as Local Feature Analysis and Diffusion Maps) through illustrative examples, and comment on important open issues regarding applications to large systems arising in molecular/atomistic simulations. (C) 2012 Elsevier Ltd. All rights reserved.

Given a legacy dynamic simulator of a chemical process plant, we construct a computational procedure that can be wrapped around the simulator to compute its steady states (both stable and unstable) and their dependence on input parameters. We apply this approach to the Tennessee Eastman (TE) challenge problem presented by Downs and Vogel, who also provided a FORTRAN process model. Using the FORTRAN simulator as a black-box input-output map, we enable it to systematically converge to isolated solutions and study their stability and parametric dependence within the equation-free framework. The presence of neutrally stable modes in TE problem (due to so-called inventories), their interplay with the problem formulation and the convergence of the solution procedure is explored and rationalized. Interestingly, our time-stepper formulation can automatically take advantage of separation of time scales, when present, to enhance computational convergence. The approach enables legacy dynamic simulators to calculate several dynamic problem characteristics useful for controller design and/or process optimization. (c) 2013 American Institute of Chemical Engineers AIChE J, 59: 3308-3321, 2013

Water emerging from similar to 100 mu m pores into millimeter-size gas flow channels forms drops that grow and become slugs which span the flow channel. Flowing gas causes the slugs to detach and move down the channel. The effect of channel geometry, surface wettability, and gravity on the formation and motion of water slugs has been analyzed using high-speed video images of the drops and differential pressure-time traces. Drops grow and appear, assuming a sequence of shapes that minimize the total interfacial energy of the gas-liquid and liquid-solid interfaces. The drops are initially spherical caps centered on the pore (the liquid contacts one wall). Above a certain size, the drops move to the corner, forming "corner drops" (the liquid contacts two walls). Corner drops grow across the channel, evolving into partial liquid bridges (drops confined by three walls), and finally the drops span the channel cross-section forming slugs (contacting all four walls). Smaller slugs are formed in channels with hydrophobic walls than in channels with hydrophilic walls. Smaller slugs are formed in channels with curved walls than in square or rectangular channels. Slugs move when the differential gas pressure overcomes the force to move the advancing and receding gas-liquid-solid contact lines of the slugs. Residual water left behind in corners by moving slugs reduces the barriers for drops to form slugs, causing the steady-state slug volumes to be smaller than those seen at start-up in dry channels.

Water emerging from micrometer-sized pores into millimeter-sized gas-flow channels forms drops. The drops grow until the force from the flowing gas is sufficient to detach the drops as either (1) slugs that completely occlude the cross section of the channel and move at the superficial gas velocity, (2) drops that partially occlude the channel and move at a velocity that is less than the gas velocity, or (3) films that flow continuously, occluding part of the channel. At steady state, small residual water droplets, similar to 100 mu m in diameter, left in corners and on surface defects from previous drops, are key in determining the shape of water drops at detachment. Slugs are formed at low-gas-phase Reynolds numbers (Re-G) in both hydrophilic and hydrophobic channels. Drops are shed in Teflon-coated hydrophobic channels for Re-G > 30. Films are formed in acrylic hydrophilic channels for Re-G > 30. Slugs form when growing drops encounter residual water droplets that nucleate the drop to slug transition. Drops are shed when the force exerted by the flowing gas on growing drops exceeds the force needed to advance the gas/liquid/solid contact line before they grow to the critical size for the drop to slug transition. Drops grow by "stick-slip" of the solid-liquid-gas contact lines and with pinned contact lines until the force on the drops results in either the downstream contact angle becoming greater than the dynamic advancing contact angle or the upstream contact angle becoming less than the dynamic receding contact angle. The upstream contact line never detaches for hydrophilic channels, which is why films form. The shape of water drops and the detachment energies are shown to be well approximated by the force balance between the force needed to advance the drop's contact lines and the force that the flowing gas exerts on a stationary drop.

We consider a coupled, heterogeneous population of relaxation oscillators used to model rhythmic oscillations in the pre-Bötzinger complex. By choosing specific values of the parameter used to describe the heterogeneity, sampled from the probability distribution of the values of that parameter, we show how the effects of heterogeneity can be studied in a computationally efficient manner. When more than one parameter is heterogeneous, full or sparse tensor product grids are used to select appropriate parameter values. The method allows us to effectively reduce the dimensionality of the model, and it provides a means for systematically investigating the effects of heterogeneity in coupled systems, linking ideas from uncertainty quantification to those for the study of network dynamics.

Despite its popularity, agent-based modeling is limited by serious barriers that constrain its usefulness as an exploratory tool. In particular, there is a paucity of systematic approaches for extracting coarse-grained, system-level information as it emerges in direct simulation. This is particularly problematic for agent-based models (ABMs) of complex urban systems in which macroscopic phenomena, such as sprawl, may manifest themselves coarsely from bottom-up dynamics among diverse agent-actors interacting across scales. Often these connections are not known, but treating them is nevertheless crucial in enabling prediction, in supporting decisions, and in facilitating the design, control, and optimization of urban systems. In this article, we describe and implement a metasimulation scheme for extracting macroscopic information from local dynamics of agent-based simulation, which allows acceleration of coarse-scale computing and which may also serve as a precursor to handle emergence in complex urban simulation. We compare direct ABM simulation, population-level equation solutions, and coarse projective integration. We apply the scheme to the simulation of urban sprawl from local drivers of urbanization, urban growth, and population dynamics. Numerical examples of the three approaches are provided to compare their accuracy and efficiency. We find that our metasimulation scheme can significantly accelerate complex urban simulations while maintaining faithful representation of the original model.

We report on a technique that utilizes an array of galvanic microreactors to guide the assembly of two-dimensional colloidal crystals with spatial and orientational order. Our system is comprised of an array of copper and gold electrodes in a coplanar arrangement, immersed in a dilute hydrochloric acid solution in which colloidal micro-spheres of polystyrene and silica are suspended. Under optimized conditions, two-dimensional colloidal crystals form at the anodic copper with patterns and crystal orientation governed by the electrode geometry. After the aggregation process, the colloidal particles are cemented to the substrate by co-deposition of reaction products. As we vary the electrode geometry, the dissolution rate of the copper electrodes is altered. This way, we control the colloidal motion as well as the degree of reaction product formation. We show that particle motion is governed by a combination of electrokinetic effects acting directly on the colloidal particles and bulk electrolyte flow generated at the copper-gold interface. (C) 2012 American Institute of Physics. [http://dx.doi.org/10.1063/1.4755807]

We consider the simplest network of coupled non-identical phase oscillators capable of displaying a "chimera" state (namely, two subnetworks with strong coupling within the subnetworks and weaker coupling between them) and systematically investigate the effects of gradually removing connections within the network, in a random but systematically specified way. We average over ensembles of networks with the same random connectivity but different intrinsic oscillator frequencies and derive ordinary differential equations (ODEs), whose fixed points describe a typical chimera state in a representative network of phase oscillators. Following these fixed points as parameters are varied we find that chimera states are quite sensitive to such random removals of connections, and that oscillations of chimera states can be either created or suppressed in apparent bifurcation points, depending on exactly how the connections are gradually removed. (C) 2012 American Institute of Physics. [http://dx.doi.org/10.1063/1.3694118]

We explore a systematic approach to studying the dynamics of evolving networks at a coarse-grained, system level. We emphasize the importance of finding good observables (network properties) in terms of which coarse-grained models can be developed. We illustrate our approach through a particular social network model: the 'rise and fall' of a networked society (Marsili M et al 2004 Proc. Natl Acad. Sci. USA 101 1439). We implement our low-dimensional description computationally using the equation-free approach and show how it can be used to (i) accelerate simulations and (ii) extract system-level stability/bifurcation information from the detailed dynamic model. We discuss other system-level tasks that can be enabled through such a computer-assisted coarse-graining approach.

Liquid water is pushed through flow channels of fuel cells, where one surface is a porous carbon electrode made up of carbon fibers. Water drops grow on the fibrous carbon surface in the gas flow channel. The drops adhere to the superficial fiber surfaces but exhibit little penetration into the voids between the fibers. The fibrous surfaces are hydrophobic, but there is. a substantial threshold force necessary to initiate water drop, motion. Once the water drops, begin to move, however, the adhesive force decreases and drops move with minimal friction, similar to motion on superhydrophobic materials. We report here studies of water wetting and water drop motion on typical porous carbon materials (carbon paper and carbon cloth) employed in fuel cells. The static coefficient of friction on these textured surfaces is comparable to that for smooth Teflon. But the dynamic coefficient of-friction-is several orders of magnitude smaller on the textured surfaces than on smooth Teflon. Carbon cloth displays a much smaller static contact angle hysteresis than carbon paper due to its two-scale roughness. The dynamic contact angle hysteresis for carbon paper is greatly reduced compared to the static contact angle hysteresis. Enhanced dynamic hydrophobicity is suggested to result from the extent to which a dynamic contact line can track topological heterogeneities of the liquid/solid interface.

Cancer is a disease that can be seen as a complex system whose dynamics and growth result from nonlinear processes coupled across wide ranges of spatio-temporal scales. The current mathematical modeling literature addresses issues at various scales but the development of theoretical methodologies capable of bridging gaps across scales needs further study. We present a new theoretical framework based on Dynamic Density Functional Theory (DDFT) extended, for the first time, to the dynamics of living tissues by accounting for cell density correlations, different cell types, phenotypes and cell birth/death processes, in order to provide a biophysically consistent description of processes across the scales. We present an application of this approach to tumor growth. Copyright 2012 Author(s). This article is distributed under a Creative Commons Attribution 3.0 Unported License. [http://dx.doi.org/10.1063/1.3699065]

The subject of this work is the development and implementation of algorithms which accelerate the simulation of early stage tumor growth models. Among the different computational approaches used for the simulation of tumor progression, discrete stochastic models (e. g., cellular automata) have been widely used to describe processes occurring at the cell and subcell scales (e. g., cell-cell interactions and signaling processes). To describe macroscopic characteristics (e. g., morphology) of growing tumors, large numbers of interacting cells must be simulated. However, the high computational demands of stochastic models make the simulation of large-scale systems impractical. Alternatively, continuum models, which can describe behavior at the tumor scale, often rely on phenomenological assumptions in place of rigorous upscaling of microscopic models. This limits their predictive power. In this work, we circumvent the derivation of closed macroscopic equations for the growing cancer cell populations; instead, we construct, based on the so-called "equation-free" framework, a computational superstructure, which wraps around the individual-based cell-level simulator and accelerates the computations required for the study of the long-time behavior of systems involving many interacting cells. The microscopic model, e. g., a cellular automaton, which simulates the evolution of cancer cell populations, is executed for relatively short time intervals, at the end of which coarse-scale information is obtained. These coarse variables evolve on slower time scales than each individual cell in the population, enabling the application of forward projection schemes, which extrapolate their values at later times. This technique is referred to as coarse projective integration. Increasing the ratio of projection times to microscopic simulator execution times enhances the computational savings. Crucial accuracy issues arising for growing tumors with radial symmetry are addressed by applying the coarse projective integration scheme in a cotraveling (cogrowing) frame. As a proof of principle, we demonstrate that the application of this scheme yields highly accurate solutions, while preserving the computational savings of coarse projective integration.

We show how the equation-free approach can be exploited to enable agent-based simulators to perform system-level computations such as bifurcation, stability analysis and controller design. We illustrate these tasks through an event-driven agent-based model describing the dynamic behaviour of many interacting investors in the presence of mimesis. Using short bursts of appropriately initialized runs of the detailed, agent-based simulator, we construct the coarse-grained bifurcation diagram of the (expected) density of agents and investigate the stability of its multiple solution branches. When the mimetic coupling between agents becomes strong enough, the stable stationary state loses its stability at a coarse turning point bifurcation. We also demonstrate how the framework can be used to design a wash-out dynamic controller that stabilizes open-loop unstable stationary states even under model uncertainty. Copyright (C) EPLA, 2012

Kinetic Monte Carlo simulations are used to study the stochastic two-species Lotka-Volterra model on a square lattice. For certain values of the model parameters, the system constitutes an excitable medium: travelling pulses and rotating spiral waves can be excited. Stable solitary pulses travel with constant (modulo stochastic fluctuations) shape and speed along a periodic lattice. The spiral waves observed persist sometimes for hundreds of rotations, but they are ultimately unstable and break-up (because of fluctuations and interactions between neighboring fronts) giving rise to complex dynamic behavior in which numerous small spiral waves rotate and interact with each other. It is interesting that travelling pulses and spiral waves can be exhibited by the model even for completely immobile species, due to the non-local reaction kinetics. (C) 2012 American Institute of Physics. [http://dx.doi.org/10.1063/1.4729141]

A mathematical model of a solid oxide fuel cell (SOFC) with a BaCe1-xSmxO3-alpha type electrolyte is developed. This class of electrolytes exhibits both proton and oxygen-anion conductivity. To develop the model, heat transfer, mass transfer and electrochemical processes are taken into account. The existence of steady-state multiplicity in this class of fuel cells is investigated under three operation modes: constant ohmic load, potentiostatic and galavanostatic. The cell has up to three steady states under the constant ohmic load and potentiostatic modes, and a unique steady state under the galvanostatic mode. This same steady state behavior has been observed in oxygen-anion conducting and proton conducting SOFCs. Interestingly, this study shows that in this class of SOFCs, thermal and concentration multiplicities can coexist; ignition in the solid temperature is accompanied by extinction in the fuel and oxygen concentrations, and ignition and extinction in concentrations of water in the anode and cathode sides, respectively.

Micro-or nano-structurally roughened solid surfaces exhibit a rich variety of wetting behavior types, ranging from superhydro- or superoleophobicity to superhydro- or superoleophilicity. Depending on their material chemistry, the scale and morphology of their roughness or even the application of external electric fields, their apparent wettability can be significantly modified giving rise to challenging technological applications by exploiting the associated capillary phenomena at the micrometer scale. Certain applications, however, are limited by hysteretic wetting transitions, which inhibit spontaneous switching between wetting states, requiring external stimuli or actuation like thermal heating. The presence of surface roughness, necessary for the manifestation of the superhydrophobicity, induces multiplicity of wetting states and the inevitable hysteresis appears due to considerable energy barriers separating the equilibrium states. Here, by using continuum as well as mesoscopic computational analysis we perform a systems level study of the mechanisms of wetting transitions on model structured solid surfaces. By tracing entire equilibrium solution families and determining their relative stability we are able to illuminate mechanisms of wetting transitions and compute the corresponding energy barriers. The implementation of our analysis to 'real world' structured or unstructured surfaces is straightforward, rendering our computational tools valuable not only for the realization of surfaces with addressable wettability through roughness design, but also for the design of suitable actuation for optimal switching between wetting states.

Agent-based modeling (ABM) constitutes a powerful computational tool for the exploration of phenomena involving emergent dynamic behavior in the social sciences. This paper demonstrates a computer-assisted approach that bridges the significant gap between the single-agent microscopic level and the macroscopic (coarse-grained population) level, where fundamental questions must be rationally answered and policies guiding the emergent dynamics devised. Our approach will be illustrated through an agent-based model of civil violence. This spatiotemporally varying ABM incorporates interactions between a heterogeneous population of citizens [active (insurgent), inactive, or jailed] and a population of police officers. Detailed simulations exhibit an equilibrium punctuated by periods of social upheavals. We show how to effectively reduce the agent-based dynamics to a stochastic model with only two coarse-grained degrees of freedom: the number of jailed citizens and the number of active ones. The coarse-grained model captures the ABM dynamics while drastically reducing the computation time (by a factor of approximately 20).

In [C. W. Gear, T. J. Kaper, I. G. Kevrekidis and A. Zagaris, Projecting to a slow manifold: Singularly perturbed systems and legacy codes, SIAM J. Appl. Dyn. Syst. 4 (2005), 711-732], we developed the family of constrained runs algorithms to find points on low-dimensional, attracting, slow manifolds in systems of nonlinear differential equations with multiple time scales. For user-specified values of a subset of the system variables parametrizing the slow manifold (which we term observables and denote collectively by u), these iterative algorithms return values of the remaining system variables v so that the point (u, v) approximates a point on a slow manifold. In particular, the m-th constrained runs algorithm (m = 0,1, ...) approximates a point (u, v(m)) that is the appropriate zero of the (m + 1) -st time derivative of v. The accuracy with which (u, v(m)) approximates the corresponding point on the slow manifold with the same value of the observables has been established in [A. Zagaris, C. W. Gear, T. J. Kaper and I.G. Kevrekidis, Analysis of the accuracy and convergence of equation-free projection to a slow manifold, ESAIM: M2AN 43(4) (2009) 757-784] for systems for which the observables u evolve exclusively on the slow time scale. There, we also determined explicit conditions under which the m-th constrained runs scheme converges to the fixed point (u, v(m)) and identified conditions under which it fails to converge. Here, we consider the questions of stability and stabilization of these iterative algorithms for the case in which the observables u are also allowed to evolve on a fast time scale. The stability question in this case is more complicated, since it involves a generalized eigenvalue problem for a pair of matrices encoding geometric and dynamical characteristics of the system of differential equations. We determine the conditions under which these schemes converge or diverge in a series of cases in which this problem is explicitly solvable. We illustrate our main stability and stabilization results for the constrained runs schemes on certain planar systems with multiple time scales, and also on a more-realistic sixth order system with multiple time scales that models a network of coupled enzymatic reactions. Finally, we consider the issue of stabilization of the m-th constrained runs algorithm when the functional iteration scheme is divergent or converges slowly. In that case, we demonstrate on concrete examples how Newton's method and Broyden's method may replace functional iteration to yield stable iterative schemes.

Oxygen transport across the cathode gas diffusion layer (GDL) in polymer electrolyte membrane (PEM) fuel cells was examined by varying the O_{2}/N_{2} ratio and by varying the area of the GDL extending laterally from the gas flow channel under the bipolar plate (under the land). As the cathode is depleted of oxygen, the current density becomes limited by oxygen transport across the GDL. Oxygen depletion from O_{2}/N_{2} mixtures limits catalyst utilization, especially under the land.The local current density with air fed PEM fuel cells falls to practically zero at lateral distances under the land more than 3 times the GDL thickness; on the other hand, catalyst utilization was not limited when the fuel cell cathode was fed with 100% oxygen. The ratio of GDL thickness to the extent of the land is thus critical to the effective utilization of the catalyst in an air fed PEM fuel cell.

The long time behavior of the dynamics of an example of a fast-slow system of ordinary differential equations is examined. The system is derived from a spatial discretization of a Korteweg-de Vries-Burgers type equation, with fast dispersion and slow diffusion. The discretization is based on a model developed by Goodman and Lax, that is composed of a fast system drifted by a slow forcing term. A difficulty to invoke available multiscale methods arises since the underlying system does not possess a natural split to fast and slow state variables. Our approach depicts the limit behavior as a Young measure with values being invariant measures of the fast contribution to the flow. The slow contribution to the dynamics causes these invariant measures to drift. We keep track of this drift via slowly evolving observables. Averaging equations for the latter lead to computation of characteristic features of the motion and the location the invariant measures. Such computations are presented in the

We present a computer-assisted approach to coarse graining the evolutionary dynamics of a system of nonidentical oscillators coupled through a (fixed) network structure. The existence of a spectral gap for the coupling network graph Laplacian suggests that the graph dynamics may quickly become low dimensional. Our first choice of coarse variables consists of the components of the oscillator states-their (complex) phase angles-along the leading eigenvectors of this Laplacian. We then use the equation-free framework, circumventing the derivation of explicit coarse-grained equations, to perform computational tasks such as coarse projective integration, coarse fixed-point, and coarse limit-cycle computations. In a second step, we explore an approach to incorporating oscillator heterogeneity in the coarse-graining process. The approach is based on the observation of fast-developing correlations between oscillator state and oscillator intrinsic properties and establishes a connection with tools developed in the context of uncertainty quantification.

We have examined the influence of coarse-graining polymer chains in dissipative particle dynamics simulations on both phase behavior and aggregation dynamics. Our coarse-graining approach involves replacing several beads of a chain with a single bead of larger size and mass, and extends a framework recently developed by Backer et al. [J. Chem. Phys. 2005, 123, 114905] The parameters governing the interactions between particles are determined on the basis of conserving the number of interactions per particle, mass density, pressure, and shear viscosity of the original reference system. Phase diagrams of coarse-grained polymer/solvent systems are conserved, aside from a simple vertical shift in the direction of a(12), the repulsion parameter between solvent and polymer particles. The dynamics of the aggregation process are well conserved upon coarse-graining in a diblock copolymer/solvent system. We find that the invariance of the phase diagrams and aggregation dynamics occurs when the molecular volume of each species is conserved upon coarse-graining. However, this suggests that our coarse-graining approach cannot be applied to a monomeric solvent, in which case each solvent molecule already has the minimum number of degrees of freedom. Our results suggest that considerable freedom exists for selection of mapping ratios from real monomer units to model beads in dissipative particle dynamics simulations of chains in a solvent.

We explore a systematic approach to studying network evolutionary models at a coarse-grained, systems level. We emphasize the importance of finding good observables (network properties) in terms of which coarse-grained models can be developed. We illustrate our approach on a particular social network model: the "rise and fall" of a networked society; we implement the low-dimensional description computationally using the equation-free approach.

We present a simple technique for the computation of coarse-scale steady states of dynamical systems with time scale separation in the form of a "wrapper" around a fine-scale simulator. We discuss how this approach alleviates certain problems encountered by comparable existing approaches, and illustrate its use by computing coarse-scale steady states of a lattice Boltzmann fine scale code. Interestingly, in the same context of multiple time scale problems, the approach can be slightly modified to provide initial conditions on the slow manifold with prescribed coarse-scale observables. The approach is based on appropriately designed short bursts of the fine-scale simulator whose results are used to track changes in the coarse variables of interest, a core component of the equation-free framework. (C) 2011 Published by Elsevier Ltd.

We have developed explicit- and implicit-solvent models for the flash nanoprecipitation process, which involves rapid coprecipitation of block copolymers and solutes by changing solvent quality. The explicit-solvent model uses the dissipative particle dynamics (DPD) method and the implicit-solvent model uses the Brownian dynamics (BD) method. Each of the two models was parameterized to match key properties of the diblock copolymer (specifically, critical micelle concentration, diffusion coefficient, polystyrene melt density, and polyethylene glycol radius of gyration) and the hydrophobic solute (aqueous solubility, diffusion coefficient, and solid density). The models were simulated in the limit of instantaneous mixing of solvent with antisolvent. Despite the significant differences in the potentials employed in the implicit- and explicit-solvent models, the polymer-stabilized nanoparticles formed in both sets of simulations are similar in size and structure; however, the dynamic evolution of the two simulations is quite different. Nanoparticles in the BD simulations have diffusion coefficients that follow Rouse behavior (D alpha M-1), whereas those in the DPD simulations have diffusion coefficients that are close to the values predicted by the Stokes-Einstein relation (D alpha R-1). As the nanoparticles become larger, the discrepancy between diffusion coefficients grows. As a consequence, BD simulations produce increasingly slower aggregation dynamics with respect to real time and result in an unphysical evolution of the nanoparticle size distribution. Surface area per polymer of the stable explicit-solvent nanoparticles agrees well with experimental values, whereas the implicit-solvent nanoparticles are stable when the surface area per particle is roughly two to four times larger. We conclude that implicit-solvent models may produce questionable results when simulating nonequilibrium processes in which hydrodynamics play a critical role. (C) 2011 American Institute of Physics. [doi:10.1063/1.3580293]

Stochastic simulation of coupled chemical reactions is often computationally intensive, especially if a chemical system contains reactions occurring on different time scales. In this paper, we introduce a multiscale methodology suitable to address this problem, assuming that the evolution of the slow species in the system is well approximated by a Langevin process. It is based on the conditional stochastic simulation algorithm (CSSA) which samples from the conditional distribution of the suitably defined fast variables, given values for the slow variables. In the constrained multiscale algorithm (CMA) a single realization of the CSSA is then used for each value of the slow variable to approximate the effective drift and diffusion terms, in a similar manner to the constrained mean-force computations in other applications such as molecular dynamics. We then show how using the ensuing Fokker-Planck equation approximation, we can in turn approximate average switching times in stochastic chemical systems. (C) 2011 American Institute of Physics. [doi:10.1063/1.3624333]

Dissipative particle dynamics simulations were used to study the effects of mixing time, solute solubility, solute and diblock copolymer concentrations, and copolymer block length on the rapid coprecipitation of polymer-protected nanoparticles. The simulations were aimed at modeling Flash NanoPrecipitation, a process in which hydrophobic solutes and amphiphilic block copolymers are dissolved in a water-miscible organic solvent and then rapidly mixed with water to produce composite nanoparticles. A previously developed model by Spaeth [J. Chem. Phys. 134, 164902 (2011)] was used. The model was parameterized to reproduce equilibrium and transport properties of the solvent, hydrophobic solute, and diblock copolymer. Anti-solvent mixing was modeled using time-dependent solvent-solute and solvent-copolymer interactions. We find that particle size increases with mixing time, due to the difference in solute and polymer solubilities. Increasing the solubility of the solute leads to larger nanoparticles for unfavorable solute-polymer interactions and to smaller nanoparticles for favorable solute-polymer interactions. A decrease in overall solute and polymer concentration produces smaller nanoparticles, because the difference in the diffusion coefficients of a single polymer and of larger clusters becomes more important to their relative rates of collisions under more dilute conditions. An increase in the solute-polymer ratio produces larger nanoparticles, since a collection of large particles has less surface area than a collection of small particles with the same total volume. An increase in the hydrophilic block length of the polymer leads to smaller nanoparticles, due to an enhanced ability of each polymer to shield the nanoparticle core. For unfavorable solute-polymer interactions, the nanoparticle size increases with hydrophobic block length. However, for favorable solute-polymer interactions, nanoparticle size exhibits a local minimum with respect to the hydrophobic block length. Our results provide insights on ways in which experimentally controllable parameters of the Flash NanoPrecipitation process can be used to influence aggregate size and composition during self-assembly. (C) 2011 American Institute of Physics. [doi: 10.1063/1.3653379]

A microfluidic device is employed to emulate water droplet emergence from a porous electrode and slug formation in the gas flow channel of a PEM fuel cell. Liquid water emerges from a 50 mu m pore forming a droplet; the droplet grows to span the entire cross-section of a microchannel and transitions into a slug which detaches and is swept downstream. Droplet growth, slug formation, detachment, and motion are analyzed using high-speed video images and pressure-time traces. Slug volume is controlled primarily by channel geometry, interfacial forces, and gravity. As water slugs move downstream, they leave residual micro-droplets that act as nucleation sites for the next droplet-to-slug transition. Residual liquid in the form of micro-droplets results in a significant decrease in slug volume between the very first slug formed in an initially dry channel and the ultimate "steady-state" slug. A physics-based model is presented to predict slug volumes and pressure drops for slug detachment and motion. (C) 2011 Elsevier B.V. All rights reserved.

Dynamic and steady-state water flux, current density, and resistance across a Nafion 115 membrane-electrode-assembly (MEA) were measured as functions of temperature, water activity, and applied potential. After step changes in applied potential, the current, MEA resistance, and water flux evolved to new values over 3000-5000 s, indicating a slow redistribution of water in the membrane. Steady-state current density initially increased linearly with increasing potential and then saturated at higher applied potentials. Water flux increases in the direction of current flow resulting from electro-osmotic drag. The coupled transport of water and protons was modeled with an explicit accounting for electro-osmotic drag, water diffusion, and interfacial water transport resistance across the vapor/membrane interface. The model shows that water is dragged inside the membrane by the proton current, but the net water flux into and out of the membrane is controlled by interfacial water transport at the membrane/vapor interface. The coupling of electro-osmotic drag and interfacial water transport redistributes the water in the membrane. Because water entering the membrane is limited by interfacial transport, an increase in current depletes Water from the anode side of the membrane, increasing the membrane resistance there, which in turn limits the current This feedback, loop between current density and membrane resistance determines the stable steady-state operation at a fixed applied potential that results in current saturation. We show that interfacial water transport resistance substantially reduces the impact of :Osmotic drag on polymer electrolyte membrane fuel cell operation.

We propose a systematic, rigorous mathematical optimization methodology for the construction, "on demand," of network structures that are guaranteed to possess a prescribed collective property: the degree-dependent clustering. The ability to generate such realizations of networks is important not only for creating artificial networks that can perform desired functions, but also to facilitate the study of networks as part of other algorithms. This problem exhibits large combinatorial complexity and is difficult to solve with off-the-shelf commercial optimization software. To that end, we also present a customized preprocessing algorithm that allows us to judiciously fix certain problem variables and, thus, significantly reduce computational times. Results from the application of the framework to data sets resulting from simulations of an acquaintance network formation model are presented.

Nonlinear dimensionality reduction techniques can be applied to molecular simulation trajectories to systematically extract a small number of variables with which to parametrize the important dynamical motions of the system. For molecular systems exhibiting free energy barriers exceeding a few k(B)T, inadequate sampling of the barrier regions between stable or metastable basins can lead to a poor global characterization of the free energy landscape. We present an adaptation of a nonlinear dimensionality reduction technique known as the diffusion map that extends its applicability to biased umbrella sampling simulation trajectories in which restraining potentials are employed to drive the system into high free energy regions and improve sampling of phase space. We then propose a bootstrapped approach to iteratively discover good low-dimensional parametrizations by interleaving successive rounds of umbrella sampling and diffusion mapping, and we illustrate the technique through a study of alanine dipeptide in explicit solvent. (C) 2011 American Institute of Physics. [doi: 10.1063/1.3574394]

Molecular simulation is an important and ubiquitous tool in the study of microscopic phenomena in fields as diverse as materials science, protein folding and drug design. While the atomic-level resolution provides unparalleled detail, it can be non-trivial to extract the important motions underlying simulations of complex systems containing many degrees of freedom. The diffusion map is a nonlinear dimensionality reduction technique with the capacity to systematically extract the essential dynamical modes of high-dimensional simulation trajectories, furnishing a kinetically meaningful low-dimensional framework with which to develop insight and understanding of the underlying dynamics and thermodynamics. We survey the potential of this approach in the field of molecular simulation, consider its challenges, and discuss its underlying concepts and means of application. We provide examples drawn from our own work on the hydrophobic collapse mechanism of n-alkane chains, folding pathways of an antimicrobial peptide, and the dynamics of a driven interface. (C) 2011 Elsevier B.V. All rights reserved.

We experiment with the the use of diffusion maps to parameterize a low-dimensional attracting manifold of a high-dimensional ODE system in order to be able to repeatedly compute trajectories on that manifold. By integrating directly on that manifold we can avoid the stiffness of the original system.

The extracellular signal-regulated kinase (ERK) controls cellular processes by phosphorylating multiple substrates. The ERK protein can use the same domains to interact with phosphatases, which dephosphorylate and deactivate ERK, and with substrates, which connect ERK to its downstream effects. As a consequence, substrates can compete with phosphatases and control the level of ERK phosphorylation. We propose that this effect can qualitatively change the dynamics of a network that controls ERK activation. On its own, this network can be bistable, but in a larger system, where ERK accelerates the degradation of a substrate competing with a phosphatase, this network can oscillate. Previous studies proposed that oscillatory ERK signaling requires a negative feedback in which active ERK reduces the rate at which it is phosphorylated by upstream kinase. We argue that oscillations can also emerge even when this rate is constant, due to substrate-dependent control of ERK phosphorylation.

The persistently fast evolutionary behavior of certain differential systems may have intrinsically slow features. We consider systems whose solution trajectories are slowly changing distributions and assume that we do not have access to the equations of the system, only to a simulator or legacy code that performs step-by-step time integration of the system. We characterize the set of all possible instantaneous solutions by points on a low-dimensional virtual slow manifold (VSM) and show how, when there is a sufficiently large gap between the time scales of the fast and slow behaviors, we can restrict fast behavior observed numerically through simulation to the VSM and perform useful computations more efficiently there.

A model fuel cell with a single transparent straight flow channel and segmented anode was constructed to measure the direct correlation of liquid water movement with the local currents along the flow channel. Water drops emerge through the largest pores of the GDL with the size of the droplets that emerge on the surface determined by the size of the pore and its location under the gas flow channel or under the land. Gravity, surface tension, and the shearing force from the gas flow control the movement of liquid in the gas flow channel. By creating a single large diameter pore in the GDL, liquid water flow emergent from the GDL was forced to be in specific locations along the length of the channel and either under the land or under the channel. The effects of gravity were amplified when the large pore was under the channel, but diminished with the large pore under the land. Current fluctuations were minimised when the dominant water transport from the GDL pore was near the cathode outlet. The results show that it is possible to engineer the water distribution in PEM fuel cells by modifying the pore sizes in the GDL.

We consider a two-layer, one-dimensional lattice of neurons; one layer consists of excitatory thalamocortical neurons, while the other is comprised of inhibitory reticular thalamic neurons. Such networks are known to support "lurching" waves, for which propagation does not appear smooth, but rather progresses in a saltatory fashion; these waves can be characterized by different spatial widths (different numbers of neurons active at the same time). We show that these lurching waves are fixed points of appropriately defined Poincare maps, and follow these fixed points as parameters are varied. In this way, we are able to explain observed transitions in behavior, and, in particular, to show how branches with different spatial widths are linked with each other. Our computer-assisted analysis is quite general and could be applied to other spatially extended systems which exhibit this non-trivial form of wave propagation.

The study of particle coagulation and sintering processes is important in a variety of research studies ranging from cell fusion and dust motion to aerosol formation applications. These processes are traditionally simulated using either Monte-Carlo methods or integro-differential equations for particle number density functions. In this paper, we present a computational technique for cases where we believe that accurate closed evolution equations for a finite number of moments of the density function exist in principle, but are not explicitly available. The so-called equation-free computational framework is then employed to numerically obtain the solution of these unavailable closed moment equations by exploiting (through intelligent design of computational experiments) the corresponding fine-scale (here, Monte-Carlo) simulation. We illustrate the use of this method by accelerating the computation of evolving moments of uni- and bivariate particle coagulation and sintering through short simulation bursts of a constant-number Monte-Carlo scheme. (C) 2010 Elsevier Inc. All rights reserved.

The anterior region of the Drosophila embryo is patterned by the concentration gradient of the homeodomain transcription factor bicoid (Bcd). The Bcd gradient was the first identified morphogen gradient and continues to be a subject of intense research at multiple levels, from the mechanisms of RNA localization in the oocyte to the evolution of the Bcd-mediated patterning events in multiple Drosophila species. Critical assessment of the mechanisms of the Bcd gradient formation requires biophysical models of the syncytial embryo. Most of the proposed models rely on reaction diffusion equations, but their formulation and applicability at high nuclear densities is a nontrivial task. We propose a straightforward alternative in which the syncytial blastoderm is approximated by a periodic arrangement of well-mixed compartments: a single nucleus and an associated cytoplasmic region. We formulate a compartmental model, constrain its parameters by experimental data, and demonstrate that it provides an adequate description of the Bcd gradient dynamics. (C) 2010 Elsevier Inc. All rights reserved.

We focus on the "trijunction" between multiscale computations, bifurcation theory and social networks. In particular, we address how the Equation-Free approach, a recently developed computational framework, can be exploited to systematically extract coarse-grained, emergent dynamical information by bridging detailed, agent-based models of social interactions on networks, with macroscopic, systems-level, continuum numerical analysis tools. For our illustrations, we use a simple dynamic agent-based model describing the propagation of information between individuals interacting under mimesis in a social network with private and public information. We describe the rules governing the evolution of the agents' emotional state dynamics and discover, through simulation, multiple stable stationary states as a function of the network topology. Using the Equation-Free approach we track the dependence of these stationary solutions on network parameters and quantify their stability in the form of coarse-grained bifurcation diagrams.

We obtain a Fokker-Planck equation describing experimental data on the collective motion of locusts. The noise is of internal origin and due to the discrete character and finite number of constituents of the swarm. The stationary probability distribution shows a rich phenomenology including nonmonotonic behavior of several order and disorder transition indicators in noise intensity. This complex behavior arises naturally as a result of the randomness in the system. Its counterintuitive character challenges standard interpretations of noise induced transitions and calls for an extension of this theory in order to capture the behavior of certain classes of biologically motivated models. Our results suggest that the collective switches of the group's direction of motion might be due to a random ergodic effect and, as such, they are inherent to group formation.

A Lattice Monte Carlo (LMC) simulation technique has been developed to describe the synthesis of a single semiconductor nanocrystal inside the droplets of a microemulsion. The LMC simulation can track the diffusion of a precursor, its irreversible reaction with a second precursor to form nuclei, and the diffusion and coalescence of the nuclei into clusters and eventually into a single particle. In this paper we compare two scenarios for forming a single nanocrystal. The first scenario involves very rapid (spontaneous) conversion of a precursor dispersed in the droplet to nuclei that diffuse and coalesce into a single particle. The second scenario involves diffusion of a precursor to the droplet interface where an irreversible reaction with a second precursor forms nuclei that subsequently diffuse into the droplet and coalesce. Simulations were performed describing the synthesis of ZnSe nanocrystals with diameters up to 7 nm, i.e., below the quantum confinement threshold of 9 nm for this material. Comparison of the time required for the formation of a single final particle in each of the two cases reveals that for particles smaller than similar to 3.5 nm the formation times are nearly equal. For particles larger than 3.5 nm, the second process is completed faster than the first one. Analysis of intermediate cluster populations indicates that the formation of a larger "sweeper" cluster accelerates the rate of coalescence and is more effective when the nuclei are supplied gradually, as in the second process, compared to spontaneous nucleation throughout the domain. The kinetics of coalescence of an initially monodisperse population of nuclei in spherical domains of finite size were studied and generalized equations were obtained that describe the evolution of the number of different sizes as function of dimensionless time; this constitutes an extension to the classical analysis of coalescence of monodisperse aerosols in an infinite domain.

A lattice Monte Carlo model has been developed to describe the formation of a single semiconductor nanocrystal (quantum dot) inside a droplet of a microemulsion. The motivation stems from the need to understand the kinetics of quantum dot formation in microemulsion templates with minimal droplet droplet coalescence. In these systems, a fixed amount of a reactant is dissolved in each droplet, and another reactant is supplied by diffusion through the interface. Nucleation is facilitated by a spontaneous reaction between the precursors at the droplet interface, and the coalescence of nuclei and clusters ultimately leads to the formation of a single particle. The size of the final particle is controlled by the concentration of the first reactant, A hard-sphere potential is used to describe cluster cluster interactions. The overall particle formation time initially increases with final particle size, quickly passes through a maximum, and subsequently decreases due to the formation of large intermediate clusters apparently acting as effective collision partners to smaller ones. Studies of the evolution of intermediate cluster sizes provided mechanistic details of the final particle formation through cluster cluster coalescence. A generalized dimensionless equation is obtained that relates the formation time of the final particle to its size for various droplet sizes and diffusivities of the first reactant and clusters. A parametric study reveals that the final particle formation time is more sensitive to changes in the cluster cluster coalescence probability than in the probability of nucleation.

This paper presents a review of recent publications on mathematical modeling, steady-state and dynamic behavior, and control of polyelectrolyte membrane and solid oxide fuel cells. We limited the scope of this review to these two fuel cell types, which have been studied more, and have been reported to be more promising, than other fuel cell types. Zero-, one-, two-, and three-dimensional models developed to describe the behavior of the fuel cells are reviewed. Essential components of these models are highlighted. Conditions under which a fuel cell exhibits steady state multiplicity are described. Stability of the steady states is discussed. Processes that take place inside the fuel cells and contribute to the existence of multiple time-scales in the fuel cells are examined. Control configurations and strategies proposed and used for the fuel cells are reviewed, and advantages and disadvantages of each are listed. At the end, in view of the current status of the research activities, topics that require further research studies are discussed.

We propose a (time) multiscale method for the coarse-grained analysis of collective motion and decision-making in self-propelled particle models of swarms comprising a mixture of 'naive' and 'informed' individuals. The method is based on projecting the particle configuration onto a single 'meta-particle' that consists of the elongation of the flock together with the mean group velocity and position. We find that the collective states can be associated with the transient and asymptotic transport properties of the random walk followed by the meta-particle, which we assume follows a continuous time random walk (CTRW). These properties can be accurately predicted at the macroscopic level by an advection-diffusion equation with memory (ADEM) whose parameters are obtained from a mean group velocity time series obtained from a single simulation run of the individual-based model. (C) 2010 Elsevier Ltd. All rights reserved,

A simple tanks-in-series model is presented, which allows for the understanding of the basic physics behind complex spatiotemporal behaviors observed in self-humidified polymer electrolyte membrane (PEM) fuel cells. Our approach is focused on how the intrinsically nonlinear dynamics of water formation couples with water transport, leading to multistability, inhomogeneous steady state current profiles through the cell and other nonlinear phenomena. We show in particular how the operating parameters determine the location of high current spots and the subsequent propagation of current waves throughout the cell during the ignition procedure. We also reproduce and explain transient current increases seen during the extinction of the cell and the unusual aspect of the polarization curves. Implications for the efficiency of self-humidified PEM fuel cells are highlighted, and possible ways to improve their performances are discussed on these bases. (C) 2009 Elsevier Ltd. All rights reserved.

Binocular rivalry occurs when two very different images are presented to the two eyes, but a subject perceives only one image at a given time. A number of computational models for binocular rivalry have been proposed; most can be categorised as either "rate" models, containing a small number of variables, or as more biophysically-realistic "spiking neuron" models. However, a principled derivation of a reduced model from a spiking model is lacking. We present two such derivations, one heuristic and a second using recently-developed data-mining techniques to extract a small number of "macroscopic" variables from the results of a spiking neuron model simulation. We also consider bifurcations that can occur as parameters are varied, and the role of noise in such systems. Our methods are applicable to a number of other models of interest.

We employ the diffusion map approach as a nonlinear dimensionality reduction technique to extract a dynamically relevant, low-dimensional description of n-alkane chains in the ideal-gas phase and in aqueous solution. In the case of C(8) we find the dynamics to be governed by torsional motions. For C(16) and C(24) we extract three global order parameters with which we characterize the fundamental dynamics, and determine that the low free-energy pathway of globular collapse proceeds by a "kink and slide" mechanism, whereby a bend near the end of the linear chain migrates toward the middle to form a hairpin and, ultimately, a coiled helix. The low-dimensional representation is subtly perturbed in the solvated phase relative to the ideal gas, and its geometric structure is conserved between C(16) and C(24). The methodology is directly extensible to biomolecular self-assembly processes, such as protein folding.

A framework for the analysis of stochastic models of chemical systems for which the deterministic mean-field description is undergoing a saddle-node infinite period (SNIPER) bifurcation is presented. Such a bifurcation occurs, for example, in the modeling of cell-cycle regulation. It is shown that the stochastic system possesses oscillatory solutions even for parameter values for which the mean-field model does not oscillate. The dependence of the mean period of these oscillations on the parameters of the model (kinetic rate constants) and the size of the system (number of molecules present) are studied. Our approach is based on the chemical Fokker-Planck equation. To gain some insight into the advantages and disadvantages of the method, a simple one-dimensional chemical switch is first analyzed, and then the chemical SNIPER problem is studied in detail. First, results obtained by solving the Fokker-Planck equation numerically are presented. Then an asymptotic analysis of the Fokker-Planck equation is used to derive explicit formulae for the period of oscillation as a function of the rate constants and as a function of the system size.

We review and discuss the use of equation-free computation in extracting coarse-grained, nonlinear dynamics information from atomistic (lattice-gas) models of surface reactions. The approach is based on circumventing the explicit derivation of macroscopic equations for the system statistics (e.g., average coverage). Short bursts of appropriately initialized computational experimentation with the lattice-gas simulator are designed "on demand" and processed in the spirit of the coarse timestepper introduced in Theodoropoulos et al. (2000) (K. Theodoropoulos, Y.-H. Qian, I.G. Kevrekidis, Proc. Natl. Acad. Sci. USA 97 (2000) 9840). The information derived from these computational experiments, processed through traditional, continuum numerical methods is used to solve the macroscopic equations without ever deriving them in closed form. The approach is illustrated through two computational examples: the CO oxidation reaction, and the NO + CO/Pt(100) reaction. (C) 2009 Elsevier B.V. All rights reserved.

Developing effective descriptions of the microscopic dynamics of many physical phenomena can both dramatically enhance their computational exploration and lead to a more fundamental understanding of the underlying physics. Previously, an effective description of a driven interface in the presence of mobile impurities, based on an Ising variant model and a single empirical coarse variable, was partially successful [M. Haataja et al., Phys. Rev. Lett. 92, 160603 (2004)]; yet it underlined the necessity of selecting additional coarse variables in certain parameter regimes. In this paper we use a data mining approach to help identify the coarse variables required. We discuss the implementation of this diffusion map approach, the selection of a similarity measure between system snapshots required in the approach, and the correspondence between empirically selected and automatically detected coarse variables. We conclude by illustrating the use of the diffusion map variables in assisting the atomistic simulations and we discuss the translation of information between fine and coarse descriptions using lifting and restriction operators.

Nonlinear independent component analysis is combined with diffusion-map data analysis techniques to detect good observables in high-dimensional dynamic data. These detections are achieved by integrating local principal component analysis of simulation bursts by using eigenvectors of a Markov matrix describing anisotropic diffusion. The widely applicable procedure, a crucial step in model reduction approaches, is illustrated on stochastic chemical reaction network simulations.

In traditional physicochemical modeling, one derives evolution equations at the (macroscopic, coarse) scale of interest; these are used to perform a variety of tasks (simulation, bifurcation analysis, Optimization) using an arsenal of analytical and numerical techniques. For many complex systems, however, although one observes evolution at a macroscopic scale of interest, accurate models ire only given at a more detailed (fine-scale, microscopic) level of description (e.g., lattice Boltzmann, kinetic Monte Carlo, molecular dynamics). Here, we review a framework for computer-aided multiscale analysis, which enables macroscopic computational tasks (over extended spatiotemporal scales) using only appropriately initialized microscopic simulation on short time and length scales. The methodology bypasses the derivation of macroscopic evolution equations when these equations conceptual), exist but are not available in closed form-hence the term equation-free. We selectively discuss basic algorithms and underlying principles and illustrate the approach through representative applications. We also discuss potential difficulties and outline areas for future research.

We describe a reverse integration approach for the exploration of low-dimensional effective potential landscapes. Coarse reverse integration initialized on a ring of coarse states enables efficient navigation on the landscape terrain: Escape from local effective potential wells, detection of saddle points, and identification of significant transition paths between wells. We consider several distinct ring evolution modes: Backward stepping in time, solution arc length, and effective potential. The performance of these approaches is illustrated for a deterministic problem where the energy landscape is known explicitly. Reverse ring integration is then applied to noisy problems where the ring integration routine serves as an outer wrapper around a forward-in-time inner simulator. Two versions of such inner simulators are considered: A Gillespie-type stochastic simulator and a molecular dynamics simulator. In these "equation-free" computational illustrations, estimation techniques are applied to the results of short bursts of inner simulation to obtain the unavailable (in closed-form) quantities (local drift and diffusion coefficient estimates) required for reverse ring integration; this naturally leads to approximations of the effective landscape. (C) 2009 American Institute of Physics. [doi:10.1063/1.3207882]

Among the most striking aspects of the movement of many animal groups are their sudden coherent changes in direction. Recent observations of locusts and starlings have shown that this directional switching is an intrinsic property of their motion. Similar direction switches are seen in self-propelled particle and other models of group motion. Comprehending the factors that determine such switches is key to understanding the movement of these groups. Here, we adopt a coarse-grained approach to the study of directional switching in a self-propelled particle model assuming an underlying one-dimensional Fokker-Planck equation for the mean velocity of the particles. We continue with this assumption in analyzing experimental data on locusts and use a similar systematic Fokker-Planck equation coefficient estimation approach to extract the relevant information for the assumed Fokker-Planck equation underlying that experimental data. In the experiment itself the motion of groups of 5 to 100 locust nymphs was investigated in a homogeneous laboratory environment, helping us to establish the intrinsic dynamics of locust marching bands. We determine the mean time between direction switches as a function of group density for the experimental data and the self-propelled particle model. This systematic approach allows us to identify key differences between the experimental data and the model, revealing that individual locusts appear to increase the randomness of their movements in response to a loss of alignment by the group. We give a quantitative description of how locusts use noise to maintain swarm alignment. We discuss further how properties of individual animal behavior, inferred by using the Fokker-Planck equation coefficient estimation approach, can be implemented in the self-propelled particle model to replicate qualitatively the group level dynamics seen in the experimental data.

We study localized modes in uniform one-dimensional chains of tightly packed and uniaxially compressed elastic beads in the presence of one or two light-mass impurities. For chains composed of beads of the same type, the intrinsic nonlinearity, which is caused by the Hertzian interaction of the beads, appears not to support localized, breathing modes. Consequently, the inclusion of light-mass impurities is crucial for their appearance. By analyzing the problem's linear limit, we identify the system's eigenfrequencies and the linear defect modes. Using continuation techniques, we find the solutions that bifurcate from their linear counterparts and study their linear stability in detail. We observe that the nonlinearity leads to a frequency dependence in the amplitude of the oscillations, a static mutual displacement of the parts of the chain separated by a defect, and for chains with two defects that are not in contact, it induces symmetry-breaking bifurcations.

Soluble peptides, susceptible to degradation and clearance in therapeutic applications, have been formulated into protected nanoparticles for the first time through the process of kinetically controlled, block copolymer directed rapid precipitation using Flash NanoPrecipitation. Complementary Brownian dynamics simulations qualitatively model the nanoparticle formation process. The simulations corroborate the hypothesis that the size of nanoparticles decreases with increasing supersaturation. Additionally, the influence of the polymer-peptide interaction energy on the efficiency of nanoparticle protection by polymer surface coverage is elucidated in both experiments and simulations.

The persistently fast evolutionary behavior of certain differential systems may have intrinsically slow features. We consider fast stochastic systems whose solution trajectories are slowly changing distributions and assume that we do not have access to the equations of the system, only to a simulator or legacy code that performs step-wise time integration of the system. We characterize the set of all possible instantaneous solutions by points on a low-dimensional virtual slow manifold (VSM) and show how, when there is a sufficiently large gap between the time scales of the fast and slow behaviors, we can restrict the fast behavior observed numerically through simulation to the VSM and perform important computations there.

We study the effects of a signalling constraint on an individual-based model of self-organizing group formation using a coarse analysis framework. This involves using an automated data-driven technique which defines a diffusion process on the graph of a sample dataset formed from a representative stationary simulation. The eigenvectors of the graph Laplacian are used to construct 'diffusion-map' coordinates which provide a geometrically meaningful low-dimensional representation of the dataset. We show that, for the parameter regime studied, the second principal eigenvector provides a sufficient representation of the dataset and use it as a coarse observable. This allows the computation of coarse bifurcation diagrams, which are used to compare the effects of the signalling constraint on the population-level behavior of the model. (C) 2008 Elsevier Inc. All rights reserved.

The thermally induced order-to-disorder transition of a monolayer of krypton (Kr) atoms adsorbed on a graphite surface is studied based on a coarse molecular-dynamics (CMD) approach for the bracketing and location of the transition onset. A planar order parameter is identified as a coarse variable, psi, that can describe the macroscopic state of the system. Implementation of the CMD method enables the construction of the underlying effective free-energy landscapes from which the transition temperature, T(t), is predicted. The CMD prediction of T(t) is validated by comparison with predictions based on conventional molecular-dynamics (MD) techniques. The conventional MD computations include the temperature dependence of the planar order parameter, the specific heat, the Kr-Kr pair correlation function, the mean square displacement and corresponding diffusion coefficient, as well as the equilibrium probability distribution function of Kr-atom coordinates. Our findings suggest that the thermally induced order-to-disorder transition at the conditions examined in this study appears to be continuous. The CMD implementation provides substantial computational gains over conventional MD.

We present a coarse-grained, implicit solvent model for polystyrene-b-poly(ethylene oxide) in aqueous solution and Study its assembly kinetics using Brownian dynamics Simulations. The polymer is modeled as a chain of freely jointed beads interacting through effective potentials. Coarse-grained force field parameters are determined by matching experimental thermodynamic quantities including radius of gyration, second virial coefficient, aggregation number, and critical micelle concentration. We investigate the influence of cooling rate (analogous to the rate of solvent quality change in rapid precipitations), polymer concentration, and friction coefficient on the assembly kinetics and compare simulation results to flash nanoprecipitation experiments. We find that assembly kinetics show a linear scaling relation with inverse friction coefficient when the friction coefficient is larger than 1. When the cooling time is less than the characteristic micellization time, stable kinetically arrested clusters are obtained; otherwise, close-to-equilibrium micelles are formed. The characteristic micellization time is estimated to be only 3-6 ms, in contrast to 30-40 ms previously determined in experiments. We suggest that previous experiments probed the formation of micellar clusters while simulations in this work studied the kinetics of a single micelle assembled from free polymer chains.

We explore a computational approach to coarse graining the evolution of the large-scale features of a randomly forced Burgers equation in one spatial dimension. The long term evolution of the solution energy spectrum appears self-similar in time. We demonstrate coarse projective integration and coarse dynamic renormalization as tools that accelerate the extraction of macroscopic information (integration in time, self-similar shapes, nontrivial dynamic exponents) from short bursts of appropriately initialized direct simulation. These procedures solve numerically an effective evolution equation for the energy spectrum without ever deriving this equation in closed form. (c) 2008 American Institute of Physics.

The concise representation of complex high dimensional stochastic systems via a few reduced coordinates is an important problem in computational physics, chemistry, and biology. In this paper we use the first few eigenfunctions of the backward Fokker-Planck diffusion operator as a coarse-grained low dimensional representation for the long-term evolution of a stochastic system and show that they are optimal under a certain mean squared error criterion. We denote the mapping from physical space to these eigenfunctions as the diffusion map. While in high dimensional systems these eigenfunctions are difficult to compute numerically by conventional methods such as finite differences or finite elements, we describe a simple computational data-driven method to approximate them from a large set of simulated data. Our method is based on de. ning an appropriately weighted graph on the set of simulated data and computing the first few eigenvectors and eigenvalues of the corresponding random walk matrix on this graph. Thus, our algorithm incorporates the local geometry and density at each point into a global picture that merges data from different simulation runs in a natural way. Furthermore, we describe lifting and restriction operators between the diffusion map space and the original space. These operators facilitate the description of the coarse-grained dynamics, possibly in the form of a low dimensional effective free energy surface parameterized by the diffusion map reduction coordinates. They also enable a systematic exploration of such effective free energy surfaces through the design of additional "intelligently biased" computational experiments. We conclude by demonstrating our method in a few examples.

We explore the effect of spatiotemporally varying substrate temperature profiles on the dynamics and the resulting reaction rate enhancement of the catalytic oxidation of CO on Pt(110). The catalytic surface is "addressed" by a focused laser beam whose motion is computer controlled. The averaged reaction rate is observed to undergo a characteristic maximum as a function of the speed of this moving laser spot. Experiments as well as modeling are used to explore and rationalize the existence of such an optimal laser speed.