Generative Learning of Separatrices

# Generative Learning of Separatrices
Source: [https://arxiv.org/html/2608.14743](https://arxiv.org/html/2608.14743)
[Ellis R\. Crabtree](https://orcid.org/0000-0000-0000-0000)Affiliation:Dept\. of Predictive AnalyticsAffiliation:Johns Hopkins All Children’s HospitalAffiliation:St\. Petersburg, FL 33701Email:[ellis@jhmi\.edu](mailto:)[Dimitris G\. Giovanis](https://orcid.org/0000-0000-0000-0000)Affiliation:Dept\. of Civil and Systems EngineeringAffiliation:Johns Hopkins UniversityAffiliation:Baltimore, MD 21218Email:[dgiovan1@jhu\.edu](mailto:)[Anastasia Georgiou](https://orcid.org/0000-0000-0000-0000)Affiliation:Dept\. of Chemical and Biomolecular EngineeringAffiliation:Johns Hopkins UniversityAffiliation:Baltimore, MD 21218Email:[ageorgi3@jhu\.edu](mailto:)[George Datseris](https://orcid.org/0000-0002-6427-2385)Affiliation:Dept\. of Mathematics and StatisticsAffiliation:University of ExeterAffiliation:Exeter, EX4 4QJEmail:[g\.datseris@exeter\.ac\.uk](mailto:)[Ioannis G\. Kevrekidis](https://orcid.org/0000-0000-0000-0000)Affiliation:Dept\. of Chemical and Biomolecular EngineeringAffiliation:Johns Hopkins UniversityAffiliation:Baltimore, MD 21218Email:[yannisk@jhu\.edu](mailto:)

###### Abstract

The identification and reconstruction of the boundaries separating basins of attraction in multistable, multidimensional dynamical systems presents a fundamental challenge in computational dynamics\. These \(often geometrically complex\) structures govern transition pathways and other important large timescale behavior, yet they remain typically under\-sampled since their neighborhood does not get routinely visited during direct simulations\. Traditional computational approaches targeted for their approximate construction \(including manifold continuation, edge tracking, and bisection algorithms, among others\) face computational limitations in high\-dimensional systems and require*a priori*knowledge of the dynamical system and its equations\. Simplistic sampling methods such as random or uniform sampling of the phase space typically fail to quantitatively approximate separatrices and their structure altogether, suffering increasingly from the curse of dimensionality as the system dimension grows\.

We introduce and implement a framework that combinessupervised classificationwithgenerative modelingto address this challenge\. Our approach first trains neural network classifiers on uniformly or randomly sampled initial conditions, importantly labeled by their corresponding basins of attraction in the system of interest\. Using uncertainty metrics of the trained classifier to quantify decision boundaries, the method then identifies these high uncertainty regions and boundaries of the classifier as preliminaryapproximate separatrices\. Subsequently, score\-based generative models are trained specifically on samples from high\-uncertainty regions, ultimately generating densities of samples consistent with the empirical density of samples on the manifold \(or regions close to this manifold\) that constitutes the separatrix between basins in the sampled region\. This approach leverages the complementary strengths of \(a\) discriminative models for global phase space partitioning and \(b\) generative models for detailed geometric sampling, resulting in a systematic, iterative, data\-driven framework that produces empirically consistent reconstructions of \(approximate\) separatrix manifolds\.

*Keywords*Generative Models⋅\\cdotUncertainty Quantification⋅\\cdotDynamical Systems

## 1Introduction

A valuable step in the modeling of multistable systems is the characterization of the boundaries of their basins of attraction, commonly known as*separatrices*\. These separatrices often possess intricate and complex geometric structures that may even be fractal in nature\. Typically, these separatrices correspond to stable manifolds of saddle\-type invariant sets, unstable \(source\-type\) invariant sets such as limit cycles, or other complex manifolds\. In multistable dynamical systems\[[14](https://arxiv.org/html/2608.14743#bib.bib25),[28](https://arxiv.org/html/2608.14743#bib.bib24)\]trajectories initialized from arbitrary, generic points in phase space tend to converge to one of several attracting sets—fixed points, limit cycles, strange attractors, or other attracting invariant manifolds\. The basin of attraction for each attractor comprises the set of initial conditions whose forward\-time trajectories asymptotically approach that attractor\[[8](https://arxiv.org/html/2608.14743#bib.bib26)\]\. Long\-time numerical simulations of said systems often experience difficulty in sampling the phase space of these systems in sufficient detail because their trajectories become effectively trapped in neighborhoods of attractors, yielding statistically stationary behavior that provides information only about the local dynamics within a single basin\.

Characterizing separatrices presents significant computational and theoretical difficulties that have motivated substantial research in computational dynamics and adjacent fields\. In principle, if one has a dense coverage of estimated basins via emergent techniques such as the*recurrences\-based*\[[9](https://arxiv.org/html/2608.14743#bib.bib27)\]or*featurize\-and\-group*\[[7](https://arxiv.org/html/2608.14743#bib.bib28)\]methods, identifying the separatrix may be tractable\. This however becomes quickly infeasible as system dimension increases, as the cost of a dense covering of the basins scales exponentially\. Sampling points near the separatrices in order to enhance the time computational orbits spent near them is equally challenging, as separatrices typically correspond to state space regions of low probability density and often have formally zero volume\. This leads to the problem that the regions that govern transitions between basins and determine the system’s response to perturbations are those least visited by typical trajectories\.

Recent advances in applied mathematics have yielded several complementary approaches for identifying and approximating separatrices in high\-dimensional dynamical systems\. These include techniques such as manifold continuation\[[24](https://arxiv.org/html/2608.14743#bib.bib11),[2](https://arxiv.org/html/2608.14743#bib.bib33),[22](https://arxiv.org/html/2608.14743#bib.bib32),[20](https://arxiv.org/html/2608.14743#bib.bib31)\], edge tracking and bisection algorithms\[[30](https://arxiv.org/html/2608.14743#bib.bib12)\], or Lyapunov exponent Monte Carlo sampling\[[1](https://arxiv.org/html/2608.14743#bib.bib30)\]\. The application of these and more rudimentary techniques spans diverse fields including chemical kinetics and reactor design \(identifying steady states and bifurcations in chemical reactions\)\[[34](https://arxiv.org/html/2608.14743#bib.bib10),[11](https://arxiv.org/html/2608.14743#bib.bib9)\], ecology \(regime shifts in population dynamics and ecosystems\)\[[29](https://arxiv.org/html/2608.14743#bib.bib8)\], climate science\[[13](https://arxiv.org/html/2608.14743#bib.bib13)\], and even comprehensive, extremely high dimensional climate models\[[3](https://arxiv.org/html/2608.14743#bib.bib29)\]\. Reconstructing separatrices provides not only a characterization of critical transition thresholds but also generates representative samples from regions of phase space that are otherwise markedly undersampled by direct numerical simulation, while the prior mentioned methods such as continuation require explicit knowledge of equations\.

In parallel, the machine learning community has developed powerful generative modeling frameworks capable of learning to sample from complex probability distributions from finite sample sets and subsequently generating new samples that reflect the statistical properties of the training data, with early work involving models such as variational autoencoders and generative adversarial networks\[[23](https://arxiv.org/html/2608.14743#bib.bib6),[18](https://arxiv.org/html/2608.14743#bib.bib7)\]\. Score\-based generative models \(SGMs\) represent a particularly promising class of deep generative models\[[33](https://arxiv.org/html/2608.14743#bib.bib1),[32](https://arxiv.org/html/2608.14743#bib.bib2)\]\. These models learn the score function∇x​log​p​\(x\)\\nabla\_\{x\}\\log p\(x\), wherep⁡\(x\)p\(x\)is the target data distribution, and use this learned score to generate samples via Langevin dynamics or by solving a reverse\-time stochastic differential equation\. Additionally, these architectures have been applied to climate modeling problems specifically targeting the generation of samples near tipping points and separatrix regions\[[31](https://arxiv.org/html/2608.14743#bib.bib5)\], demonstrating that deep generative models can be specialized to focus on dynamically critical regions of phase space\. Crucially, recent work has identified connections between generative models such as SGMs and classical techniques in statistical physics and enhanced sampling\[[4](https://arxiv.org/html/2608.14743#bib.bib3),[6](https://arxiv.org/html/2608.14743#bib.bib14)\]\. Furthermore, SGMs and other generative models have been shown to generate density\-consistent data on manifolds\[[16](https://arxiv.org/html/2608.14743#bib.bib15),[5](https://arxiv.org/html/2608.14743#bib.bib16)\], and have demonstrated state\-of\-the\-art performance in high\-dimensional generative tasks, and they exhibit a remarkable ability to recover low\-dimensional manifold structure embedded in high\-dimensional ambient spaces\[[27](https://arxiv.org/html/2608.14743#bib.bib17)\]\. This has enabled hybrid approaches that combine the representational power of deep neural networks with physics\-informed constraints and sampling strategies\.

As a continuation/extension of these efforts, this work presents a framework that combines \(a\) supervised classification of basins of attraction with \(b\) sampling of densities on a manifold performed by a generative model, to achieve comprehensive characterization and reconstruction of separatrices in multistable dynamical systems\. The methodology proceeds as follows:

1. 1\.Basin classification:Train neural network classifiers to predict basin membership from provided phase space coordinates with known \(labeled\) basins\.
2. 2\.Decision boundary identification:Extract the decision boundaries of the trained classifier with the help of uncertainty metrics; these boundaries will approximate the true separatrices under appropriate conditions\.
3. 3\.Generative modeling on boundaries:Train score\-based or other generative models specifically on data of high uncertainty, which lie near the identified decision boundaries, enabling targeted sampling of separatrix regions\.
4. 4\.Refinement and visualization:Use generated samples to refine separatrix approximations iteratively and produce high\-fidelity visual representations of basin boundary geometry\.

This integrated approach leverages the complementary strengths of discriminative and generative models: classification provides global phase space partitioning and basin boundary localization, while generative models guide \(and thus enhance\) detailed sampling and reconstruction of the geometric structures comprising those boundaries\. The framework is particularly well\-suited to high\-dimensional systems where traditional continuation or bisection methods become computationally prohibitive, and provides a data\-driven \(and density\-consistent with the empirical data\) pathway to extracting geometric structures of separatrices from simulation or experimental time series data\.

## 2Proposed Mathematical Framework

Consider an autonomous dynamical system inℝn\\mathbb\{R\}^\{n\}governed by a flowΦt​\(𝐱\)=𝐱​\(t\)\\Phi^\{t\}\(\\mathbf\{x\}\)=\\mathbf\{x\}\(t\)\. Assume that this system containskkdistinct attractorsA1,A2,…,AkA\_\{1\},A\_\{2\},\\ldots,A\_\{k\}with corresponding basins of attractionB1,B2,…,BkB\_\{1\},B\_\{2\},\\ldots,B\_\{k\}\. The basin boundaries∂Bi\\partial B\_\{i\}are\(n−1\)\(n\-1\)\-dimensional sets \(generically\) that coincide with stable manifoldsWs​\(S\)W^\{s\}\(S\)of saddle\-type invariant setsSS\. These saddle points, unstable periodic orbits, or chaotic saddles play a critical role in organizing the phase space topology and mediating transitions between basins under perturbation\.

Our proposed approach uses a classification neural network to classify \(assign the correct\) basins of attractionBiB\_\{i\}for given initial conditions\. Beyond the assignment of the right basin for a given point, the classifier also reportsan uncertainty metric\(i\.e\. how sure the network is that the predicted basin is correct\)\. The uncertainty metric used in this work is based on Monte Carlo dropout\[[15](https://arxiv.org/html/2608.14743#bib.bib4)\], which is represented by an entropy\. This entropy,HH, is a measure of how much the network’s predictions vary when it is run multiple times with different random subsets of neurons temporarily switched off \(dropout\); if the predictions consistently agree, H is low and the classifier is confident, but if they disagree, H is high, signaling that the point likely lies in an uncertain region—such as near a boundary between basins\.

A generative model is then given initial conditions as training data with the uncertainty metric as a label; this then allows for the generation of new points at high uncertainty phase space areas, which will naturally lie at or near boundaries between the basins of the system\. These high uncertainty points function as proxy sampling points for local and global separatrices and/or their surrounding regions of phase space, depending on sampling density and model parameters\. Algorithms 1 and 2 outline the framework\.

Algorithm 1Basin Classifier Network Training and InferenceInput:Training data𝒟1=\{\(x​\(t0\)i,Bi\)\}i=1N\\mathcal\{D\}\_\{1\}=\\\{\(x\(t\_\{0\}\)\_\{i\},B\_\{i\}\)\\\}\_\{i=1\}^\{N\}\(initial conditions from a dynamical system,x⁡\(t0\)x\(t\_\{0\}\), and discrete labels given by the initial condition’s basin of attraction,BiB\_\{i\}that is reached from each initial condition after some timett\), neural networkOPENf⁡\(𝒟1\);θ\)f\(\\mathcal\{D\}\_\{1\}\);\\theta\)with dropout layers, loss functionℒ\\mathcal\{L\}, number of epochsEE, batch sizeBB, number of MC samplesMM

1:— Training Phase —

2:forepoch

=1=1to

EEdo

3:Shuffle training data

4:foreach batch

\{\(xb,yb\)\}b=1B\\\{\(x\_\{b\},y\_\{b\}\)\\\}\_\{b=1\}^\{B\}do

5:Enable dropout

6:Compute predictions:

y^b=f⁡\(xb,θ\)\\hat\{y\}\_\{b\}=f\(x\_\{b\};\\theta\)
7:Compute loss:

ℒbatch=ℒ⁡\(y^b,yb\)\\mathcal\{L\}\_\{\\text\{batch\}\}=\\mathcal\{L\}\(\\hat\{y\}\_\{b\},y\_\{b\}\)
8:Update parameters:

θ←θ−η​∇θℒbatch\\theta\\leftarrow\\theta\-\\eta\\nabla\_\{\\theta\}\\mathcal\{L\}\_\{\\text\{batch\}\}
9:endfor

10:endfor

11:— Inference Phase \(Monte Carlo Dropout\) —

12:New input sample

x⁡\(t0\)x\(t\_\{0\}\)
13:Initialize average probabilities,

p¯\\bar\{p\}as a zero vector

14:for

m=1m=1to

MMdo

15:Enable dropout at inference

16:Run forward pass of network:

y^m\\hat\{y\}\_\{m\}←\\leftarrowf⁡\(x​\(t0\)m\)f\(x\(t\_\{0\}\)\_\{m\}\)
17:Apply softmax to get probabilities:

pmp\_\{m\}←\\leftarrowsoftmax\(

y^m\\hat\{y\}\_\{m\}\)

18:Get average probabilities:

p¯←p¯\+1M​pm\\bar\{p\}\\leftarrow\\bar\{p\}\+\\frac\{1\}\{M\}p\_\{m\}
19:endfor

20:Compute uncertainty \(entropy\):

H\(p¯\)=−∑mp¯klogp¯kH\(\\bar\{p\}\)=\-\\sum\_\{m\}\\bar\{p\}\_\{k\}\\log\\bar\{p\}\_\{k\}
21:return

y^,p¯,H⁡\(p¯\)\\hat\{y\},\\bar\{p\},H\(\\bar\{p\}\)

Algorithm 2Training a score\-based generative model \(SGM\) to approximate the separatricesInput:Training data𝒟2=\{\(x​\(t0\)j\)\}j=1N\\mathcal\{D\}\_\{2\}=\\\{\(x\(t\_\{0\}\)\_\{j\}\)\\\}\_\{j=1\}^\{N\}\(initial condition samples from a dynamical system,x⁡\(t0\)x\(t\_\{0\}\), with corresponding distributionP⁡\(𝒟2\)P\(\\mathcal\{D\}\_\{2\}\), where𝒟2\\mathcal\{D\}\_\{2\}is chosen by selecting an entropy threshold from the Monte Carlo dropout of the classifier network,H⁡\(p¯CLOSEH\(\\bar\{p\}\)\. For practical notation purposes we will call the uncertainty thresholdzz\.

1:Retain high entropy

𝒟2\\mathcal\{D\}\_\{2\}data from Algorithm[1](https://arxiv.org/html/2608.14743#alg1)\. In this work, all initial conditions with

zzabove the 90th percentile were retained\.

2:Train a SGM on

𝒟2\\mathcal\{D\}\_\{2\}to model the distributions of high entropy samples

3:Use the trained SGM to sample from

P^​\(𝒟2\)\\hat\{P\}\(\\mathcal\{D\}\_\{2\}\)— approximated measures in ambient space at the prescribed high entropy

4:return

P^​\(𝒟2\)\\hat\{P\}\(\\mathcal\{D\}\_\{2\}\)for high

zz

The identification of separatrices with high\-uncertainty regions from the trained classifier represents an approximation whose validity depends critically on the relationship between the classifier’s learned decision boundary and the true dynamical basin boundary\. Under ideal conditions—sufficient training data uniformly distributed across phase space, adequate network capacity, and successful optimization—the classifier learns to assign high probability to the correct basin membership; classification uncertainty \(quantified via entropy of the MC dropout predictive distribution\) should peak near the true separatrix∂Bi\\partial B\_\{i\}, where infinitesimal changes in initial conditions lead to different asymptotic outcomes\. However, several sources of systematic bias affect this approximation\. First, the classifier decision boundary depends on the training data distribution\. Regions of phase space with sparse coverage of sampled initial conditions will exhibit elevated uncertainty regardless of proximity to the true separatrix, since the network has insufficient information to make confident predictions\. This creates a fundamental confound: high uncertainty may indicate either proximity to a basin boundary or absence of training data\. If the model is not trained on sufficiently numerous examples near the true decision boundary, and that boundary is not geometrically simple, the model may confidently misclassify instances by erroneously extrapolating learned patterns beyond their valid range\.

Furthermore, the geometry of the classifier decision boundary is constrained by the neural network’s inductive bias\. Standard feedforward networks with smooth activation functions \(e\.g\., ReLU, tanh\) produce piecewise smooth decision boundaries whose complexity scales with network depth and width\. For separatrices with fractal structure or those involving fine\-scale geometric detail, finite networks necessarily produce smoothed approximations\. The decision boundary may capture large\-scale topology while missing fine structure below the effective resolution determined by network architecture and training procedure\. Additionally, training data imbalance between basins systematically biases decision boundaries\. If basinB1B\_\{1\}is overrepresented relative to basinB2B\_\{2\}in the training set, standard cross\-entropy loss could incentivize the classifier to expand the predicted region forB1B\_\{1\}at the expense ofB2B\_\{2\}, shifting the decision boundary away from its dynamically correct location\. This is particularly problematic when the dynamics itself creates sampling bias\.

The approximation is most reliable when: \(1\) training trajectories provide relatively uniform coverage of phase space, perhaps through strategic initial condition selection, \(2\) sufficient trajectory data exists near the separatrix, which can be achieved iteratively by generating initial conditions \(and subsequent trajectories\) from high\-uncertainty regions, \(3\) the network architecture is sufficiently expressive to capture the separatrix geometry at the scale of interest, and \(4\) class balance is maintained or explicitly corrected through loss function weighting\. These considerations motivate the iterative refinement approach and the integration of generative models: by producing samples near/along provisional boundary estimates, the generative model enables targeted generation of new trajectories that reduce sampling bias, while validation against dynamical properties \(e\.g\., checking that generated points lie near saddle points/manifolds or exhibit long transient times\) provides quality control on the reconstruction\. To this point, the entropy threshold selected can offer some fine control on how "wide or narrow" the end result approximation of the separatrix becomes\. In this work, a threshold of ninety percent was decided to be quantitative and qualitatively sufficient for the examples that follow\.

## 3Numerical Examples

We demonstrate the proposed framework on three dynamical systems that exhibit multistability with progressively complex basin boundary structures\. The Newton method fractal boundary provides a two\-dimensional benchmark with intricate, self\-similar basin boundaries whose fractal geometry is well\-characterized, allowing us to assess the method’s ability to reconstruct fine\-scale geometric detail\. The continuous stirred tank reactor \(CSTR\) model represents a canonical problem in chemical process engineering where separatrices delineate operationally critical boundaries between stable regimes, and where traditional continuation methods have been extensively applied\. Finally, the three\-dimensional Lorenz system exemplifies chaotic dynamics where basin boundaries, formed by two\-dimensional stable manifolds of saddle equilibria, separate trajectories leading to qualitatively different asymptotic behaviors, possibly on a strange attractor\. Together, these examples demonstrate the framework’s capability to handle systems with varying dimensionality, geometric complexity, and practical significance, while demonstrating advantages over traditional methods by sidestepping computationally intensive integration/continuation and avoiding any required*a priori*knowledge of the dynamical system\.

### 3\.1The Newton Method Fractal Boundary

The Newton fractal is generated by applying Newton’s method for root finding to the polynomialf⁡\(z\)=z3−1f\(z\)=z^\{3\}\-1in the complex plane\. The iterative map is given by

zn\+1=zn−f⁡\(zn\)f′​\(zn\)=zn−zn3−13​zn2=2​zn3\+13​zn2,z\_\{n\+1\}=z\_\{n\}\-\\frac\{f\(z\_\{n\}\)\}\{f^\{\\prime\}\(z\_\{n\}\)\}=z\_\{n\}\-\\frac\{z\_\{n\}^\{3\}\-1\}\{3z\_\{n\}^\{2\}\}=\\frac\{2z\_\{n\}^\{3\}\+1\}\{3z\_\{n\}^\{2\}\},\(1\)
wherez∈ℂz\\in\\mathbb\{C\}\. Each initial conditionz0∈ℂz\_\{0\}\\in\\mathbb\{C\}converges to one of three roots:r1=1r\_\{1\}=1,r2=ei​2​π/3r\_\{2\}=e^\{i2\\pi/3\}, andr3=ei​4​π/3r\_\{3\}=e^\{i4\\pi/3\}\. The basins of attraction for these roots are separated by fractal boundaries\. This system provides an interesting test case for validating separatrix reconstruction accuracy due to its well\-studied geometric properties and the self\-similar structure of its basin boundaries across multiple spatial scales due to its fractal nature\. Figure[1](https://arxiv.org/html/2608.14743#S3.F1)demonstrates the algorithmic framework presented in section[2](https://arxiv.org/html/2608.14743#S2)in three steps: basin classification, using UQ metrics to identify the classifier decision boundary given by the high entropy initial conditions, and generative sampling of the separatrix by an SGM trained on the high entropy initial conditions\. Forty thousand uniformly sampled initial conditions \(a 200x200 mesh grid\) on the complex plane were integrated 50 steps each using Newton’s method using equation[1](https://arxiv.org/html/2608.14743#S3.E1)\. This grid of initial conditions labeled by the endpoints of their trajectories were used to train the neural network classifier\. Once trained, the classifier evaluated a new finer mesh grid of 160,000 \(400x400\) samples\. The top ten percent of the evaluated samples with respect to their entropy predicted by the classifier \(16,000 samples\) were retained for SGM training\. Once trained on the high entropy points, the SGM then generated 16,000 new samples for comparison with the high entropy training dataset that effectively approximates the decision boundary \(and thus the separatrix\) of the system\.

Once the framework was implemented and executed, the marginal densities of the high entropy points corresponding to the decision boundary of the classifier were compared with the marginal densities of the SGM generated samples\. Two metrics were used to quantitatively compare the two sets: Wasserstein distance, used generally to measure the "work" required to transform one density into another\[[17](https://arxiv.org/html/2608.14743#bib.bib18)\], and the Chamfer distance, the average of nearest\-neighbor distances from each point in one set to the closest point in the other set\[[12](https://arxiv.org/html/2608.14743#bib.bib19)\]\. These quantities effectively measure the difference in density and difference in coverage of the two sets, respectively, and the averages of the distances between the two marginals \(real and complex\) are reported in Table[1](https://arxiv.org/html/2608.14743#S3.T1)\. Figure[2](https://arxiv.org/html/2608.14743#S3.F2)displays the marginal distributions plotted against each other and qualitatively compared\.

![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/diffusion_vs_high_entropy_main.png)Figure 1:Newton fractal basin classification for the complex polynomialf⁡\(z\)=z3−1f\(z\)=z^\{3\}\-1, illustrating the uncertainty quantification and generative sampling framework across four subfigures\.Top\-left:uniformly sampled initial conditions in the complex plane, labeled and colored by their basin of attraction under Newton’s method iterationzn\+1=zn−f⁡\(zn\)/f′​\(zn\)z\_\{n\+1\}=z\_\{n\}\-f\(z\_\{n\}\)/f^\{\\prime\}\(z\_\{n\}\)\.Top\-right:high\-entropy points identified by the UQ metric derived from the classifier neural network, highlighting regions of fractal boundary ambiguity where basin membership is uncertain\.Bottom\-left:generated samples produced by a score\-based generative model \(SGM\) trained exclusively on the high\-entropy initial conditions, learning to reproduce the complex boundary structure\.Bottom\-right:overlay comparison in the complex plane contrasting the true high\-entropy data against the SGM\-generated samples, demonstrating the model’s ability to capture the fractal geometry near basin boundaries\.![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/diffusion_marginal_distributions.png)Figure 2:For the Newton fractal, marginal distributions in histogram form of \(a\) the high entropy samples corresponding to the decision boundary of the classifier plotted along with \(b\) the marginal distributions of the SGM generated samples\.Table 1:Marginal Distribution Comparison Metrics for the Newton Fractal
### 3\.2Sensitivity Analysis of the Entropy Threshold with a Continuous Stirred Tank Reactor \(CSTR\) Example

A critical consideration in the proposed framework is the selection of the entropy percentile threshold used to identify points near the separatrix\. The 90th percentile threshold is chosen in the analysis throughout this work as a reasonable illustration, but depending on the desired quantity and accuracy of the output, other threshold values may be warranted\. We performed a sensitivity analysis to examine how variations in the entropy threshold affect the extracted high\-entropy point set and to identify regions of stability in the threshold parameter space\.

Consider the non\-isothermal continuous stirred tank reactor \(CSTR\) model studied by Uppal, Ray, and Poore\[[34](https://arxiv.org/html/2608.14743#bib.bib10)\], given by the functional dimensionless equations

d​x1d​t\\displaystyle\\frac\{dx\_\{1\}\}\{dt\}=−x1\+D​a​\(1−x1\)​exp⁡\(x2\)=f1​\(x1,x2\)\\displaystyle=\-x\_\{1\}\+Da\(1\-x\_\{1\}\)\\exp\(\{x\_\{2\}\}\)=f\_\{1\}\(x\_\{1\},x\_\{2\}\)\(2\)d​x2d​t\\displaystyle\\frac\{dx\_\{2\}\}\{dt\}=−x2\+B​D​a​\(1−x1\)​exp⁡\(x2\)−β⁡\(x2\)=f2​\(x1,x2\)\\displaystyle=\-x\_\{2\}\+BDa\(1\-x\_\{1\}\)\\exp\(\{x\_\{2\}\}\)\-\\beta\(x\_\{2\}\)=f\_\{2\}\(x\_\{1\},x\_\{2\}\)\(3\)which exhibit multiple steady states arising from the interplay between reaction kinetics dependent on temperature, concentration, and other factors\. Characterizing this separatrix is important for reactor design and control, as it can delineate safe operating regions and determine startup trajectories that avoid extinction or runaway reactions\. We investigate this system at the parameter valuesB,D​a,β=16\.2,0\.12823,3\.0B,Da,\\beta=16\.2,0\.12823,3\.0, which exhibit a stable steady state surrounded by an unstable limit cycle \(the separatrix\), surrounded by a larger stable limit cycle\. Following the framework, 40,000 uniformly sampled initial conditions in the range ofx1,x2=\[0,1\],\[0,10\]\.x\_\{1\},x\_\{2\}=\[0,1\],\[0,10\]\.The initial conditions were then integrated using equations[2](https://arxiv.org/html/2608.14743#S3.E2)and[3](https://arxiv.org/html/2608.14743#S3.E3), and because the stable steady state lies close to the unstable limit cycle, the basins of attraction were labeled bytime to convergencerather than location in phase space\. Small perturbations near the steady state can lead to qualitatively different transient behaviors, despite originating from nearly identical regions in phase space, making this separatrix vanishingly small and difficult for the network to identify when using the traditional labels\. This grid of initial conditions labeled by basin based on time to convergence of their trajectories were used to train the neural network classifier\. Once trained, the classifier evaluated a new finer mesh grid of 160,000 \(400x400\) samples\.

The sensitivity analysis evaluates thresholds ranging from the 75th to 99th percentile of the entropy distribution\. For each candidate threshold, we compute several diagnostic quantities: the number of points exceeding the threshold, the spatial coverage of the resulting point set as a fraction of the domain, the centroid location of the extracted points, and the spatial spread characterized by the standard deviation in each coordinate direction\. From these quantities, we derive a combined instability score that measures the sensitivity of the results to small perturbations in the threshold\. Specifically, the instability score incorporates the normalized rate of change in point count with respect to the percentile threshold and the magnitude of centroid displacement as the threshold varies\. Regions where this score is minimized correspond to threshold values where the extracted separatrix approximation is most stable with respect to the threshold hyperparameter\. Mathematically, the score is defined as follows:

LetN⁡\(p\)N\(p\)denote the number of grid points whose entropy exceeds thepp\-th percentile threshold, and letμx1​\(p\)\\mu\_\{x\_\{1\}\}\(p\),μx2​\(p\)\\mu\_\{x\_\{2\}\}\(p\)denote the coordinates of the centroid of that high\-entropy point set, each regarded as continuous functions of the percentilepp\. The combined instability scoreS⁡\(p\)S\(p\), which quantifies how sensitive the detected separatrix is to the choice of threshold, is defined as:

S⁡\(p\)=\|d​N​\(p\)d​p\|maxp⁡\|d​N​\(p\)d​p\|\+\(d​μx1​\(p\)d​p\)2\+\(d​μx2​\(p\)d​p\)2maxp⁡\(d​μx1​\(p\)d​p\)2\+\(d​μx2​\(p\)d​p\)2S\(p\)\\;=\\;\\frac\{\\left\|\\dfrac\{dN\(p\)\}\{dp\}\\right\|\}\{\\displaystyle\\max\_\{p\}\\left\|\\dfrac\{dN\(p\)\}\{dp\}\\right\|\}\\;\+\\;\\frac\{\\sqrt\{\\left\(\\dfrac\{d\\mu\_\{x\_\{1\}\}\(p\)\}\{dp\}\\right\)^\{2\}\+\\left\(\\dfrac\{d\\mu\_\{x\_\{2\}\}\(p\)\}\{dp\}\\right\)^\{2\}\}\}\{\\displaystyle\\max\_\{p\}\\sqrt\{\\left\(\\dfrac\{d\\mu\_\{x\_\{1\}\}\(p\)\}\{dp\}\\right\)^\{2\}\+\\left\(\\dfrac\{d\\mu\_\{x\_\{2\}\}\(p\)\}\{dp\}\\right\)^\{2\}\}\}\(4\)
where the first term captures the normalized rate of change of the number of detected high\-entropy points with respect to the threshold percentile \(sensitivity of the point count\), and the second term captures the normalized magnitude of the rate of change of the centroid position \(sensitivity of the spatial location of the detected region\)\. Both derivatives are normalized by their respective maxima over the sampled percentile range so that each term lies in\[0,1\]\[0,1\], and the two are summed to produce a single scalar measure of overall instability at each thresholdpp\. Lower values ofS⁡\(p\)S\(p\)indicate percentile choices for which the separatrix detection is more robust to small perturbations in the threshold\.

Figure I[4](https://arxiv.org/html/2608.14743#S3.F4)presents the results of this analysis\. The top row displays the high\-entropy regions identified at three representative threshold values \(85th, 90th, and 95th percentiles\), overlaid on the entropy field computed by the MC Dropout classifier\. As the threshold increases, the extracted region contracts to a tighter estimation of the separatrix where classification uncertainty is highest\. The bottom panel shows the instability score as a function of the percentile threshold\. The vertical dashed lines indicate the thresholds corresponding to the top row visualizations, while the green dashed line marks the most stable threshold identified by the analysis\. The 97th percentile of points gives the optimal instability score, indicating that it would be the appropriate threshold for the most exact point\-wise approximation of the separatrix in this example, trading off less SGM training data for estimation accuracy\.

Maintaining consistency with the rest of this work, rather than optimizing each example’s cutoff distance, we demonstrate the method with a uniform entropy cutoff\. The top ten percent of the evaluated samples with respect to their entropy predicted by the classifier \(16,000 samples\) were retained for SGM training\. Once trained on the high entropy points, the SGM then generated 16,000 new samples for comparison with the high entropy training dataset that effectively approximates the decision boundary \(and thus the separatrix\) of the system\. Figure[3](https://arxiv.org/html/2608.14743#S3.F3)displays the training data given by the classifier, and Figure[5](https://arxiv.org/html/2608.14743#S3.F5)displays the scatter plot of the high entropy points and SGM generated samples as well as their marginal distributions\. The Wasserstein and Chamfer distance metrics were again computed to quantify the difference in density and difference in coverage of the two sets and the averages of the distances between the two marginals \(real and complex\)\. Their values are0\.05720\.0572and0\.00880\.0088, respectively\.

![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/cstr_data_overview.png)Figure 3:Uniformly sampled initial conditions of the CSTR model equations defined by equations[2](https://arxiv.org/html/2608.14743#S3.E2)and[3](https://arxiv.org/html/2608.14743#S3.E3)colored by binary label of the corresponding attractor identified by time to convergence \(purple for the stable steady state and yellow for the stable limit cycle\)\.![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/cstr_threshold_sensitivity_combined.png)Figure 4:Sensitivity analysis of the entropy threshold selection for the CSTR system\.Top row:High\-entropy regions identified at the 85th, 90th, and 95th percentile thresholds of the instability score defined by equation[4](https://arxiv.org/html/2608.14743#S3.E4)colored as purple, cyan, and light green, respectively, plotted with the entropy field computed by the MC Dropout classifier shaded from dark red to white\.Bottom:Instability score as a function of the percentile threshold; lower values indicate greater stability with respect to threshold perturbations\. The green dashed line indicates the most stable threshold, while colored dotted lines correspond to the thresholds visualized in the top row\.![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/cstr_scatter_and_marginals.png)Figure 5:Comparison of high\-entropy points and the SGM\-generated samples for the CSTR system\.Left:Scatter plot overlay of the high\-entropy training points \(blue\) extracted from the MC Dropout classifier at the 90th percentile threshold and the samples \(red\) generated by the trained score\-based generative model\. The spatial agreement between the two distributions indicates that the SGM has successfully learned to sample from the separatrix\.Middle:Marginal density distributions along thex1x\_\{1\}coordinate \(dimensionless concentration\), showing close agreement between the high\-entropy and the generated samples\.Right:Marginal density distributions along thex2x\_\{2\}coordinate\. The strong overlap in both marginal distributions demonstrates that the generative model accurately captures the statistical structure of the classifier’s decision boundary corresponding to the dynamical separatrix\. The “peaks" in the high entropy marginals are a product of their uniform sampling\.
### 3\.3A Comparison of Methods with the Lorenz System

The Lorenz system, originally derived as a simplified model of atmospheric convection, is a three\-dimensional autonomous system given by

d​xd​t\\displaystyle\\frac\{dx\}\{dt\}=σ⁡\(y−x\),\\displaystyle=\\sigma\(y\-x\),\(5\)d​yd​t\\displaystyle\\frac\{dy\}\{dt\}=x⁡\(ρ−z\)−y,\\displaystyle=x\(\\rho\-z\)\-y,\(6\)d​zd​t\\displaystyle\\frac\{dz\}\{dt\}=x​y−β​z,\\displaystyle=xy\-\\beta z,\(7\)where\(x,y,z\)∈ℝ3\(x,y,z\)\\in\\mathbb\{R\}^\{3\}represent the system state andσ\\sigma,ρ\\rho, andβ\\betaare parameters\. We consider the parameter valuesσ=10\\sigma=10,ρ=20\\rho=20, andβ=8/3\\beta=8/3\. At these parameters, which lie below the critical valueρc≈24\.74\\rho\_\{c\}\\approx 24\.74for the onset of chaos, the system exhibits bistability between two stable fixed pointsC±=\(±β⁡\(ρ−1\),±β⁡\(ρ−1\),ρ−1\)C^\{\\pm\}=\(\\pm\\sqrt\{\\beta\(\\rho\-1\)\},\\pm\\sqrt\{\\beta\(\\rho\-1\)\},\\rho\-1\), while the origin remains an unstable saddle with a two\-dimensional stable manifold and a one\-dimensional unstable manifold\. The stable manifoldWs​\(𝟎\)W^\{s\}\(\\mathbf\{0\}\)of the origin forms the separatrix, a two\-dimensional surface that divides the phase space into two basins of attraction corresponding toC\+C^\{\+\}andC−C^\{\-\}\. This configuration provides an interesting test case for three\-dimensional separatrix reconstruction in a system with well\-defined multistability, a geometrically complex basin boundary structure, and chaotic transient behavior\. This separatrix is the 2D stable manifold of the origin of the 3D Lorenz system described by the above parameter values, which we will briefly call the "2D Lorenz Manifold" from now on\.

![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/Lorenz_metrics_comparison.png)Figure 6:For the Lorenz System, Comparison of high entropy points from three different methods \(MC dropout, deep ensemble averaging, and Laplace approximation\) versus a sampled manifold obtained via SCIGMA\. The three metrics used are Chamfer distance, coverage, and Hausdorff distance\.To assess the effectiveness of different uncertainty quantification approaches for separatrix identification, we compared three methods for estimating predictive uncertainty from the trained classifier: Monte Carlo \(MC\) dropout, deep ensemble averaging\[[25](https://arxiv.org/html/2608.14743#bib.bib22)\], and Laplace approximation\[[26](https://arxiv.org/html/2608.14743#bib.bib21)\]\. For each method, we identified high\-uncertainty regions and evaluated their geometric correspondence to the ground truth separatrix computed using the SCIGMA\[[21](https://arxiv.org/html/2608.14743#bib.bib23)\]continuation software\.111The ground truth separatrix was computed using SCIGMA, a software package for invariant\-manifold computation and visualization\[[21](https://arxiv.org/html/2608.14743#bib.bib23)\]\. We used the December 2025 macOS ARM release with themstablecontinuation commands, initialized from the unstable saddle at\(0,0,0\)\(0,0,0\)\. Numerical integration usedd​t=0\.005dt=0\.005; AUTO continuation usedd​s=0\.01ds=0\.01,d​smin=10−6ds\_\{\\min\}=10^\{\-6\}, andd​smax=1ds\_\{\\max\}=1, with all other parameters set to their defaults\.Performance was quantified using three metrics: coverage \(the fraction of true separatrix points within a specified distance of high\-uncertainty regions\), Chamfer distance \(measuring average nearest\-neighbor distances between predicted and true boundary points in both directions, described and cited in previous sections\), and Hausdorff distance \(capturing worst\-case deviations\)\[[19](https://arxiv.org/html/2608.14743#bib.bib20)\]\. MC dropout outperformed the alternatives in coverage and Chamfer distance while staying comparatively equal or better with the other methods in Hausdorff distance, as detailed in Figure[6](https://arxiv.org/html/2608.14743#S3.F6)\. The superior performance of MC dropout likely stems from its ability to capture local uncertainty variations more effectively than the global uncertainty estimates provided by methods such as Laplace approximation, while requiring substantially less computational cost than training full ensembles\.

![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/Lorenz_manifold_comparison_3d.png)Figure 7:Three\-dimensional plots of the approximated 2D Lorenz Manifold\. Left: An overlaid comparison of the generated approximation of the manifold using the generative method plotted as red points and the sampled manifold found via SCIGMA in blue points\. Right: the SCIGMA calculated manifold points, colored by distance to the approximated manifold \(effectively the Chamfer distance\)\.Using MC dropout as the uncertainty quantification method of choice, we applied the complete framework \(classification, boundary identification via high\-uncertainty regions, and generative modeling on boundary data\) to reconstruct the separatrix of the Lorenz system\. The final output from the score\-based generative model, comprising samples along the reconstructed boundary, was compared against the ground truth manifold\. The generative model substantially improved the quality of the reconstruction relative to the classifier alone: the Chamfer distance decreased from 3\.05 to 1\.58 \(a 48\.2% reduction\) and coverage increased from 96\.7% to 99%\. The Hausdorff distance improved slightly \(from 20\.9 to 20\.58\), consistent with the low variation seen between the UQ methods\. Qualitatively, the generated samples captured the fine\-scale geometric features of the separatrix\. While this is demonstrated in the improved metrics and not readily apparent from the three dimensional point clouds in Figure[7](https://arxiv.org/html/2608.14743#S3.F7), this can be seen when looking at two dimensional projections of the ground truth manifold versus the approximate manifold from the generative model, demonstrated in Figure[8](https://arxiv.org/html/2608.14743#S3.F8)\. These projections do not exactly approximate the true calculated separatrix due to the large sampled region and inherent complexity of the manifold, but with iterations on specific regions, they would approach an exact approximation\. These results demonstrate that the integrated classification\-generation approach produces separatrix reconstructions with fidelity comparable to specialized continuation methods, while offering advantages in computational efficiency and applicability to systems where continuation methods face challenges, such as when the equations to a system are not known*a priori*\.

![Refer to caption](https://arxiv.org/html/2608.14743v1/Figures/Lorenz_projections.png)Figure 8:Sets of two dimensional projections from various sections of the separatrix calculated using traditional methods \(plotted in blue\), versus the generated points that approximate the separatrix \(maked in red\)\. Inset metrics describing the Chamfer distance, coverage, and sliced Wasserstein distance of each projection are included in each plot\.

## 4Conclusion

This work has introduced a data\-driven framework for computationally approximating separatrices in multistable dynamical systems by combining supervised classification with generative modeling\. The methodology addresses fundamental sampling challenges inherent in characterizing basin boundaries: these structures govern transitions between stable states yet remain poorly sampled in standard numerical simulations due to their repelling behavior\. We have demonstrated the framework’s effectiveness on three representative dynamical systems\. Across these examples, the method successfully identified and reconstructed basin boundaries with fidelity comparable to traditional techniques while offering advantages in computational efficiency, applicability to high\-dimensional systems, and circumvention of required system knowledge\. The framework provides particular value for systems where traditional continuation and bisection methods face limitations\. By learning directly from trajectory data rather than requiring known equations, it accommodates systems known only through simulation or experiment\. The approach naturally extends to high\-dimensional phase spaces where visualization fails and manifold continuation becomes prohibitively expensive\. Moreover, the generative component addresses a fundamental need: producing representative samples from regions of phase space that standard simulations systematically undersample\. This framework demonstrates that modern machine learning, when appropriately integrated with dynamical systems theory, provides powerful capabilities for characterizing global geometric structures in complex nonlinear systems by enabling a data\-driven understanding of critical transitions that govern system behavior\.

### 4\.1Future Work

The future work will focus on augmenting this method with an iterative and active learning framework to mitigate its existing limitations\.Isolated high\-uncertainty regionscan emerge in sparsely sampled areas far from any true separatrix, leading to false positive identification of spurious boundaries\. These arise when the network, lacking local training data, defaults to uncertain predictions while still being far from any dynamically significant structure\. Several promising directions exist for extending this methodology and mitigating this limitation\. First, systematic evaluation of the framework’s performance in sparse data regimes would address a critical practical constraint\. Understanding how reconstruction quality degrades with reduced training data and identifying minimum data requirements for reliable separatrix identification would guide deployment in data\-limited contexts\. Developing uncertainty quantification measures that indicate when additional data is needed would further enhance practical utility, and some readily available metrics already exist, e\.g\., basin entropy\[[10](https://arxiv.org/html/2608.14743#bib.bib34)\]\. Additionally, active learning strategies offer substantial potential for improving sample efficiency\. Current implementation treats initial condition sampling as uniform or predetermined, but an adaptive approach could iteratively select the most informative regions for new trajectory generation\. In future work we plan to utilize functionality offered by the software DynamicalSystems\.jl, which provides generic functions for labeling trajectories according to their basin membership\[[9](https://arxiv.org/html/2608.14743#bib.bib27)\]\. This will allow for an iterative refinement of separatrices, iteratively probing high uncertainty regions with denser sampling\. Finally, integration with optimal experimental design principles could formalize selection criteria to maximize information gain about basin boundary geometry per simulation\. Further tests will be performed for the improved methods, including multi\-seed robustness tests and improvements on our existing threshold sensitivity studies\.

## References

- \[1\]\(2014\)A novel method to identify boundaries of basins of attraction in a dynamical system using lyapunov exponents and monte carlo techniques\.Nonlinear Dynamics79\(1\),pp\. 275–293\.External Links:ISSN 1573\-269X,[Link](http://dx.doi.org/10.1007/s11071-014-1663-z),[Document](https://dx.doi.org/10.1007/s11071-014-1663-z)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[2\]D\. G\. Aronson, M\. A\. Chory, G\. R\. Hall, and R\. P\. McGehee\(1982\)Bifurcations from an invariant circle for two\-parameter families of maps of the plane: a computer\-assisted study\.Communications in Mathematical Physics83\(3\),pp\. 303 – 354\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[3\]R\. Börner, O\. Mehling, J\. von Hardenberg, and V\. Lucarini\(2025\)Global stability of the atlantic overturning circulation: edge state, long transients and boundary crisis under co2\{\}\_\{2\}forcing\.arXiv\.External Links:[Document](https://dx.doi.org/10.48550/ARXIV.2504.20002),[Link](https://arxiv.org/abs/2504.20002)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[4\]E\. R\. Crabtree, J\. M\. Bello\-Rivas, A\. L\. Ferguson, and I\. G\. Kevrekidis\(2023\)Gans and closures: micro\-macro consistency in multiscale modeling\.Multiscale modeling & simulation21\(3\),pp\. 1122–1146\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[5\]E\. R\. Crabtree, D\. G\. Giovanis, N\. Evangelou, J\. M\. Bello\-Rivas, and I\. G\. Kevrekidis\(2025\)Generative learning for slow manifolds and bifurcation diagrams\.Computers & Chemical Engineering,pp\. 109544\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[6\]E\. Crabtree, J\. Bello\-Rivas, and I\. Kevrekidis\(2024\)Micro\-macro consistency in multiscale modeling: score\-based model assisted sampling of fast/slow dynamical systems\.Chaos: An Interdisciplinary Journal of Nonlinear Science34\(5\)\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[7\]G\. Datseris, K\. Luiz Rossi, and A\. Wagemakers\(2023\)Framework for global stability analysis of dynamical systems\.Chaos: An Interdisciplinary Journal of Nonlinear Science33\(7\),pp\. 073151\.External Links:ISSN 1054\-1500, 1089\-7682,[Document](https://dx.doi.org/10.1063/5.0159675)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p2.1)\.
- \[8\]G\. Datseris and U\. Parlitz\(2022\)Nonlinear Dynamics: A Concise Introduction Interlaced with Code\.Undergraduate Lecture Notes in Physics,Springer International Publishing,Cham\.External Links:ISBN 978\-3\-030\-91031\-0 978\-3\-030\-91032\-7,[Document](https://dx.doi.org/10.1007/978-3-030-91032-7)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p1.1)\.
- \[9\]G\. Datseris and A\. Wagemakers\(2022\)Effortless estimation of basins of attraction\.Chaos: An Interdisciplinary Journal of Nonlinear Science32\(2\),pp\. 023104\.External Links:ISSN 1054\-1500, 1089\-7682,[Document](https://dx.doi.org/10.1063/5.0076568)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p2.1),[§4\.1](https://arxiv.org/html/2608.14743#S4.SS1.p1.1)\.
- \[10\]A\. Daza, A\. Wagemakers, B\. Georgeot, D\. Guéry\-Odelin, and M\. A\. F\. Sanjuán\(2016\)Basin entropy: a new tool to analyze uncertainty in dynamical systems\.Scientific Reports6\(1\)\.External Links:ISSN 2045\-2322,[Link](http://dx.doi.org/10.1038/srep31416),[Document](https://dx.doi.org/10.1038/srep31416)Cited by:[§4\.1](https://arxiv.org/html/2608.14743#S4.SS1.p1.1)\.
- \[11\]S\.S\. Elnashaie, M\.E\. Abashar, and F\.A\. Teymour\(1993\)Bifurcation, instability and chaos in fluidized bed catalytic reactors with consecutive exothermic chemical reactions\.Chaos, Solitons and Fractals3\(1\),pp\. 1–33\.External Links:ISSN 0960\-0779,[Document](https://dx.doi.org/https%3A//doi.org/10.1016/0960-0779%2893%2990037-2),[Link](https://www.sciencedirect.com/science/article/pii/0960077993900372)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[12\]H\. Fan, H\. Su, and L\. J\. Guibas\(2017\)A point set generation network for 3d object reconstruction from a single image\.InProceedings of the IEEE conference on computer vision and pattern recognition,pp\. 605–613\.Cited by:[§3\.1](https://arxiv.org/html/2608.14743#S3.SS1.p4.1)\.
- \[13\]U\. Feudel and C\. Grebogi\(1997\)Multistability and the control of complexity\.Chaos: An Interdisciplinary Journal of Nonlinear Science7\(4\),pp\. 597–604\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[14\]U\. Feudel, A\. N\. Pisarchik, and K\. Showalter\(2018\)Multistability and tipping: From mathematics and physics to climate and brain—Minireview and preface to the focus issue\.Chaos: An Interdisciplinary Journal of Nonlinear Science28\(3\),pp\. 033501\.External Links:ISSN 1054\-1500, 1089\-7682,[Document](https://dx.doi.org/10.1063/1.5027718)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p1.1)\.
- \[15\]Y\. Gal and Z\. Ghahramani\(2016\)Dropout as a bayesian approximation: representing model uncertainty in deep learning\.Ininternational conference on machine learning,pp\. 1050–1059\.Cited by:[§2](https://arxiv.org/html/2608.14743#S2.p2.1)\.
- \[16\]D\. G\. Giovanis, E\. Crabtree, R\. G\. Ghanem, and I\. G\. Kevrekidis\(2025\)Generative learning of densities on manifolds\.Computer Methods in Applied Mechanics and Engineering446,pp\. 118266\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[17\]C\. R\. Givens and R\. M\. Shortt\(1984\)A class of wasserstein metrics for probability distributions\.\.Michigan Mathematical Journal31\(2\),pp\. 231–240\.Cited by:[§3\.1](https://arxiv.org/html/2608.14743#S3.SS1.p4.1)\.
- \[18\]I\. J\. Goodfellow, J\. Pouget\-Abadie, M\. Mirza, B\. Xu, D\. Warde\-Farley, S\. Ozair, A\. Courville, and Y\. Bengio\(2014\)Generative Adversarial Networks\.arXiv:1406\.2661 \[cs, stat\]\.Note:arXiv: 1406\.2661External Links:[Link](http://arxiv.org/abs/1406.2661)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[19\]D\. P\. Huttenlocher, G\. A\. Klanderman, and W\. J\. Rucklidge\(2002\)Comparing images using the hausdorff distance\.IEEE Transactions on pattern analysis and machine intelligence15\(9\),pp\. 850–863\.Cited by:[§3\.3](https://arxiv.org/html/2608.14743#S3.SS3.p2.1)\.
- \[20\]M\. E\. Johnson, M\. S\. Jolly, and I\. G\. Kevrekidis\(1997\)Two\-dimensional invariant manifolds and global bifurcations: some approximation and visualization studies\.Numerical Algorithms14\(1\),pp\. 125–140\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[21\]I\. G\. Kevrekidis, M\. S\. Jolly, M\. A\. Taylor, M\. E\. Johnson, and R\. Hölzel\(2025\)SCIGMA: stability computations and interactive graphics for invariant manifold analysis\.Note:[https://scigma\.org/](https://scigma.org/)Cited by:[§3\.3](https://arxiv.org/html/2608.14743#S3.SS3.p2.1),[footnote 1](https://arxiv.org/html/2608.14743#footnote1)\.
- \[22\]I\. G\. Kevrekidis\(1987\)A numerical study of global bifurcations in chemical dynamics\.AIChE Journal33\(11\),pp\. 1850–1864\.External Links:[Document](https://dx.doi.org/https%3A//doi.org/10.1002/aic.690331112),[Link](https://aiche.onlinelibrary.wiley.com/doi/abs/10.1002/aic.690331112),https://aiche\.onlinelibrary\.wiley\.com/doi/pdf/10\.1002/aic\.690331112Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[23\]D\. P\. Kingma and M\. Welling\(2013\)Auto\-encoding variational bayes\.arXiv\.External Links:[Document](https://dx.doi.org/10.48550/ARXIV.1312.6114),[Link](https://arxiv.org/abs/1312.6114)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[24\]B\. Krauskopf, H\. M\. Osinga, E\. J\. Doedel, M\. E\. Henderson, J\. Guckenheimer, A\. Vladimirsky, M\. Dellnitz, and O\. Junge\(2005\)A survey of methods for computing \(un\) stable manifolds of vector fields\.International Journal of Bifurcation and Chaos15\(03\),pp\. 763–791\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[25\]B\. Lakshminarayanan, A\. Pritzel, and C\. Blundell\(2017\)Simple and scalable predictive uncertainty estimation using deep ensembles\.Advances in neural information processing systems30\.Cited by:[§3\.3](https://arxiv.org/html/2608.14743#S3.SS3.p2.1)\.
- \[26\]D\. J\. MacKay\(1992\)A practical bayesian framework for backpropagation networks\.Neural computation4\(3\),pp\. 448–472\.Cited by:[§3\.3](https://arxiv.org/html/2608.14743#S3.SS3.p2.1)\.
- \[27\]J\. Pidstrigach\(2022\)Score\-based generative models detect manifolds\.arXiv\.External Links:[Document](https://dx.doi.org/10.48550/ARXIV.2206.01018),[Link](https://arxiv.org/abs/2206.01018)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[28\]A\. N\. Pisarchik and A\. E\. Hramov\(2022\)Multistability in physical and living systems\.Springer International Publishing\.External Links:[Document](https://dx.doi.org/10.1007/978-3-030-98396-3),ISBN 978\-3\-030\-98395\-6,[Link](https://link.springer.com/10.1007/978-3-030-98396-3)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p1.1)\.
- \[29\]M\. Scheffer, S\. Carpenter, J\. A\. Foley, C\. Folke, and B\. Walker\(2001\)Catastrophic shifts in ecosystems\.Nature413\(6856\),pp\. 591–596\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[30\]J\. D\. Skufca, J\. A\. Yorke, and B\. Eckhardt\(2006\)Edge of chaos in a parallel shear flow\.Physical review letters96\(17\),pp\. 174101\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1)\.
- \[31\]J\. Sleeman, D\. Chung, A\. Gnanadesikan, J\. Brett, Y\. Kevrekidis, M\. Hughes, T\. Haine, M\. Pradal, R\. Gelderloos, C\. Ashcraft,et al\.\(2023\)A generative adversarial network for climate tipping point discovery \(tip\-gan\)\.arXiv preprint arXiv:2302\.10274\.Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[32\]Y\. Song and S\. Ermon\(2020\)Improved techniques for training score\-based generative models\.arXiv\.External Links:[Document](https://dx.doi.org/10.48550/ARXIV.2006.09011),[Link](https://arxiv.org/abs/2006.09011)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[33\]Y\. Song, J\. Sohl\-Dickstein, D\. P\. Kingma, A\. Kumar, S\. Ermon, and B\. Poole\(2020\)Score\-based generative modeling through stochastic differential equations\.arXiv\.External Links:[Document](https://dx.doi.org/10.48550/ARXIV.2011.13456),[Link](https://arxiv.org/abs/2011.13456)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p4.1)\.
- \[34\]A\. Uppal, W\.H\. Ray, and A\.B\. Poore\(1974\)On the dynamic behavior of continuous stirred tank reactors\.Chemical Engineering Science29\(4\),pp\. 967–985\.External Links:ISSN 0009\-2509,[Document](https://dx.doi.org/https%3A//doi.org/10.1016/0009-2509%2874%2980089-8),[Link](https://www.sciencedirect.com/science/article/pii/0009250974800898)Cited by:[§1](https://arxiv.org/html/2608.14743#S1.p3.1),[§3\.2](https://arxiv.org/html/2608.14743#S3.SS2.p2.1)\.