Intermediate

Let us consider complex quadratic polynomials that can be expressed as 

$$\phi_c(z)=z^2+c.$$

The filled Julia set gathers $z$ such that the sequence of their iterates by $\phi_c$ is bounded. 

Fix some $R>0$ large enough that $R^2-R\geq |c|$, the filled Julia set for this system is the subset of the complex plane given by

$$K(\phi_c)=\{z\in C\; :\; \forall n \in N, \;|\phi_c^n(z)|\leq R\},$$ 

where $\phi_c^n(z)$ is the nth iterate of $\phi_c(z)$. The Julia set $J(\phi_c)$ of this function is the boundary of $K(\phi_c)$.

In 1982, a deep theorem by Adrien Douady and Hubbard states that Julia sets are path-connected if and only if the sequence of the iterates of $0$ is bounded. 

Moreover, $K(\phi_c)$ and $J(\phi_c)$ are either path-connected (there is a path within the set that connects any two given points in the set) or path-disconnected (for any two given points in the set, it is impossible to find a path within the set that connects the pair). The set of complex numbers $c$ that can form a path-connected Julia set is called the Mandelbrot set. It is defined in the complex plane as the complex numbers $c$ for which the function $\phi_c$ does not diverge to infinity when iterated starting at $z=0$.

The set Mandelbrot set is bounded, included in $B(0,2)$. Different components of the Mandelbrot set correspond to different dynamical behaviors of the corresponding filled Julia sets.

From

Aubin, J.-P., Bayen, A., & Saint-Pierre, P. Viability Theory: New Directions. Springer. 2011.

• Obviously, the subset $K_c:=Viab_{\phi_c}(B(0,2))$ is the filled Julia set for the function $\phi_c$ whenever $|c|\leq 2$ (which is true when $c$ belongs to the Mandelbrot set) and its boundary $J_u:=\partial K_c$ is the Julia set.

We can approximate these Julia sets thanks to the Viability Algorithm.

Two famous Julia sets : 

• Douady Rabbit when $c$ is a complex root of $c^3+2c^2+c+$,  $c\approx -0.123 \pm 0.745i $ which implies that the orbit of 0 is periodic with period 3. Douady Rabbit is then connected. See below an approximation obtained by the ViabLab software :

You can see below an approximation of the filled Julia and the Julia set for $c=???$ obtained using the VIABLAB software:.

This system describes a farm and a restaurant belonging to a same project. Consequently, they function in full cooperation.

Description of the system : 

de Lapparent, A., Martin, S. & Sabatier, R. Using System Modularity to Simplify Viability Studies: An Application to a Farm-Restaurant Interaction. Environ Model Assess (2024). https://doi.org/10.1007/s10666-024-10014-w

The object computed is a viability kernel. The model is discrete in states, controls and time.

The computation takes a moment (4088s on my computer), be patient...

Model

States and controls

State variables

upper limit of $x_1$ can be relaxed.

Control variables

Dynamics

Overall dynamics are:

\begin{equation}
   \mathcal{S}_U
   \begin{cases}
   x_{1}^{t+1} = x_{1}^t + G(x_{2}^t,u_{3}^t,R(x_3^t,u_1^t,u_2^t)) - E(u_{1}^t,u_{2}^t)\\
   x_{2}^{t+1} = \alpha(x_{2}^t,u_{3}^t,R(x_3^t,u_1^t,u_2^t))\\
   x_{3}^{t+1} = \Phi (x_{3}^t ,u_{1}^t,u_{2}^t) \\
   \end{cases}
\end{equation}

with the following functions:

Some dynamics require to use grid parameters. Consequently, a function has been implemented into the source file to get these values.

Constraints

