Research, Projects and Interests

Research Interests

Projects

Publications

• Modeling Focal Epileptic Activity in the Wilson-Cowan model with Depolarization Block. H.G.E. Meijer, T.L. Eissa, B. Kiewiet, J.F. Neuman, C.A. Schevon, R.G. Emerson, R.R. Goodman, G.M. McKhann Jr, C.J. Marcuccilli, A.K. Tryba, J.D. Cowan, S.A. van Gils and W. van Drongelen, Journal of Mathematical Neuroscience, Accepted, 2015 Abstract
  Description:

  Measurements of neuronal signals during human seizure activity and evoked epileptic activity in experimental models suggest that, in these pathological states, the individual nerve cells experience an activity driven depolarization block, i.e. they saturate. We examined the effect of such a saturation in the Wilson-Cowan formalism by adapting the nonlinear activation function; we substituted the commonly applied sigmoid for a Gaussian function. We discuss experimental recordings during a seizure that support this substitution. Next we perform a bifurcation analysis on the Wilson-Cowan model with a Gaussian activation function. The main effect is an additional stable equilibrium with high excitatory and low inhibitory activity. Analysis of coupled local networks then shows that such high activity can stay localized or spread. Specifically, in a spatial continuum we show a wavefront with inhibition leading followed by excitatory activity. We relate our model simulations to observations of spreading activity during seizures.

  DOI
• Tracking pattern evolution through extended center manifold reduction and singular perturbations. L. Sewalt, A. Doelman, H.G.E. Meijer, V. Rottschäfer, A. Zagaris, Physica D 2015, Vol. 298-299, pp. 48-67 Abstract
  Description:

  In this paper we develop an extended center manifold reduction method: a methodology to analyze the formation and bifurcations of small-amplitude patterns in certain classes of multi-component, singularly perturbed systems of partial differential equations. We specifically consider systems with a spatially homogeneous state whose stability spectrum partitions into eigenvalue groups with distinct asymptotic properties. One group of successive eigenvalues in the bifurcating group are widely interspaced, while the eigenvalues in the other are stable and cluster asymptotically close to the origin along the stable semi-axis. The classical center manifold reduction provides a rigorous framework to analyze destabilizations of the trivial state, as long as there is a spectral gap of sufficient width. When the bifurcating eigenvalue becomes commensurate to the stable eigenvalues clustering close to the origin, the center manifold reduction breaks down. Moreover, it cannot capture subsequent bifurcations of the bifurcating pattern. Through our methodology, we formally derive expressions for low-dimensional manifolds exponentially attracting the full flow for parameter combinations that go beyond those allowed for the (classical) center manifold reduction, i.e. to cases in which the spectral gap condition no longer can be satisfied. Our method also includes an explicit description of the flow on these manifolds and thus provides an analytical tool to study subsequent bifurcations. Our analysis centers around primary bifurcations of transcritical type -- that can be either of co-dimension 1 or 2 -- in two- and three-component PDE systems. We employ our method to study bifurcation scenarios of small-amplitude patterns and the possible appearance of low-dimensional spatio-temporal chaos. We also exemplify our analysis by a number of characteristic reaction-diffusion systems with disparate diffusivities.

  DOI
• Exploiting pallidal plasticity for stimulation in Parkinson's disease, M.A.J. Lourens, B. Schwab, J.A. Nirody, H.G.E. Meijer and S.A. van Gils, Journal of Neural Engineering, 2015, 12: 026005 Abstract
  Description:

  Objective. Continuous application of high-frequency deep brain stimulation (DBS) often effectively reduces motor symptoms of Parkinson's disease patients. While there is a growing need for more effective and less traumatic stimulation, the exact mechanism of DBS is still unknown. Here, we present a methodology to exploit the plasticity of GABAergic synapses inside the external globus pallidus (GPe) for the optimization of DBS. Approach. Assuming the existence of spike-timing-dependent plasticity (STDP) at GABAergic GPe-GPe synapses, we simulate neural activity in a network model of the subthalamic nucleus and GPe. In particular, we test different DBS protocols in our model and quantify their influence on neural synchrony. Main results. In an exemplary set of biologically plausible model parameters, we show that STDP in the GPe has a direct influence on neural activity and especially the stability of firing patterns. STDP stabilizes both uncorrelated firing in the healthy state and correlated firing in the parkinsonian state. Alternative stimulation protocols such as coordinated reset stimulation can clearly profit from the stabilizing effect of STDP. These results are widely independent of the STDP learning rule. Significance. Once the model settings, e.g., connection architectures, have been described experimentally, our model can be adjusted and directly applied in the development of novel stimulation protocols. More efficient stimulation leads to both minimization of side effects and savings in battery power.

  DOI