There are two cconstraints in this system: the global system has to be profitable and a minimal soil quality has to be preserved in order to address sustainability concerns. These constraints take the form of thresholds on the cumulative cash flow ($x_{1} \geq x_{1min}$) and on soil quality ($x_3 \geq x_{3min}$), respectively. In other words, $(x_1^t,x_2^t,x_3^t)$ must remain in $K$ for all $t\in \mathbb{N}$ with : 
\begin{equation}
K:=\{(x_1,x_2,x_3)\in \mathbb{R}^+\times [0;1]^2 \; |\; x_1\geq x_{1min} \text{ and }x_3\geq x_{3min}\}.
\end{equation}

Implementation parameters

Time horizon

The time horizon (for trajectory computations) is 20 years.

Algorithm parameters

Default parameters are used.

System parameters

We used the parameters for a low-hypotheses computation.

   "SYSTEM_PARAMETERS": {
       "DYNAMICS_TYPE": 2,
       "DYN_BOUND": 1,
       "DYN_BOUND_COMPUTE_METHOD": 2,
       "IS_TIMESTEP_GLOBAL": 0,
       "LIPSCHITZ_CONSTANT": 1,
       "LIPSCHITZ_CONSTANT_COMPUTE_METHOD": 2,
       "TIME_DISCRETIZATION_SCHEME": 4
   }

Viability kernel computed using ViabLab

This model represents the effect of grazing multi-species herds on feed resources, as well as its effects on farmers' income and self-sufficiency. Due to the highly extensive nature of the herds considered, only culling allows for herd management (no breeding or supplemental feeding is practiced).

These pastoral systems are particularly susceptible to extreme bioclimatic hazards known as Dzuds, which will be central to the modeling approach we will undertake.

Dzuds are episodes of sudden and massive mortality in herds, generally occurring at the end of winter when stocking rates are too high relative to forage resources and climatic conditions are particularly harsh. This binary phenomenon therefore depends on both farmers' practices and climatic conditions, as it is the combination of these two factors that triggers it.

The model we have developed links the dynamics of herds of 5 herbivore species (cattle, sheep, goats, horses and camels) to the dynamics of the forage resource and posits the Dzud as a central mechanism in the dynamics of the system.

Methodologically, this case study allowed us to evaluate the three properties of robustness, adaptability, and resilience within a single case.

Applied to practical application, the long-term work we conducted on this case study enabled us to gain a detailed understanding of the system's dynamics. We were thus able to distinguish two viable operating regimes for the system: small herds that prevent the onset of Dzuds by preserving resources, and very large herds that can withstand Dzuds when they occur.

The detailed model can be found in the following article:
Rodolphe Sabatier, Frédéric Joly, Bernard Hubert. Assessing both ecological and engineering resilience of a steppe agroecosystem using the viability theory. Agricultural Systems, Elsevier Masson, 2017, 157, pp.146-156. ⟨10.1016/j.agsy.2017.07.009⟩. ⟨hal-02627849⟩

Introduction $k$

Le problème d'eutrophisation du lac décrit comme le problème de viabilité du lac et des exploitations riveraines est traité à un niveau global et suppose un décideur unique. Dans le présent exemple les parties prenantes sont réunies dans un comité, et les membres du comité (appelés individus dans la suite) ne sont pas nécessairement d'accord sur la dynamique du lac. Cet exemple illustre les cas où la gestion est collective et où il n'y a pas de consensus sur la dynamique. Ce problème est un problème de viabilité garantie.

Ce problème est décrit en détail dans :

Alvarez, I., Zaleski, L., Briot, J.-P., Irving, M. de A. (2023). Collective management of environmental commons with multiple usages: A guaranteed viability approach. Ecological Modelling, 475, Article 110186

Version auteur disponible ici

Les caractéristiques 

Les caractéristiques du problème modélisé sont les suivantes :  