• Bifurcation analysis of a model for atherosclerotic plaque evolution. M.A.K. Bulelzai, J.L.A. Dubbeldam, H.G.E. Meijer, Physica D 2014, Vol. 278-9, pp. 31-43 Abstract
  Description:

  We analyze two ordinary differential equation (ODE) models for atherosclerosis. The ODE models describe long time evolution of plaques in arteries. We show how the dynamics of the first atherosclerosis model (model A) can be understood using codimension-two bifurcation analysis. The Low-Density Lipoprotein (LDL) intake parameter (d) is the first control parameter and the second control parameter is either taken to be the conversion rate of macrophages (b) or the wall shear stress (σ). Our analysis reveals that in both cases a Bogdanov-Takens (BT) point acts as an organizing center. The bifurcation diagrams are calculated partly analytically and to a large extent numerically using AUTO07 and MATCONT. The bifurcation curves show that the concentration of LDL in the plaque as well as the monocyte and the macrophage concentration exhibit oscillations for a certain range of values of the control parameters. Moreover, we find that there are threshold values for both the cholesterol intake rate dcritdcrit and the conversion rate of the macrophages bcritbcrit, which depend on the values of other parameters, above which the plaque volume increases with time. It is found that larger conversion rates of macrophages lower the threshold value of cholesterol intake and vice versa. We further argue that the dynamics for model A can still be discerned in the second model (model B) in which the slow evolution of the radius of the artery is coupled self-consistently to changes in the plaque volume. The very slow evolution of the radius of the artery compared to the other processes makes it possible to use a slow manifold approximation to study the dynamics in this case. We find that in this case the model predicts that the concentrations of the plaque constituents may go through a period of oscillations before the radius of the artery will start to decrease. These oscillations hence act as a precursor for the reduction of the artery radius by plaque growth.

  DOI
• Improved homoclinic predictor for Bogdanov-Takens bifurcation. Yu.A. Kuznetsov, H.G.E. Meijer, B. Al Hdaibat and W. Govaerts, Internat. J. Bifur. Chaos Appl. Sci. Engrg., 2014, Vol. 24, No. 4: 1450057 Abstract
  Description:

  An improved homoclinic predictor at a generic codim 2 Bogdanov-Takens (BT) bifucation is derived. We use the classical "blow-up" technique to reduce the canonical smooth normal form near a generic BT bifurcation to a perturbed Hamiltonian system. With a simple perturbation method, we not only rederive the known approximation for the bifurcation parameter values but also obtain an explicit first-order correction of the unperturbed homoclinic orbit. To obtain the normal form on the center manifold, we apply the standard parameter-dependent center manifold reduction combined with the normalization, that is based on the Fredholm solvability of the homological equation. However, by systematically solving all linear systems appearing from the homological equation, we remove an ambiguity in the parameter transformation existing in the literature. The actual implementation of the improved predictor in MatCont and numerical examples illustrating its effciency are discussed.

  DOI
• Travelling waves in models of neural tissue: from localised structures to periodic waves, H.G.E. Meijer and S. Coombes, EPJ Nonlinear Biomedical Physics 2014, 2:3 Abstract
  Description:

  We consider travelling waves (fronts, pulses and periodics) in spatially extended one dimensional neural field models. We demonstrate for an excitatory field with linear adaptation that, in addition to an expected stable pulse solution, a stable anti-pulse can exist. Varying the adaptation strength we unravel the organizing centers of the bifurcation diagram for fronts and pulses, with a mixture of exact analysis for a Heaviside firing rate function and novel numerical schemes otherwise. These schemes, for non-local models with space-dependent delays, further allow for the construction and continuation of periodic waves. We use them to construct the dispersion curve -- wave speed as a function of period -- and find that they can be oscillatory and multi-valued, suggesting bistability of periodic waves. A kinematic theory predicts the onset of wave instabilities at stationary points in the dispersion curve, leading to period doubling behaviour, and is confirmed with direct numerical simulations. We end with a discussion of how the construction of dispersion curves may allow a useful classification scheme of neural field models for epileptic waves.

  DOI (open access) Poster version (CNS, 2013 Paris) see DOI . Follow the link for documented simulation and continuation Code .
• Travelling waves in a neural field model with refractoriness, H.G.E. Meijer and S. Coombes, Journal of Mathematical Biology 2014, Vol. 68, No. 5, pp. 1249-1268 Abstract
  Description:

  At one level of abstraction neural tissue can be regarded as a medium for turning local synaptic activity into output signals that propagate over large distances via axons to generate further synaptic activity that can cause reverberant activity in networks that possess a mixture of excitatory and inhibitory connections. This output is often taken to be a firing rate, and the mathematical form for the evolution equation of activity depends upon a spatial convolution of this rate with a fixed anatomical connectivity pattern. Such formulations often neglect the metabolic processes that would ultimately limit synaptic activity. Here we reinstate such a process, in the spirit of an original prescription by Wilson and Cowan [1972, ``Excitatory and inhibitory interactions in localised populations of model neurons'', Biophysical Journal 12, 1-24.], using a term that multiplies the usual spatial convolution with a moving time average of local activity over some refractory time-scale. This modulation can substantially affect network behaviour, and in particular give rise to periodic travelling waves in a purely excitatory network (with exponentially decaying anatomical connectivity), which in the absence of refractoriness would only support travelling fronts. We construct these solutions numerically as stationary periodic solutions in a co-moving frame (of both an equivalent delay differential model as well as the original delay integro-differential model). Continuation methods are used to obtain the dispersion curve for periodic travelling waves (speed as a function of period), and found to be reminiscent of those for spatially extended models of excitable tissue. A kinematic analysis (based on the dispersion curve) predicts the onset of wave instabilities, which are confirmed numerically.

  DOI (open access)
• Tracking of nociceptive thresholds using adaptive psychophysical methods, R.J. Doll, J.R. Buitenweg, H.G.E. Meijer and P. Veltink, Behavior Research Methods 2014, Vol. 46, No. 1, pp. 55-66 Abstract
  Description:
  Tracking psychophysical thresholds over time is beneficial for the investigation of dynamic behavior of underlying mechanisms. For example, noxious events can activate endogenous analgesic mechanisms reflected in an increased nociceptive threshold. Psychophysical techniques for tracking these thresholds should include efficient methods for stimulus selection and threshold estimation. This study compares the performance of different adaptive stimulus selection procedures and estimation methods in simulation and human psychophysical experiments.
  Monte Carlo simulations were performed to compare the bias and precision of threshold estimates between three stimulus selection procedures (simple staircase, random staircase and a minimum entropy procedure) and two estimation methods (logistic regression and Bayesian estimation). Logistic regression and Bayesian estimations resulted in similar precision only when prior probability distributions were appropriately chosen. The minimum entropy (with an appropriate prior distribution) and simple staircase procedures showed the highest precision but the random staircase procedure was less sensitive to different procedure specific settings than the other procedures.
  The simple staircase and random staircase stimulus selection procedures were compared in a human subject study (n=30), using logistic regression for threshold estimation. Electrocutaneous stimulation was used to track the nociceptive perception threshold before, during and after a cold pressor task as a conditioning stimulus. The two stimulus selection procedures were able to detect habituation as well as changes induced by the conditioning stimulus, but the random staircase procedure achievedthe highest precision.
  Based on these results, we recommend the random staircase procedure in combination with logistic regression for nonstationary threshold tracking experiments.
  DOI
• Analysis of bifurcations of limit cycles with Lyapunov exponents and numerical normal forms, V. De Witte, W. Govaerts, Yu.A. Kuznetsov, H.G.E. Meijer, Physica D 2014, Vol. 269, pp. 126-141 Abstract
  Description:

  Explicit computational formulas for coefficients of the periodic normal forms of the three most complex codim 2 bifurcations of limit cycles with dimension of the center manifold equal to 4 or to 5 in generic autonomous ODEs are derived. The resulting formulas are independent of the dimension of the phase space and involve solutions of certain boundary-value problems as well as multilinear functions from the Taylor expansion of the ODE right-hand side near the cycle. The formulas allow one to distinguish between the complicated bifurcation scenarios which can happen near these codim 2 bifurcations of limit cycles, where 3-tori and 4-tori can be present. We apply our techniques to the study of a known laser model, a novel model from population biology, and one for mechanical vibrations; these models exhibit Limit Point--Neimark-Sacker, Period-Doubling--Neimark-Sacker and double Neimark-Sacker bifurcations. Lyapunov exponents are computed to numerically confirm the results of the normal form analysis, in particular with respect to the existence of stable invariant tori of various dimensions and chaos.

  DOI