1.  $N$ individus utilisent les même variables pour décrire le système. Dans cet exemple du lac, ce sont les mêmes variables $L$ (apports de phosphore) et $P$ (concentration totale de phosphore) que dans le cas du décideur unique.
2. Les contrôles admissibles retenus dans les modèles sont communs à tous les individus. On garde donc la même variable de contrôle et les mêmes bornes : $u\in U=[u_{min},u_{max}]$.
3. Tous les individus $i$ peuvent définir un ensemble d'états souhaitables $K_i$ dont l'intersection est non vide. On garde donc un ensemble de contraintes similaires au cas du décideur unique :  $x=\left(L(t),P(t)\right) \in K=[L_{min}, +\infty[ \times [0,P_{max}]$
4. La dynamique du système est modélisée pour chaque individu par un ensemble d'équations (ou d'inclusion différentielles) $S_i$
5. L'objectif de chaque individu est de maintenir le système modélisé par leur dynamique $S_i$ dans leur ensemble d'états souhaitables $K_i$

6. Les individus acceptent de partager leurs informations personnelles avec un tiers de confiance.

   Dans cet exemple, on considère $N=4$ individus. Trois d'entre eux ($i\in\{1,2,3\}$) adoptent la dynamique classique, mais avec des paramètres différents :

   \begin{equation}
   (P_i \text{ modèle 1})\left\{
   \begin{array}{l}
   \frac{dL}{dt}=u\in U=\left[ u_{min},u_{max}\right] \\
   \frac{dP}{dt}=-b_i P(t) + L(t) +r_i\frac{P(t)^{q_i}}{m_i^{q_i} + P(t)^{q_i}} \\
   \left(L(t),P(t)\right) \in K=[L_{min}, +\infty[ \times [0,P_{max}]
   \end{array}
   \right.
   \end{equation}

   Un individu ($i=4$) considère que le processus de relargage peut se produire différemment, et utilise une formule différente pour la pseudo-sygmoïde :

   \begin{equation}
   (S' \text{ modèle 0}) \quad \frac{dP}{dt}=-b_i P(t) + L(t) + r_i \frac{P(t)}{P(t) + m_i e^{(-\lambda_i(P(t)-m_i)) }}
   \end{equation}

   Les valeurs numériques pour les calculs sont les suivantes : $u_{min}=-u_{max}/2$ ; $u_{max}\approx 3,15\; \mu g.l^{-1}.\text{an}^{-1}$, $L_{min}\approx6,94 \; \mu g.l^{-1}$,  $P_{max}=24.76 \; \mu g.l^{-1}$, 

Le  tableau suivant rassemble les valeurs des paramètres pour chaque individu (b,r,m sont en $\mu g.l^{-1}$).

Les incertitudes sont traitées comme des variables prenant leur valeur dans des ensembles. Le problème est un problème de viabilité garantie, dans lequel le vecteur $v$ des "tyches" regroupent les paramètres incertains : $v_1$ représente le paramètre $b_i$, $v_2$ : $\alpha_i$, $v_3$ : $q_i$ et $v_4$ représente $\lambda_i$. La dynamique du modèle avec incertitudes devient :

\begin{equation}
f_v(x=(L,P),u,v)=\left( 
\begin{array}{l}
u\\
- v_1 P + L +  r \left( (1-v_2) \frac{P^{v_3}}{m^{v_3} + P^{v_3}} + v_2 \frac{P}{P + m e^{(-v_4(P-m)) }}\right)  \;
\end{array}
\right)
\end{equation}

Problèmes

Le problème de viabilité garantie qui doit être résolu est le suivant :

\begin{equation}
(P_v)\left\{
\begin{array}{lcl}
(L,P)'(t)&=&f_B((L,P)(t),u(t),v(t))\\
u(t)&\in& U=\left[ u_{min},u_{max}\right] \\
v(t)&\in & V=[2.2,2.3]\times[0,1]\times[2.2,2.3]\times[1/19,1/18]\\
(L,P)(t) &\in& K \; 
\end{array}
\right.
\end{equation}

Le noyau garanti peut être approximé avec ViabLab (version 2.2). Le code correspondant est ici : https://forgemia.inra.fr/isabelle.alvarez/emlake