• Epileptic spikes are phase-locked to rhythmic delta activity in non-convulsive status epilepticus. R. Hindriks, H.G.E. Meijer, S.A. van Gils and M.J.A.M van Putten, Frontiers in Systems Neuroscience 2013, Vol. 7: 111 Abstract
  Description:

  The EEG of patients in non-convulsive status epilepticus (NCSE) often displays delta oscillations or generalized spike-wave discharges. In some patients, these delta oscillations coexist with intermittent epileptic spikes. In this study we verify the prediction of a computational model of the thalamo-cortical system that these spikes are phase-locked to the delta oscillations. We subsequently describe the physiological mechanism underlying this observation as suggested by the model. It is suggested that the spikes reflect inhibitory stochastic fluctuations in the input to thalamo-cortical relay neurons and phase-locking is a consequence of differential excitability of relay neurons over the delta cycle. Further analysis shows that the observed phase-locking can be regarded as a stochastic precursor of generalized spike-wave discharges. This study thus provides an explanation of intermittent spikes during delta oscillations in NCSE and might be generalized to other encephathologies in which delta activity can be observed.

  DOI
• Functional neuronal activity and connectivity within the subthalamic nucleus in Parkinson's disease, M.A.J. Lourens, H.G.E. Meijer, M.F. Contarino, P. van den Munckhof, P.R. Schuurman, S.A. van Gils and L.J. Bour, Clinical Neurophysiology 2013, Vol. 124, No. 5, pp 967-981 Abstract
  Description:

  Objective: Characterization of the functional neuronal activity and connectivity within the subthalamic nucleus (STN) in patients with Parkinson's disease (PD).
  Methods: Single units were extracted from micro-electrode recording (MER) of 18 PD patients who underwent STN deep brain stimulation (DBS) surgery. The firing rate and pattern of simultaneously recorded spike trains and their coherence were analyzed. To provide a precise functional assignment of position to the observed activities, for each patient we mapped its classified multichannel STN MERs to a generic atlas representation with a sensorimotor part and a remaining part.
  Results: Within the sensorimotor part we found significantly higher mean firing rate (P < 0.05) and significantly more burst-like activity (P < 0.05) than within the remaining part. The proportion of significant coherence in the beta band (13-30 Hz) is significantly higher in the sensorimotor part of the STN than elsewhere (P=0.015).
  Conclusions: The STN sensorimotor part distinguishes itself from the remaining part with respect to beta coherence, firing rate and burst-like activity and postoperatively was found as the preferred target area. Significance: Our firing behavior analysis may help to discriminate the STN sensorimotor part for the placement of the DBS electrode.

  DOI
• Analysis of Stability and Bifurcations of Fixed Points and Periodic Solutions of a Lumped Model of Neocortex with Two Delays. S. Visser S, H.G.E. Meijer, M.J.A.M. van Putten and S.A. van Gils, The Journal of Mathematical Neuroscience 2012, 2:8 Abstract
  Description:

  A lumped model of neural activity in neocortex is studied to identify regions of multi-stability of both steady states and periodic solutions. Presence of both steady states and periodic solutions is considered to correspond with epileptogenesis. The model, which consists of two delay differential equations with two fixed time lags, is mainly studied for its dependency on varying connection strength between populations. Equilibria are identified, and using linear stability analysis, all transitions are determined under which both trivial and non-trivial fixed points lose stability. Periodic solutions arising at some of these bifurcations are numerically studied with a two-parameter bifurcation analysis.

  DOI
• The Hopf-van der Pol System: Failure of a Homotopy Method, H.G.E. Meijer and T. Kalmar-Nagy, Differential Equations and Dynamical Systems: Vol. 20, No. 3, pp 323-328, 2012 Abstract
  Description:

  The purpose of this article to provide an explicit example where continuation based on the homotopy method fails. The example is a one-parameter homotopy for periodic orbits between two well-known nonlinear systems, the normal form of the Hopf bifurcation and the van der Pol system. Our analysis shows that various types of obstructions can make approximation over the whole range of the homotopy parameter impossible. The Hopf-van der Pol system demonstrates that homotopy methods may fail even for seemingly innocent systems.

  DOI
• From parkinsonian thalamic activity to restoring thalamic relay using Deep Brain Stimulation: new insights from computational modeling, H.G.E. Meijer, M. Krupa, H. Cagnan, M.A.J. Lourens, T. Heida, H.C.F. Martens, L.J. Bour and S.A. van Gils, Journal of Neural Engineering, 2011, 8: 066005 Abstract
  Description:

  We present a computational model of a thalamocortical relay neuron for exploring basal ganglia thalamocortical loop behavior in relation to Parkinson's disease and deep brain stimulation (DBS). Previous microelectrode, single-unit recording studies demonstrated that oscillatory interaction within and between basal ganglia nuclei is very often accompanied by synchronization at parkinsonian rest tremor frequencies (3-10 Hz). These oscillations have a profound influence on thalamic projections and impair the thalamic relaying of cortical input by generating rebound action potentials. Our model describes convergent inhibitory input received from basal ganglia by the thalamocortical cells based on characteristics of normal activity, and/or low-frequency oscillations (activity associated with Parkinson's disease). In addition to simulated input, we also used microelectrode recordings as inputs for the model. In the resting state, and without additional sensorimotor input, pathological rebound activity is generated for even mild parkinsonian input. We have found a specific stimulation window of amplitudes and frequencies for periodic input, which corresponds to highfrequency DBS, and which also suppresses rebound activity for mild and even more prominent parkinsonian input. When low-frequency pathological rebound activity disables the thalamocortical cell's ability to relay excitatory cortical input, a stimulation signal with parameter settings corresponding to our stimulation window can restore the thalamocortical cell's relay functionality.

  DOI
• The pedunculopontine nucleus as an additional target for deep brain stimulation. M.A.J. Lourens, H.G.E. Meijer, T. Heida, E. Marani and S.A. van Gils, Neural networks, 2011, Vol. 24, No. 6, pp 617-630 Abstract
  Description:

  The pedunculopontine nucleus has been suggested as a target for DBS. In this paper we propose a single compartment computational model for a PPN Type I cell and compare its dynamic behavior with experimental data. The model shows bursts after a period of hyperpolarization and spontaneous firing at 8 Hz. Bifurcation analysis of the single PPN cell shows bistability of fast and slow spiking solutions for a range of applied currents. A network model for STN, GPe and GPi produces basal ganglia output that is used as input for the PPN cell. The conductances for projections from the STN and the GPi to the PPN are determined from experimental data. The resulting behavior of the PPN cell is studied under normal and Parkinsonian conditions of the basal ganglia network. The effect of high frequency stimulation of the STN is considered as well as the effect of combined high frequency stimulation of the STN and the PPN at various frequencies. The relay properties of the PPN cell demonstrate that the combined high frequency stimulation of STN and low frequency (10 Hz, 25 Hz, 40 Hz) stimulation of PPN hardly improves the effect of exclusive STN stimulation. Moreover, PPN-DBS at low stimulation amplitude has a better effect than at higher stimulation amplitude. The effect of PPN output on the basal ganglia is investigated, in particular the effect of STN-DBS and/or PPN-DBS on the pathological firing pattern of STN and GPe cells. PPN-DBS eliminates the pathological firing pattern of STN and GPe cells, whereas STN-DBS and combined STN-DBS and PPN-DBS eliminate the pathological firing pattern only from STN cells.

  DOI
• Comparing Epileptiform Behavior of Mesoscale Detailed Models and Population Models of Neocortex. S. Visser, H.G.E. Meijer, H.C. Lee, W. van Drongelen, M.J.A.M. van Putten and S.A. van Gils, Journal of Clinical Neurophysiology, 2010, Vol. 27, No. 6, pp 471-478 Abstract
  Description:

  Two models of the neocortex are developed to study normal and pathologic neuronal activity. One model contains a detailed description of a neocortical microcolumn represented by 656 neurons, including superficial and deep pyramidal cells, four types of inhibitory neurons, and realistic synaptic contacts. Simulations show that neurons of a given type exhibit similar, synchronized behavior in this detailed model. This observation is captured by a population model that describes the activity of large neuronal populations with two differential equations with two delays. Both models appear to have similar sensitivity to variations of total network excitation. Analysis of the population model reveals the presence of multistability, which was also observed in various simulations of the detailed model.

  DOI
• Frequency-selectivity of a thalamocortical relay neuron during Parkinson's disease and deep brain stimulation: a computational study. H. Cagnan, H.G.E. Meijer, S.A. van Gils, M. Krupa, T. Heida, M. Rudolph, W.J. Wadman, and H.C.F. Martens, European Journal of Neuroscience, 2009, Vol. 30, No. 7, pp 1306-1317 Abstract
  Description:

  In this computational study, we investigated (i) the functional importance of correlated basal ganglia (BG) activity associated with Parkinson's disease (PD) motor symptoms by analysing the effects of globus pallidus internum (GPi) bursting frequency and synchrony on a thalamocortical (TC) relay neuron, which received GABAergic projections from this nucleus; (ii) the effects of subthalamic nucleus (STN) deep brain stimulation (DBS) on the response of the TC relay neuron to synchronized GPi oscillations; and (iii) the functional basis of the inverse relationship that has been reported between DBS frequency and stimulus amplitude, required to alleviate PD motor symptoms [A. L. Benabid et al. (1991)Lancet, 337, 403-406]. The TC relay neuron selectively responded to and relayed synchronized GPi inputs bursting at a frequency located in the range 2-25 Hz. Input selectivity of the TC relay neuron is dictated by low-threshold calcium current dynamics and passive membrane properties of the neuron. STN-DBS prevented the TC relay neuron from relaying synchronized GPi oscillations to cortex. Our model indicates that DBS alters BG output and input selectivity of the TC relay neuron, providing an explanation for the clinically observed inverse relationship between DBS frequency and stimulus amplitude.

  DOI
• Emergence and bifurcations of slow manifolds in nonlinear wave equations. T. Bakri, H.G.E. Meijer and F. Verhulst, Journal of Nonlinear Science, 2009, Vol. 19, No. 5, pp 571-596 Abstract
  Description:

  Persistence and bifurcations of Lyapunov manifolds can be studied by a combination of averaging-normalization and numerical bifurcation methods. This can be extended to infinite-dimensional cases when using suitable averaging theorems. The theory is applied to the case of a parametrically excited wave equation. We find fast dynamics in a finite, resonant part of the spectrum and slow dynamics elsewhere. The resonant part corresponds with an almost-invariant manifold and displays bifurcations into a wide variety of phenomena among which are 2- and 3-tori.

  DOI (open access)
• Numerical Bifurcation Analysis. H.G.E. Meijer, F. Dercole and B. Oldeman, in Meyers, Robert (Ed.) Encyclopedia of Complexity and Systems Science, Part 14, Pages 6329-6352, Springer New York, 2009. Abstract
  Description:

  The theory of dynamical systems studies the behavior of solutions of systems, like nonlinear ordinary differential equations (ODEs), depending upon parameters. Using qualitative methods of bifurcation theory, the behavior of the system is characterized for various parameter combinations. In particular, the catalog of system behaviors showing qualitative differences can be identified, together with the regions in parameter space where the different behaviors occur. Bifurcations delimit such regions. Symbolic and analytical approaches are in general infeasible, but numerical bifurcation analysis is a powerful tool that aids in the understanding of a nonlinear system. When computing power became widely available, algorithms for this type of analysis matured and the first codes were developed. With the development of suitable algorithms, the advancement in the qualitative theory has found its way into several software projects evolving over time. The availability of software packages allows scientists to study and adjust their models and to draw conclusions about their dynamics.

  DOI
• Numerical continuation of connecting orbits of maps. W. Govaerts, R. Khoshsiar Ghaziani, Yu.A. Kuznetsov, and H.G.E. Meijer, Journal of Difference Equations and Applications, 2009, Vol. 15, No. 8, pp 849-875 Abstract
  Description:

  We present new or improved methods to continue heteroclinic and homoclinic orbits to fixed points in iterated maps and to compute their fold bifurcation curves, corresponding to the tangency of the invariant manifolds. The proposed methods are applicable to general n-dimensional maps and are implemented in matlab. They are based on the continuation of invariant subspaces (CIS) algorithm, which is presented in a novel way. The systems of defining equations include the Riccati equations appearing in CIS for bases of the generalized stable and unstable eigenspaces. We use the bordering techniques to continue the folds, and provide full algorithmic details on how to treat the Jacobian matrix of the defining system as a sparse matrix in matlab. For a special -- but important in applications -- case n=2 we describe the first matlab implementation of known algorithms to grow one-dimensional stable and unstable manifolds of the fixed points of noninvertible maps. The methods are applied to study heteroclinic and homoclinic connections in the generalized Hénon map.

  DOI
• Switching to nonhyperbolic cycles from codim 2 bifurcations of equilibria in ODEs. Yu.A. Kuznetsov, H.G.E. Meijer, W. Govaerts and B. Sautois, Physica D 2008, Vol. 237, pp 3061-3068. Abstract
  Description:

  The paper provides full algorithmic details on switching to the continuation of all possible codim 1 cycle bifurcations from generic codim 2 equilibrium bifurcation points in n-dimensional ODEs. We discuss the implementation and the performance of the algorithm in several examples, including an extended Lorenz-84 model and a laser system.

  DOI
• New features of the software MatCont for bifurcation analysis of dynamical systems. A. Dhooge, and W. Govaerts, Yu.A. Kuznetsov, H.G.E. Meijer and B. Sautois, Mathematical and Computer Modelling of Dynamical Systems 2008, Vol. 14, No. 2, pp 147-175 Abstract
  Description:

  Bifurcation software is an essential tool in the study of dynamical systems. From the beginning (the first packages were written in the 1970's) it was also used in the modelling process, in particular to determine the values of critical parameters. More recently, it is used in a systematic way in the design of dynamical models and to determine which parameters are relevant. MatCont and Cl_MatCont are freely available matlab numerical continuation packages for the interactive study of dynamical systems and bifurcations. MatCont is the GUI-version, Cl_MatCont is the command-line version. The work started in 2000 and the first publications appeared in 2003. Since that time many new functionalities were added. Some of these are fairly simple but were never before implemented in continuation codes, e.g. Poincaré maps. Others were only available as toolboxes that can be used by experts, e.g. continuation of homoclinic orbits. Several others were never implemented at all, such as periodic normal forms for codimension 1 bifurcations of limit cycles, normal forms for codimension 2 bifurcations of equilibria, detection of codimension 2 bifurcations of limit cycles, automatic computation of phase response curves and their derivatives, continuation of branch points of equilibria and limit cycles. New numerical algorithms for these computations have been published or will appear elsewhere; here we restrict to their software implementation. We further discuss software issues that are in practice important for many users, e.g. how to define a new system starting from an existing one, how to import and export data, system descriptions, and computed results.

  DOI
• Numerical Methods for Two-Parameter Local Bifurcation Analysis of Maps. W. Govaerts, R. Khoshsiar Ghaziani, Yu.A. Kuznetsov and H.G.E. Meijer, SIAM Journal on Scientific Computing 2007, Vol. 29, No. 6, pp 2644-2667 Abstract
  Description:

  We discuss new and improved algorithms for the bifurcation analysis of fixed points and periodic orbits (cycles) of maps and their implementation in matcont, a MATLAB toolbox for continuation and bifurcation analysis of dynamical systems. This includes the numerical continuation of fixed points of iterates of the map with one control parameter, detecting and locating their bifurcation points (i.e., limit point, period-doubling, and Neimark-Sacker) and their continuation in two control parameters, as well as detection and location of all codimension 2 bifurcation points on the corresponding curves. For all bifurcations of codim 1 and 2, the critical normal form coefficients are computed, both numerically with finite directional differences and using symbolic derivatives of the original map. Using a parameter-dependent center manifold reduction, explicit asymptotics are derived for bifurcation curves of double and quadruple period cycles rooted at codim 2 points of cycles with arbitrary period. These asymptotics are implemented into the software and allow one to switch at codim 2 points to the continuation of the double and quadruple period bifurcations. We provide two examples illustrating the developed techniques: a generalized Henon map and a juvenile/adult competition model from mathematical biology.

  DOI
• Remarks on interacting Neimark-Sacker bifurcations. Yu.A. Kuznetsov and H.G.E. Meijer, Journal of Difference Equations and Applications 2006, Vol. 12, No. 10, pp 1009-1035, 2006 Abstract
  Description:

  We study codimension-2 bifurcations of fixed points of dissipative diffeomorphisms with a pair of complex eigenvalues together with either an eigenvalue $-1$ or another such a pair. In the previous studies only cubic normal forms were considered. However, in some cases the unfolding requires higher order terms and these are investigated here. We (re)derive the normal forms and reduce them to a single amplitude map. This map is similar to the amplitude system for the double-Hopf bifurcation for vector fields. We show how the critical normal form coefficients determine the general bifurcation picture for this amplitude map. Representative nonsymmetric perturbations of the normal forms are studied numerically. Our case studies show a detailed picture near various bifurcation curves, which is richer than known theoretical predictions. For arbitrary maps with these bifurcations we give explicit formulas for critical normal form coefficients on center manifolds. We apply them to an example from robotics where we are able to demonstrate the existence of a bubble-structure, which was only observed in perturbed normal forms before.

  DOI
• Generalized Hénon map and bifurcations of homoclinic tangencies. V.S. Gonchenko, Yu.A. Kuznetsov and H.G.E. Meijer, SIAM Journal on Applied Dynamical Systems 2005, Vol. 4, No. 2, pp 407-436 Abstract
  Description:

  The Generalized Hénon Map(GHM) appears in codim 2 bifurcation of homoclinic orbits for maps, for example a neutral saddle with a homoclinic tangency (NHT). The standard Hénon map with the standard parameters itself is interesting for the strange attractor, the homoclinic orbit and the resulting chaos. However, from the bifurcation point of view it turns out that many appearing bifurcations of the standard Hénon map are degenerate, for example, the strong resonances. This is not the case for GHM. One interesting feature is that the GHM has the NHT-bifurcation, so we show that the existence of infinitely many of such bifurcations.

  DOI
• Numerical normal forms for codim 2 bifurcations of fixed points with at most two critical eigenvalues. Yu.A. Kuznetsov and H.G.E Meijer, SIAM Journal on Scientific Computing 2005, Vol. 26, No. 6, pp 1932-1954 Abstract
  Description:

  We derive the formulas for the center manifold reduction and the critical normal form coefficients for the codimension two bifurcations for diffeomorphisms of the plane. Some formulas involve third or fifth order derivatives. This poses no problems for maps, but for limit cycles of ODE's it does. A RungeKuttaFehlberg-78 integrator adapted with automatic differentiation(ADOL-C 1.8) is used to obtain the required data.

  DOI
• The Fold-Flip bifurcation. Yu.A. Kuznetsov, H.G.E Meijer and L. van Veen, Internat. J. Bifur. Chaos Appl. Sci. Engrg. 2004, Vol. 14, No. 7, pp 2253-2282 Abstract
  Description:

  The fold-flip bifurcation occurs if a map has a fixed point with multipliers +1 and -1 simultaneously. In this paper the normal form of this singularity is calculated explicitly. Both local and global bifurcations of the unfolding are analysed by exploring a close relationship between the derived normal form and the truncated amplitude system for the fold-Hopf bifurcation of ODEs. Two examples are presented, the generalized H\'enon map and an extension of the Lorenz-84 model. In the latter example the first, second and third order derivatives of the Poincar\'e map are computed to find the normal form coefficients.

  DOI

Preprints/Conference Proceedings/Posters

• Computational modeling of nociceptive detection thresholds. H. Yang, H.G.E. Meijer, R.J. Doll, J.R. Buitenweg and S.A. van Gils, submitted, 2015
• A study of resonance tongues near a Chenciner bifurcation using MatcontM. W. Govaerts, Y.A. Kuznetsov, H.G.E. Meijer and N. Neirynck, ENOC 2011, Rome DOI
• On phase-locking of oscillators with delay coupling. H.G.E. Meijer and T. Kalmar-Nagy, ENOC 2011, Rome DOI
• A dynamic neural model of localization of brief successive stimuli in saltation. Meijer H.G.E., Trojan J., Kleinböhl D., Hölzl R. and Buitenweg J.R., Poster presentation at CNS 2009, Berlin DOI
• A transcritical-flip bifurcation in a model for a robot-arm. H.G.E. Meijer, SPT 2004: Symmetry and perturbation theory (Cala Gonone), pp. 209-215, World Sci. Publishing, River Edge, NJ, 2004 Abstract
  Description:

  We study a model for a robot-arm due to Steindl, Troger and Lindtner near a bifurcation of a fixed point which has multipliers +1 and -1 simultaneously. There is a degeneracy in the parameter dependent part, so that the multiplier +1 is associated to a transcritical bifurcation instead of an ordinary fold. The critical center-manifold reduction is given and the nonstandard unfolding as well.

  DOI

Ph.D Thesis

Master Thesis

• Mathematics:
  Title: On a codimension two bifurcation of fold-flip type
  Supervisor: Yu. A. Kuznetsov
  Info: See publications above
• Physics:
  Title: On the diffraction of high-amplitude picosecond acoustic pulses
  Supervisors:prof.dr. H. van Beijeren, prof. dr. J.I. Dijkhuis and dr. O.L. Muskens.
  Info: See the Soliton Project.