# Sensitivity and optimality analysis of breathing scenarios for 1D or 0D models of gas diffusion in the lung

Céline Grandmont, Cyril Karamaoun, Sébastien Martin, Frédérique Noël

## ► To cite this version:

Céline Grandmont, Cyril Karamaoun, Sébastien Martin, Frédérique Noël. Sensitivity and optimality analysis of breathing scenarios for 1D or 0D models of gas diffusion in the lung. 2025. hal-05019689v1

HAL Id: hal-05019689

https://hal.science/hal-05019689v1

Preprint submitted on 4 Apr 2025 (v1), last revised 19 Nov 2025 (v2)

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# SENSITIVITY AND OPTIMALITY ANALYSIS OF BREATHING SCENARIOS FOR 1D OR 0D MODELS OF GAS DIFFUSION IN THE LUNG

CÉLINE GRANDMONT¹, CYRIL KARAMAOUN², SÉBASTIEN MARTIN³, AND FRÉDÉRIQUE NOËL⁴

**ABSTRACT.** In the present work we propose a new nonlinear coupled 1D model to describe lung ventilation and the transport and diffusion of both oxygen and carbon dioxide in the bronchial tree through the blood. It takes into account the *so called* Bohr-Haldane effect and is driven by the applied pleural pressure. The ability of this model to recover standard acknowledged values in healthy situations is provided. One key aspect is that, contrary to its 0D counterpart, it naturally takes into account mixing of gases along the tree and a time delay as the gases have to be transported before reaching the alveoli. We further investigate the sensitivity of both the 1D model and its 0D counterpart with respect to the breathing pattern by considering two types of pleural applied pressure: a piecewise constant one and a piecewise exponential one for various values of the breathing period, inspiratory ratio and pressure amplitude. We finally explore which cost functions the observed stereotypical breathing scenario may optimize, underlying the fact that it should be a combination of several criteria: low effort and small lung distension while maintaining carbon dioxide arterial partial pressure at a given level.

## CONTENTS

<table><tr><td>1. Motivations</td><td>1</td></tr><tr><td>2. Modelling and numerical approaches</td><td>4</td></tr><tr><td>2.1. Morphometric 1D model of the lung</td><td>4</td></tr><tr><td>2.2. 0D mechanical model of ventilation</td><td>6</td></tr><tr><td>2.3. Transport of the gases in the tree</td><td>6</td></tr><tr><td>2.4. Exchange of respiratory gases between alveoli and blood</td><td>7</td></tr><tr><td>2.5. Full 1D coupled model</td><td>8</td></tr><tr><td>2.6. A 0D model</td><td>8</td></tr><tr><td>2.7. Inputs of the models</td><td>9</td></tr><tr><td>2.8. Outputs of the models</td><td>10</td></tr><tr><td>2.9. Numerical methods</td><td>11</td></tr><tr><td>3. Numerical Experiments – Results</td><td>12</td></tr><tr><td>3.1. 1D model</td><td>12</td></tr><tr><td>3.2. Models comparison, sensitivity and optimality</td><td>13</td></tr><tr><td>4. Conclusion and perspectives</td><td>28</td></tr><tr><td>Appendix A. Numerical data</td><td>29</td></tr><tr><td>References</td><td>29</td></tr></table>

## 1. MOTIVATIONS

Mathematical modelling approaches of the respiratory system may bring major insights in the understanding of the lung functioning, in health and disease. Depending on the studied physiological phenomena and the available data, different models can be developed from macroscopic ones describing, for instance, the time evolution of global quantities [1, 3, 4, 7, 10, 30] to multiphysics and multiscale ones describing the multidimensional lung mechanics [5, 13, 14, 26, 27, 31] or airflow
[2, 8, 9] and gas diffusion [18]. The common trait between all the models is that they shall offer a strong generalizing framework or good predictive values [20].

The lung function is to transport, through breathing cycles, the air from the outside into the acini where diffusion in the blood of oxygen and carbon dioxide takes place. This diffusion process couples the dynamics of oxygen and carbon dioxide as the affinity of hemoglobin depends on pH and on carbon dioxide partial pressure in the blood. A large panel of models has been developed in the literature to describe ventilation and gases diffusion depending on the desired outputs and applications.

First, the so-called 0D models are based on a simplified representation of the lungs, with time dependency and no explicit spatial dimension. The simplest model to describe ventilation is the single balloon model where a rigid tube – the bronchial tree – is connected to an inflatable elastic balloon i.e., the alveoli [3]. More complex 0D models have been proposed, taking into account various physiological phenomena [1, 4, 7, 10, 30] and thus possibly involving a large number of parameters making their fitting and use for understanding preponderant phenomena challenging. Moreover because they do not include any spatial dimension, they are not well suited, *per se*, for studying the transport phenomena as the transport of gases in the bronchial and alveolar compartments. Note that transport induces time delay that has to be modeled and calibrated for 0D models to be accurate in particular for modelling the control of respiration [1, 7, 10, 30].

Dimensional models have been explored to tackle these issues. Full 3D models of the lungs can provide numerical simulations of air transport [2, 8, 9], lung mechanics [26] or multiscaled coupled airflow and parenchyma deformation models [5, 13, 14, 27, 31]. Yet, simulations in 3D geometries require medical images data and imply high computational costs. Moreover they are adapted to patient specific issues but may be less adaptable for global functioning understanding.

An intermediate approach is possible. One or two-dimensional models are based on key hypotheses of symmetry of the lung structure. This allows for introducing spatial dimensions in the models and consider some heterogeneities while keeping the complexity and the simulation cost low. Among these models, the 1D approach [18, 23, 25] is based on the approximation that each bronchus is a perfect cylinder presenting the same dimensions at each level of bifurcation of the bronchial tree. In essence, this strong hypothesis allows to aggregate all bronchi of the same level of subdivision i.e., generation in a single tube whose transverse surface corresponds to the sum of the transverse surface of each bronchus in this generation (see Section 2.1 for detailed modelling). This leads to a characteristic *trumpet-shaped* model of the bronchi and alveoli, with a dominant axial dimension and increased cross-sectional surface along the generations (see Fig. 1) as the number of bronchi increases exponentially at each level of subdivision. These models enable to describe the airflow as well as the transport and diffusion of oxygen [18] and carbon dioxide [23]. They can be used to simulate gas mixing along the bronchial tree, understand the impact of the breathing scenario and why one observes stereotypical scenarios for normal breathing [23].

Mechanically, the primary objective of breathing is to reduce the work involved and the mechanical elongation in the bronchi while, physiologically, it aims to sustain and regulate the amount of oxygen and carbon dioxide in the blood. Hence, the regulation of ventilation is governed, thanks to chemical sensors, by the levels of arterial pressure for oxygen ($P_{a,O_2}$) and carbon dioxide ($P_{a,CO_2}$) with $P_{a,CO_2}$ serving as the most vital regulator. To keep $P_{a,CO_2}$ near 40 mmHg, both the breathing rate and the depth of each breath are carefully controlled [17]. Different models of the control of breathing have been developed in the literature in order to recover the stereotypical breathing parameters observed in the physiology. Among those, two different approaches have appeared: a feedback control and an optimal control.

The feedback control approach often relies on complex ODE systems that take into account not only the ventilation and gas diffusion but also other parts of the body: in [10] gas concentrations in the tissues and the brain are considered, in [30] the consumption of muscles is further added and in [7] the respiratory system, the central neural system, and the cardiovascular system are modeled.
<figure><img src="image_2.png" /><figcaption>FIGURE 1. Trumpet. Light grey: alveolar part. Adapted from [16]</figcaption></figure>

In these studies, one crucial issue is to accurately model the time delay coming in particular from transport phenomena and thus to fit the various delay times appearing in the whole system.

On the other hand, the optimal control approach is based on the minimizing of the work of breathing [23, 24] or the average respiratory muscle force [21] (or even both [12]) while satisfying a constraint on the amount of oxygen exchanged with blood [23] or the alveolar ventilation rate [21, 24]. It focuses only on the respiratory system and on the exchanges made with the blood.

Our aim here is to propose a dynamical model able to describe the ventilation and gases diffusion processes that takes into account transport phenomena along the bronchial tree. One main issue is to describe the *so called* Bohr-Haldane effect [35] which induces a strong coupling between oxygen and carbon dioxide blood partial pressures. With this model at hand, our second goal is to explore which quantities control the breathing cycle: indeed, as already stated the role of the lung is to ensure the oxygen transfer from the outside air to the blood and the carbon dioxide transfer from the blood to the outside air while maintaining oxygen and carbon dioxide blood partial pressure at given levels. To achieve this, the breathing pattern has to adapt and a key question is what quantities does it optimize.

The same issues have been tackled in [23], yet with some main differences compared to what we proposed herein. Here, the breathing cycle is driven by an applied pressure (not a given flux), the proposed 1D model is based on morphometric data, includes mass-flow conservation of the gases, and the diffusion process includes the Bohr-Haldane effect. Consequently, the proposed model enables to recover standard breathing quantities without any fitting of the physiological parameters.

More precisely, we take gas transport into account by connecting a 1D geometrical model of the bronchial tree [18] to the breathing dynamics characteristic of 0D models [20]. This connected model is enriched by a recent dynamical model of alveolar gas exchange taking into account the interaction between oxygen and carbon dioxide in the blood [6]. Then, our new full model is compared to 0D model developed in [6] pointing out the transport time delay phenomena. Next both models are used
<figure><img src="image_3.png" /><figcaption>FIGURE 2. Modeling of the human lung</figcaption></figure>

to explore the impact of the breathing scenario (modifying the applied pressure), and to investigate the dependencies with respect to breathing inspiratory pressure, frequency, inspiratory ratio and understand which criteria the breathing pattern is optimizing. We investigate whether standard breathing minimizes work or the variations of the lung volume compared to the functional residual capacity $V_{FRC}$ while maintaining the carbon dioxide partial pressure in the arterial blood $P_{a,CO_2}$ as close as possible to 40 mmHg. Finally, we evaluate the sensibility of the model with respect to the breathing pattern.

## 2. MODELLING AND NUMERICAL APPROACHES

### 2.1. Morphometric 1D model of the lung.
The lung is modelled as a symmetric dichotomic bifurcating tree composed of $N$ generations (see Fig. 2). Typically there are around $N = 24$ generations in the lung of a human, numbered from $g = 0$ (trachea) to $g = 23$. The lung can then be divided into two parts, the bronchial tree which consists of the first 16 generations and the acini where the exchange with the blood is performed which consists of the last 8 generations [25, 33, 34].

Each bronchus of a generation $g$ is then represented by a cylinder with a length $l_g$ and a radius $r_g$. We can then define what we call the length of the lung $L$. It corresponds to the length of a path going from the trachea to the last generation of the tree

$$ (1) \qquad L = \sum_{g=0}^{N-1} l_g. $$

We introduce $L_k := \sum_{g=0}^k l_g$ the path length from the entrance of the trachea to the end of a bronchus of generation $k$. It means that a bronchus of generation $g$ is located along the path length by $\{x \in (L_{g-1}, L_g)\}$, with $L_{-1} = 0$. Let us now introduce lineic quantities related to the path length.

* **Bronchial volume and lineic bronchial volume.** The airway cross-section and volume of the $2^g$ cylinders of generation $g$ are $S_g := 2^g \pi r_g^2$ and $V_g := 2^g \pi r_g^2 l_g$ respectively. This allows us to define the lineic volume of the bronchial tree as

$$ (2) \qquad \tilde{V}_b(x) := \sum_{g=0}^{23} \frac{V_g}{l_g} \mathbf{1}_{(L_{g-1}, L_g)}(x), $$
which corresponds to the airway cross-section along the path: as the volumes are made of ducts, the lineic volume distribution simply reads as the piecewise constant function (on each generation), equal to the airway cross section and thus $\tilde{V}_b(x) = S(x)$ with

$$S(x) = \sum_{g=0}^{23} S_g \mathbf{1}_{(L_{g-1}, L_g)}(x).$$

The total bronchial volume $V_b$ is equal to $\int_0^L \tilde{V}_b(x) dx = \sum_{g=0}^{23} V_g$.

* **Alveolar surface area and lineic alveolar surface area.** The total alveolar surface at rest is set to $\Sigma = 120 \text{ m}^2$ and, as in [35, 29], we dispatch the corresponding surface in the acinar region only, i.e. from generation 16 to 23, as follows: the alveolar surface area is distributed to the different generations in proportion to the duct surface $2^g \times 2\pi r_g \times l_g$. This provides the distribution of alveolar area $g \mapsto \Sigma_g$ such that, on the one hand, $\Sigma_g = 0$ for $g \in \{0, \dots, 15\}$ as the generations from 0 to 15 are not alveolated and, on the other hand, $\Sigma = \sum_{g=0}^{23} \Sigma_g$. This allows us to define the *lineic surface area of alveoli* as

$$(3) \qquad \tilde{\Sigma}(x) = \sum_{g=0}^{23} \frac{\Sigma_g}{l_g} \mathbf{1}_{(L_{g-1}, L_g)}(x).$$

Noteworthily the multiplication by the function $x \mapsto \frac{\tilde{\Sigma}(x)}{\Sigma}$ allows us to convert any surface area at time $t$ into a lineic surface area, in accordance to the alveolar surface area distribution determined at rest.

* **Alveolar volume and lineic bronchial volume.** As for the alveolar surface area, we define $\Omega_g := (\Sigma_g)^{\frac{3}{2}}$ and $\Omega := \sum_{g=0}^{23} \Omega_g$. This allows us to define

$$\tilde{\Omega}(x) = \sum_{g=0}^{23} \frac{\Omega_g}{l_g} \mathbf{1}_{(L_{g-1}, L_g)}(x).$$

Noteworthily the multiplication by the function $x \mapsto \frac{\tilde{\Omega}(x)}{\Omega}$ allows us to convert any total alveolar volume at time $t$ into a lineic alveolar volume, in accordance to the alveolar volume distribution determined at rest: for a given alveolar volume at time $t$ $V_a(t)$, the lineic alveolar volume is defined as

$$\tilde{V}_a(t, x) := \frac{\tilde{\Omega}(x)}{\Omega} V_a(t).$$

By construction, $\tilde{V}_a$ and $V_a$ are related by $V_a(t) = \int_0^L \tilde{V}_a(t, x) dx$.

Let us now consider the lung volume $V(t)$ at time $t$: this volume will be the solution of a 0D mechanical model. The total volume divides into two contributions: the bronchial volume, denoted $V_b$, and defined by the morphometry and assumed to be time independent and the alveolar volume, denoted $V_a(t)$ and which corresponds to $V(t) - V_b$. Then lineic volumes are given by

$$(4) \qquad \tilde{V}_a(t, x) := \frac{\tilde{\Omega}(x)}{\Omega} (V(t) - V_b), \qquad \tilde{V}(t, x) := \tilde{V}_b(x) + \tilde{V}_a(t, x),$$

with $\tilde{V}_b$ defined by (2). By construction, we have

$$V_b = \int_0^L \tilde{V}_b(x) dx, \qquad V_a(t) = \int_0^L \tilde{V}_a(t, x) dx, \qquad V(t) = \int_0^L \tilde{V}(t, x) dx.$$
## 2.2. 0D mechanical model of ventilation.
In order to compute the evolution of the volume of the lung and the airflow through the lung, we use a well-known mechanical model that assumes that the lung can be considered as a resistive tube representing the bronchial tree linked to an elastic balloon representing the acinus [6, 20]. This model writes

$$ (5) \qquad -P(t) = E(V(t) - V_{\text{FRC}}) + R\dot{V}(t), $$

where $P(t)$ is the pleural pressure during spontaneous breathing, $E$ is the elastance of the lung, $V$ is the volume of the lung, $V_{\text{FRC}}$ is the functional residual capacity (namely, the lung volume at the end of a passive expiration) and $R$ is the resistance of the lung (that accounts for the resistance of bronchial tree and tissue viscosity).

Thanks to this model, we can compute the lung volume $V(t)$ and, as a consequence, the alveolar lineic volume $\tilde{V}_a(t, x)$ and the total lineic volume $\tilde{V}(t, x)$. Note that, as we assumed that the bronchial tree does not expand in time, the lineic bronchial volume is computed only once and does not evolve in time: $\tilde{V}_b(x) = S(x)$.

## 2.3. Transport of the gases in the tree.
To model the transport of the respiratory gases in the tree, we assume that the displacement of each gas is homogeneous in a cross section $S_g$. It allows us to model the transport of the respiratory gases with a one dimensional model that goes from $x = 0$ the entrance of the trachea to $x = L$ the end of the last generation of the acinus.

In our model, considering that air in the lung is incompressible as in [18], the time variation of the total lineic volume of the lung is equal to the flux of the fluid velocity across a section. It writes

$$ (6) \qquad \partial_t \tilde{V} + \partial_x (Su) = 0, \text{ in } (0, L), $$

where $\tilde{V}$ is the total lineic volume defined by (4) and $u$ is the velocity of the air. In the equation of the conservation of air volume, our only unknown is $u$ the velocity of the air along the tree (the volume $\tilde{V}$ defined by (4) is known from $V$ obtained by solving the mechanical equation (5)). This velocity is obtained, by integration in space of (6), using the boundary condition $S(0)u(t, 0) = \dot{V}$ during inspiration (or equivalently $S(L)u(t, L) = 0$ during expiration).

Next, we assume that the transport of the respiratory gases in the tree is mainly due to convection and diffusion. A reaction term models the exchange made with the blood in the acinus. This results in the following coupled equations for the oxygen and carbon dioxide

$$ (7) \qquad \partial_t(\tilde{V}\chi_{O_2}) + \partial_x(Su\chi_{O_2}) - D_{O_2}\partial_x(S\partial_x\chi_{O_2}) = -\frac{\tilde{\Sigma}}{\Sigma}\dot{V}_{O_2}(\chi_{O_2}, \chi_{CO_2}), \text{ in } (0, L), $$

$$ (8) \qquad \partial_t(\tilde{V}\chi_{CO_2}) + \partial_x(Su\chi_{CO_2}) - D_{CO_2}\partial_x(S\partial_x\chi_{CO_2}) = -\frac{\tilde{\Sigma}}{\Sigma}\dot{V}_{CO_2}(\chi_{CO_2}, \chi_{O_2}), \text{ in } (0, L), $$

where $\chi_{gas}$ corresponds to the mole fraction of the considered gas, $D_{gas}$ is the diffusion coefficient of the considered gas in the air and $\dot{V}_{gas}(\chi_{gas}, \cdot)$ corresponds to the amount of gas exchanged with blood: note that due to Bohr-Haldane effect, gas exchanges of $O_2$ and $CO_2$ at the alveolar-capillary membrane involve a coupling in the dynamics of $O_2$ and $CO_2$ in the bloodstream, hence in the flux of gases that goes through the membrane. Due to this, equations (7) and (8) are strongly coupled. This was not the case in [18] that considers only the dynamics of $O_2$ nor in [23] where the source term of the oxygen part was independent of the carbon dioxide. We complete this set of equations with initial conditions $\chi_{gas}^0$ and with boundary conditions depending also on breathing cycle, namely

• At inspiration,

$$ \chi_{\text{gas}}(t, 0) = \bar{\chi}_{\text{gas}}, \qquad \partial_x \chi_{\text{gas}}(t, L) = 0. $$

Here the boundary condition at $x = 0$ implies that the bronchial tree is supplied with fresh air in which the considered gas has a fraction $\bar{\chi}_{gas}$, namely 0.21 for $O_2$, 0.04 for $CO_2$ (these values could be modified according to the specific environment, for instance in the case of artificial oxygen supply). The boundary condition at $x = L$ means that the total flux
(advective+diffusive) is zero at the end of the tree (note that the advective flux is zero also as the air velocity at the end of the tree is zero).

* At expiration,

$$ \partial_x \chi_{\text{gas}}(t, 0) = 0, \quad \partial_x \chi_{\text{gas}}(t, L) = 0. $$

It means that the flux that goes out from the tree at mouth is purely advective (whereas the boundary condition at the end of the tree is the same as in the inspiration regime).

2.4. Exchange of respiratory gases between alveoli and blood. Let us first describe the transfer dynamics between the blood and the alveoli as in [6]. We assume that a volume $V_c$ of blood is brought near the acinus where the exchange can be made. It remains in the neighbourhood during a transient time $\tau_b$. The capillary volume of blood $V_c$ is then removed to be replaced by a new volume $V_c$ of blood. This results in the following dynamics for the gases concentrations, for $\theta \in [0, \tau_b]$,

$$ (9) \quad \begin{cases} V_c \frac{d}{d\theta} (C_{O_2} (P_{O_2}(\chi_{O_2}, \theta), P_{CO_2}(\chi_{CO_2}, \theta))) = D_{m,O_2} (\chi_{O_2}(x,t)P_{\text{atm}} - P_{O_2}(\chi_{O_2}, \theta)), \\ V_c \frac{d}{d\theta} (C_{CO_2} (P_{O_2}(\chi_{O_2}, \theta), P_{CO_2}(\chi_{CO_2}, \theta))) = D_{m,CO_2} (\chi_{CO_2}(x,t)P_{\text{atm}} - P_{CO_2}(\chi_{CO_2}, \theta)), \end{cases} $$

where $C_{\text{gas}}$ is the concentration of the considered gas in the blood, $P_{\text{gas}}$ is the partial pressure of the gas in the blood during a cycle $[0, \tau_b]$, $D_{m,\text{gas}}$ is the diffusive capacity of the gas through the alveolar membrane and $P_{\text{atm}}$ is the atmospheric pressure. The dynamics of the gas concentration is then driven by the difference of partial pressures in the alveoli and in the blood. Note that (9) is solved at each time $t$ and at each point $x$. Moreover note that the gas concentrations $C_{\text{gas}}$ depend on the blood partial pressures of each species and the specific forms of these functions model the Bohr-Haldane effect. Here we consider the functions used in [6]. We complete this dynamical system with the following initial conditions

$$ P_{O_2}(\chi_{O_2}, 0) = P_{v,O_2}, \quad P_{CO_2}(\chi_{CO_2}, 0) = P_{v,CO_2}, $$

where $P_{v,\text{gas}}$ corresponds to the deoxygenated blood partial pressure of the considered gas corresponding to a blood which is poor in oxygen and rich in carbon dioxide.

Finally, we define the oxygen $\dot{V}_{O_2}$ and carbon dioxide $\dot{V}_{CO_2}$ fluxes that go to and from the bloodstream. They are obtained by computing the difference between the concentration of the gas in the blood at time $\tau_b$ (oxygenated blood) and the concentration of the gas in the blood at the initial time (deoxygenated blood). They write

$$ (10) \quad \begin{cases} \dot{V}_{O_2}(\chi_{O_2}, \chi_{CO_2}) = \frac{V_c}{\tau_b} (C_{O_2}(P_{O_2}(\chi_{O_2}, \tau_b), P_{CO_2}(\chi_{CO_2}, \tau_b)) - C_{O_2}(P_{v,O_2}, P_{v,CO_2})), \\ \dot{V}_{CO_2}(\chi_{CO_2}, \chi_{O_2}) = \frac{V_c}{\tau_b} (C_{CO_2}(P_{O_2}(\chi_{O_2}, \tau_b), P_{CO_2}(\chi_{CO_2}, \tau_b)) - C_{CO_2}(P_{v,O_2}, P_{v,CO_2})). \end{cases} $$

Let us mention that the oxygenated partial pressure of the gas $P_{a,\text{gas}}$ at each time $t$ is defined as the partial pressure of the gas obtained at the final time $\tau_b$ and writes

$$ P_{a,\text{gas}}(t, x) = P_{\text{gas}}(\chi_{\text{gas}}(t, x), \tau_b). $$

In what follows, the arterial partial blood pressure shall denote the blood partial pressure in oxygenated blood.
2.5. **Full 1D coupled model.** Let us rewrite the whole model

$$
\left\{
\begin{array}{ll}
\text{Mechanical model} & \\
\quad -P(t) = E(V(t) - V_{\text{FRC}}) + R\dot{V}(t), & \\
\text{Conservation of air volume along the tree} & \\
\quad \partial_t \tilde{V} + \partial_x (Su) = 0, \quad \text{in } (0, L), & \\
\text{Gas transport in the lung} & \\
\quad \partial_t (\tilde{V}\chi_{\text{gas}}) + \partial_x (Su\chi_{\text{gas}}) - D_{\text{gas}}\partial_x (S\partial_x \chi_{\text{gas}}) = -\frac{\tilde{\nabla}}{\Sigma}\dot{V}_{\text{gas}}, \quad \text{in } (0, L), & \\
\text{Instantaneous gas fluxes through the alveolar membrane} & \\
\quad \dot{V}_{\text{O}_2}(t, x) = \frac{V_c}{\tau_b}(C_{\text{O}_2}(P_{\text{O}_2}(\chi_{\text{O}_2}(t, x), \tau_b), P_{\text{CO}_2}(\chi_{\text{CO}_2}(t, x), \tau_b)) - C_{\text{O}_2}(P_{v,\text{O}_2}, P_{v,\text{CO}_2})), & \\
\quad \dot{V}_{\text{CO}_2}(t, x) = \frac{V_c}{\tau_b}(C_{\text{CO}_2}(P_{\text{O}_2}(\chi_{\text{O}_2}(t, x), \tau_b), P_{\text{CO}_2}(\chi_{\text{CO}_2}(t, x), \tau_b)) - C_{\text{CO}_2}(P_{v,\text{O}_2}, P_{v,\text{CO}_2})), & \\
\text{Instantaneous partial pressures in arterial blood at each time } t \text{ and } x \in [0, L] & \\
\quad \text{for } \theta \in [0, \tau_b], & \\
\left\{
\begin{array}{ll}
V_c \frac{d}{d\theta} (C_{\text{O}_2}(P_{\text{O}_2}(\chi_{\text{O}_2}(t, x), \theta), P_{\text{CO}_2}(\chi_{\text{CO}_2}(t, x), \theta))) & = D_{\text{m,O}_2}(\chi_{\text{O}_2}(t, x)P_{\text{atm}} - P_{\text{O}_2}(\chi_{\text{O}_2}, \theta)), \\
V_c \frac{d}{d\theta} (C_{\text{CO}_2}(P_{\text{O}_2}(\chi_{\text{O}_2}, \theta), P_{\text{CO}_2}(\chi_{\text{CO}_2}(t, x), \theta))) & = D_{\text{m,CO}_2}(\chi_{\text{CO}_2}(t, x)P_{\text{atm}} - P_{\text{CO}_2}(\chi_{\text{CO}_2}(t, x), \theta)), \\
P_{\text{O}_2}(\chi_{\text{O}_2}, 0) = P_{v,\text{O}_2}, \quad P_{\text{CO}_2}(\chi_{\text{CO}_2}, 0) = P_{v,\text{CO}_2}.
\end{array}
\right.
\end{array}
\right.
$$

Since one of our aims is to understand the time delay induced by the transport of air along the respiratory tree we shall compare numerically the considered model to its 0D analogue introduced and studied in [6], which we recall below for the sake of completeness.

2.6. **A 0D model.** The 0D analogue of our 1D model can be described as follows: we consider the same 0D mechanical model that provides the evolution of the lung volume $t \mapsto V(t)$. We assume that this space is occupied by air, including a specific gas (e.g. $O_2$ and $CO_2$) with uniform molar fraction $t \mapsto \chi_{\text{gas}}(t)$. The available volume for gas exchange is the volume $V$ to which we subtract the dead volume, denoted $V_D$. This “dead space” volume does not take part in the exchange process and can be split in a physiological part and an anatomic part. The volume $V_D$ is defined as the sum of a reference dead space and $\max(0, V_{\text{endexpi}} - V_{\text{FRC}})$, where $V_{\text{endexpi}}$ denotes the end expiration volume at each breathing cycle. Since we are considering herein non-pathological situations, the dead volume shall be close to the reference anatomical dead volume. Taking into account this dead space volume, it means that the volume of the specific gas in the lung volume is $\chi_{\text{gas}}(V - V_D)$. During inspiration, the rate of gas which enters the balloon is $\bar{\chi}_{\text{gas}}\dot{V} > 0$, which corresponds to a mixture of fresh air and lung air. Consequently $\bar{\chi}_{\text{gas}}$ is defined by $\bar{\chi}_{\text{gas}} = (1 - \phi(t))\chi_{\text{gas}}^{\text{atm}} + \phi(t)\chi_{\text{gas}}$, with $\phi(t) = V_D/(V(t) - V_{\text{endexpi}} + V_D)$, with $\chi_{\text{O}_2}^{\text{atm}} = 0.21$ and $\chi_{\text{CO}_2}^{\text{atm}} = 0.04$. During expiration, the rate of gas writes $\chi_{\text{gas}}\dot{V} < 0$ as air is rejected (including gas at current fraction $\chi_{\text{gas}}$). Now denoting by $\dot{V}_{\text{gas}}$ the flux of gas which diffuses through the alveolar-blood membrane, the dynamics of the gas fraction in the lung is governed by

$$
\frac{d}{dt} (\chi_{\text{gas}}(V - V_D)) = H(\dot{V}) \bar{\chi}_{\text{gas}} \dot{V} + (1 - H(\dot{V})) \chi_{\text{gas}} \dot{V} - \dot{V}_{\text{gas}},
$$
where $H(\cdot)$ is the Heaviside graph (characteristic function of $\mathbb{R}^+$). One obtains

$$ \dot{\chi}_{\text{gas}} = \frac{1}{V - V_D} \left( \dot{V}(1 - \phi)(\chi_{\text{gas}}^{\text{atm}} - \chi_{\text{gas}})H(\dot{V}) - \dot{V}_{\text{gas}} \right). $$

The flux of oxygen $\dot{V}_{O_2}$ and carbon dioxide $\dot{V}_{CO_2}$ are defined as in the 1D model by (10) and the full reduced 0D model thus writes

$$
\left\{
\begin{array}{l}
\text{Mechanical model} \\
- P(t) = E(V(t) - V_{\text{FRC}}) + R\dot{V}(t), \\
\\
\text{Gas exchange} \\
\dot{\chi}_{\text{gas}}(t) = \frac{1}{V(t)-V_D} \left( \dot{V}(t)(1-\phi(t))(\chi_{\text{gas}}^{\text{atm}} - \chi_{\text{gas}}(t))H(\dot{V}(t)) - \dot{V}_{\text{gas}}(t) \right), \\
\\
\text{Instantaneous gas fluxes through the alveolar membrane,} \\
\dot{V}_{O_2}(t) = \dot{V}_{O_2}(\chi_{O_2}(t), \chi_{CO_2}(t)) = \frac{V_c}{\tau_b} (C_{O_2}(P_{O_2}(\chi_{O_2}(t), \tau_b), P_{CO_2}(\chi_{CO_2}(t), \tau_b)) - C_{O_2}(P_{v,O_2}, P_{v,CO_2})), \\
\dot{V}_{CO_2}(t) = \dot{V}_{CO_2}(\chi_{CO_2}(t), \chi_{O_2}(t)) = \frac{V_c}{\tau_b} (C_{CO_2}(P_{CO_2}(\chi_{CO_2}(t), \tau_b), P_{O_2}(\chi_{O_2}(t), \tau_b)) - C_{CO_2}(P_{v,CO_2}, P_{v,O_2})), \\
\\
\text{Instantaneous partial pressures in arterial blood at time } t, \theta \in [0, \tau_b], \\
\left\{
\begin{array}{l}
V_c \frac{d}{d\theta} (C_{O_2}(P_{O_2}(\chi_{O_2}(t), \theta), P_{CO_2}(\chi_{CO_2}(t), \theta))) = D_{m,O_2}(\chi_{O_2}(t)P_{\text{atm}} - P_{O_2}(\chi_{O_2}(t), \theta)), \\
\\
V_c \frac{d}{d\theta} (C_{CO_2}(P_{O_2}(\chi_{O_2}(t), \theta), P_{CO_2}(\chi_{CO_2}(t), \theta))) = D_{m,CO_2}(\chi_{CO_2}(t)P_{\text{atm}} - P_{CO_2}(\chi_{CO_2}(t), \theta)), \\
\\
P_{O_2}(\chi_{O_2}, 0) = P_{v,O_2}, \qquad P_{CO_2}(\chi_{CO_2}, 0) = P_{v,CO_2}.
\end{array}
\right.
\end{array}
\right.
$$

Note that in this 0D model, the mixing of fresh air and lung air and the dead space have to be taken into account to recover standard physiological behavior (see Remark 8 in [6]). In the 1D model this mixing and the dead space are taken into account by the “geometry”. Indeed, the mixing takes place during air transport along the bronchial tree in which there is no gas exchange, and this bronchial region represents the dead space region.

2.7. **Inputs of the models.** Let us note that our model is totally controlled by the pleural pressure $t \mapsto P(t)$. This function is the driving force of the ventilation. It is regulated to achieve good gases exchanges. Our aim is to explore numerically the ability of our 1D model to reproduce normal breathing scenarios. We shall consider periodic breathing scenarios with two different shapes: in the first scenario, the applied pressure on the lung is piecewise constant and writes

$$ (11) \qquad P(t) = \begin{cases} P_{\text{insp}}, & t \in [0, \text{IE} \times T[, \\ P_{\text{exp}}, & t \in [\text{IE} \times T, T[, \end{cases} $$

where $T$ is the respiratory cycle period and $\text{IE}$ is the inspiratory ratio. Typical values for normal breathing are given by

$$ T = 5 \text{ s}, \quad \text{IE} = 0.4, \quad P_{\text{insp}} = -2 \text{ cmH}_2\text{O}, \quad P_{\text{exp}} = 0 \text{ cmH}_2\text{O}. $$

In the second scenario, we propose another pressure pattern: we consider increasing (resp. decreasing) exponential shape at inspiration (resp. at expiration). This scenario corresponds to a piecewise constant muscle activity. This signal is transformed into an exponential applied pressure through a resistance-compliance (RC) model also called Voigt model of viscoelastic behavior [22]. It writes

$$ (12) \qquad P_m(t) = \tau \dot{P}(t) + P(t), $$
<figure><img src="image_4.png" /><figcaption>FIGURE 3. Pleural pressure for two scenarios: Timeslot profile defined by (11) with $T = 5$ s, $P_{\text{insp}} = -2.0$ cmH$_2$O, $P_{\text{exp}} = 0$ cmH$_2$O, $\text{IE} = 0.4$, vs. Exponential profile defined by (12) with $T = 5$ s, $P_{\text{insp}} = -2.4$ cmH$_2$O, $P_{\text{exp}} = 0$ cmH$_2$O, $\text{IE} = 0.4$. Both scenarios correspond to the same tidal volume $V_T = 0.49$ L.</figcaption></figure>

where $P_m(t)$ is the muscle activity pressure, $\tau = RC$ is the characteristic time and $P(t)$ is the pressure applied on the lung. The pressure applied on the respiratory muscles is defined as follows

$$ (13) \qquad P_m(t) = \begin{cases} P_{\text{insp}}, & t \in [0, \text{IE} \times T[, \\ P_{\text{exp}}, & t \in [\text{IE} \times T, T[. \end{cases} $$

The reference exponential pattern has been chosen with characteristic time $\tau = 0.7$ s and typical values for normal breathing are given by

$$ T = 5 \text{ s, } \text{IE} = 0.4, \quad P_{\text{insp}} = -2.4 \text{ cmH}_2\text{O, } P_{\text{exp}} = 0 \text{ cmH}_2\text{O.} $$

This breathing scenario leads to a tidal volume of 0.49 L, as obtained with the piecewise constant pattern. Standard breathing patterns are illustrated in Fig. 3.

## 2.8. Outputs of the models.
For comparison purposes and sensitivity analysis, we choose to consider the following outputs of our models.

* The mean amount of oxygen and carbon dioxide exchanged with the blood during a respiratory cycle. It is computed thanks to the following expression

$$ (14) \qquad \bar{V}_{\text{gas}} = \frac{1}{T} \int_{0}^{T} \left( \int_{0}^{L} \frac{\bar{\Sigma}(x)}{\Sigma} \dot{V}_{\text{gas}}(t, x) dx \right) dt. $$

* The mean arterial partial pressure of oxygen and carbon dioxide during a respiratory cycle. It writes

$$ (15) \qquad \bar{P}_{a,\text{gas}} = \frac{1}{T} \int_{0}^{T} \left( \int_{0}^{L} \frac{\bar{\Sigma}(x)}{\Sigma} P_{a,\text{gas}}(t, x) dx \right) dt. $$

* The tidal volume $V_T$, which is the amount of air inhaled during a respiratory cycle.

* The power of external forces $W$, which is defined as the opposite of the product between the applied pleural pressure and the airflow going through the lung [20] and writes

$$ (16) \qquad W = \frac{1}{T} \int_{0}^{T} -P(t) \dot{V}(t) dt. $$
Thanks to equation (5), and to the fact that $V(T) = V(0)$ we can rewrite it as

$$W = \frac{1}{T} \int_{0}^{T} R \dot{V}^{2}(t) dt \geq 0.$$

* The variation of the maximal volume with respect to a given $V_{\text{FRC}}$ which corresponds to the amount of distension of the lung during a respiratory cycle and writes

$$\mathcal{D}_V = (V_{\text{max}} - V_{\text{FRC}})^2.$$

## 2.9. Numerical methods.
We use the same numerical method developed in [18]. First, we discretize the mechanical 0D model (5) with an implicit Euler scheme. It gives us the volume $V$ and the flux $\dot{V}$ which corresponds also to the derivative of the volume at each time.

Next, our PDEs are approximated considering a finite volume scheme. We discretize the time interval $[0, T]$ with $N$ points $(t^n)_{1 \le i \le N}$ and the space domain $[0, L]$ with $M$ points $(x^i)_{1 \le i \le M}$. Let us define at each nodes $(x^i, t^n)_{n=1,...,N}$, the derivative of the lineic volume $\partial_t \tilde{V}^{i,n}$ which corresponds to the airflow $\dot{V}^{i,n}$ per unit of length. The conservation of air volume can then be discretized with an upwind scheme which reads

$$(17) \qquad \left( \frac{(Su)^{i,n} - (Su)^{i-1,n}}{\Delta x_g} \right) \mathbb{1}_{\dot{V}^n > 0} + \left( \frac{(Su)^{i+1,n} - (Su)^{i,n}}{\Delta x_g} \right) \mathbb{1}_{\dot{V}^n \le 0} = -\partial_t \tilde{V}^{i,n},$$

where $\Delta x_g$ is an adaptive space step depending of the $g$th generation. Since the length of the bronchi decrease generation per generation this adaptive space step is needed to have enough discretization points in the last generations of the respiratory tree without refining in the whole domain. We choose $\Delta x_g$ such that we have the same number of points per generation. We complete our equation with the following boundary conditions: $(Su)^{1,n} = \dot{V}^n$ during inspiration ($\dot{V}^n > 0$) and $(Su)^{M,n} = 0$ at expiration ($\dot{V}^n < 0$).

Let us now denote $\chi_{\text{gas}}^{i,n}$ the approximations of the mole fractions at each point $(x^i, t^n)$. The transport equation of the respiratory gases is then discretized as

$$(18) \qquad \frac{(\tilde{V}\chi_{\text{gas}})^{i,n+1} - (\tilde{V}\chi_{\text{gas}})^{i,n}}{\Delta t} + \frac{(\Phi_{\text{gas}}^{i+1/2,n+1} - \Phi_{\text{gas}}^{i-1/2,n+1})}{\Delta x_g} = -\frac{\tilde{\Sigma}(x^i)}{\Sigma} \dot{V}_{\text{gas}}^{i,n},$$

where $\Phi_{\text{gas}}^{i+1/2,n+1}$ and $\Phi_{\text{gas}}^{i-1/2,n+1}$ are approximations of the total flux which can be expressed as

$$\left\{ \begin{array}{l} \Phi_{\text{gas}}^{i+1/2,n+1} = (Su)^{i,n+1} \chi_{\text{gas}}^{i,n+1} \mathbb{1}_{\dot{V}^n > 0} + (Su)^{i+1,n+1} \chi_{\text{gas}}^{i+1,n+1} \mathbb{1}_{\dot{V}^n \le 0} - D_{\text{gas}} S^{i+1/2,n+1} \frac{\chi_{\text{gas}}^{i+1,n+1} - \chi_{\text{gas}}^{i,n+1}}{\Delta x_g}, \\ \\ \Phi_{\text{gas}}^{i-1/2,n+1} = (Su)^{i-1,n+1} \chi_{\text{gas}}^{i-1,n+1} \mathbb{1}_{\dot{V}^n > 0} + (Su)^{i,n+1} \chi_{\text{gas}}^{i,n+1} \mathbb{1}_{\dot{V}^n \le 0} - D_{\text{gas}} S^{i-1/2,n+1} \frac{\chi_{\text{gas}}^{i,n+1} - \chi_{\text{gas}}^{i-1,n+1}}{\Delta x_g}. \end{array} \right.$$

We complete the equations with boundary conditions that depend on the breathing cycle. During inspiration, we assume that the molar fraction at the entrance of the trachea is the same as the one in the air, $\chi_{\text{gas}}^{1,n} = \chi_{\text{gas}}^{\text{atm}}$, and during expiration, we assume that there is a free flux at the entrance of the trachea, $\chi_{\text{gas}}^{1,n} = \chi_{\text{gas}}^{2,n}$. For the whole breathing cycle, we impose at the end of the tree that $\Phi_{\text{gas}}^{M+1/2,n} = 0$.

Note that the flux $\dot{V}_{\text{gas}}^{i,n}$ is treated in an explicit way in (18) and computed from the previous time step. More precisely, the approximation of the oxygen and carbon dioxide fluxes $\dot{V}_{\text{gas}}$ are computed by resolving the differential equations (9) using a 4th order Runge-Kutta scheme. Thanks to this resolution, we have access to the approximation of the arterial partial pressure of the each gas at each node $P_{a,\text{gas}}^{n,i} = P_{\text{gas}}(\chi_{\text{gas}}^{i,n}, \tau_b)$. We can then compute the mean value of the partial pressure of
the gases in the arterial blood during a respiratory cycle with

$$ (19) \qquad \bar{P}_{a,gas} = \frac{1}{N} \sum_{n=1}^{N} \sum_{i=1}^{M} \frac{\tilde{\Sigma}(x^i)}{\Sigma} \Delta x_g P_{a,gas}^{n,i} $$

where $\tilde{\Sigma}$ corresponds to the lineic alveolar surface defined previously in equation (3) and $\Sigma$ corresponds to the total alveolar surface.

Likewise, we can compute the mean amount of oxygen and carbon dioxide exchanged with the blood during a respiratory cycle by

$$ (20) \qquad \bar{V}_{gas} = \frac{1}{N} \sum_{n=1}^{N} \sum_{i=1}^{M} \frac{\tilde{\Sigma}(x^i)}{\Sigma} \Delta x_g \dot{V}_{gas}^{i,n} $$

### 3. NUMERICAL EXPERIMENTS – RESULTS

In this section, we numerically investigate the behavior of the two models. First we consider the 1D model for standard breathing scenario in the normal breathing case and then we compare the 0D model and 1D model by performing simulations with different applied pressure making $T$, $IE$, $P_{insp}$ vary while $P_{exp}$ remains constant and equal to 0 cmH$_2$O since $P_{insp}+P_{exp}$ is an invariant of the linear mechanical model (5) with the piecewise constant pressure. The objectives are multiple: to assess the ability of the 1D model to recover standard physiological behavior for normal breathing, to investigate the importance of the time delay induced by the transport along the tree, to explore the sensitivity of the models with respect to the applied pressure, to analyse which are the criteria the breathing scenario may optimize and understand the ability of each model to accurately reproduce physiological behavior. Note that in order to reach convergence towards the periodic solution for each simulation, we perform several cycles (around 100).

In all the simulations, the same physiological parameters are used and correspond to standard values that can be found in the literature for normal breathing and an healthy situation (see section A and Table 12). Moreover the geometry parameters are the physiological ones provided in [34].

We shall emphasize that there are three main differences between the chosen parameters and the one considered in [23]. The geometry is different since in [23] a homothetic geometry is chosen (each generation length and radius is obtained by a geometric law from the morphometry of the trachea) which gives a different length and radius distribution along the tree and thus different surface exchange. The second main difference in the parameters is the gas diffusion parameter through the alveolo-capillary membrane. In [23], its value was chosen much higher than in our simulations, leading to, possibly, overestimate the gas exchanges. As a consequence the given deoxygenated blood partial pressures of gases $P_{v,gas}$ have to be fitted in order to recover good orders of magnitude of the output data. These pressures were seen as the pressures of a mixture of deoxygenated blood and oxygenated one, surely to counterbalance the high diffusion coefficient between alveolar and blood partial pressures that tends to overestimate the gas fluxes.

Note that with the considered values of the parameters, we recover the standard values of the outputs without any fitting.

3.1. **1D model.** Let us first present numerical simulations of the 1D model for two different external forcing terms corresponding to Fig. 3. The associated periodic volume and flux are given by Fig. 4. In particular, the airflow instantaneously reaches its maximum at the beginning of inhalation and exhalation for the timeslot pattern while the dynamics of the air flow is smoother for the exponential pressure pattern. Yet, both scenarios correspond to a similar maximal inhaled volume and a similar tidal volume $V_T$. Note that, in the exponential scenario, the volume does not return to $V_{FRC}$ and the airflow has a smoother dynamics. Despite these differences, the average quantities by generation over a time period are strongly similar (see Fig. 5) as well as the average fluxes with respect to
<figure><img src="image_5.png" /><figcaption>FIGURE 4. Airflow $t \mapsto \dot{V}(t)$ and lung volume $t \mapsto V(t)$ for two scenarios of pleural pressure: Timeslot profile defined by (11) with $T = 5$ s, $P_{\text{insp}} = -2.0$ cmH$_2$O, $P_{\text{exp}} = 0$ cmH$_2$O, $\text{IE} = 0.4$ vs. Exponential profile defined by (12) with $T = 5$ s, $P_{\text{insp}} = -2.4$ cmH$_2$O, $P_{\text{exp}} = 0$ cmH$_2$O, $\text{IE} = 0.4$. Each scenario corresponds to the same tidal volume $V_T = 0.49$ L.</figcaption></figure>

time. We also observed that the fluxes of gas are important at the end of the tree. As in [18, 32] and to the contrary of [23, 28], no screening effect is observed in our model.

In Fig. 6, the magnitudes and time evolution of the average oxygen and carbon dioxide fluxes are really similar in both scenarios and a time delay can be observed in both cases: smaller in the timeslot case and surely due to the transport along the tree. A larger time delay coming from the exponential pattern of the pressure appears in the other case. The same time delay is observed in Fig. 7 for the time evolution of the alveolar partial pressure of carbon dioxide in different regions of the bronchial tree. Observe that, at the end of the acinus, the partial pressure is almost constant over time and equal to 40 mmHg, which is the value of the oxygenated blood partial pressure.

The main gas exchanges take place in the last generations of the tree and in particular 73% of the exchanges take place in the last two generations, enlightening that no screening effect is present in our model. One explanation is that in [23] the coefficient of diffusion was chosen much higher than the acknowledged value chosen herein. Thus in [23] the gases exchanges take place upper in the respiratory tree and require a smaller surface of exchange. Moreover we see that, despite smaller airflow in the exponential scenario, the exchanges are sufficient to reproduce in a really accurate way the well-documented physiological quantities as shown in Table 5.

The observable quantities that are commonly referred to in the literature are the tidal volume $V_T$, the mean oxygen and carbon dioxide blood partial pressures $\bar{P}_{a,gas}$, the mean oxygen and carbon dioxide fluxes $\bar{V}_{gas}$ through the alveolar-capillary membrane, and the respiratory quotient

$$RQ = \frac{-\bar{V}_{CO_2}}{\bar{V}_{O_2}}.$$

In our simulations, the standard values:

$$\bar{V}_{O_2} = 250 \text{ mL} \cdot \text{min}^{-1}, \bar{V}_{CO_2} = -200 \text{ mL} \cdot \text{min}^{-1}, RQ = 0.8,$$

$$\bar{P}_{a,O_2} = 100 \text{ mmHg}, \bar{P}_{a,CO_2} = 40 \text{ mmHg},$$

are recovered (as seen in Table 1) without any fitting of the input parameters.

3.2. Models comparison, sensitivity and optimality. The 1D model has been developed in order to take into account the heterogeneity of transport and gas exchange along the bronchial tree. The 1D structure has an impact on the gas exchange, due to the time that is required for fresh air
<figure><img src="image_6.png" /></figure>

FIGURE 5. Time average by generation of $P_{alv,O_2}$, $P_{alv,CO_2}$, $\dot{V}_{O_2}$, $\dot{V}_{CO_2}$ for both scenarios: piecewise constant applied pressure defined by (11) with $T = 5$ s, $P_{insp} = -2$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$ and piecewise exponential defined by (12) with $T = 5$ s, $P_{insp} = -2.4$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$.

<figure><img src="image_7.png" /></figure>

FIGURE 6. Time evolution of $\dot{V}_{O_2}$, $\dot{V}_{CO_2}$, for both scenarios: piecewise constant applied pressure, defined by (11) with $T = 5$ s, $P_{insp} = -2$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$ and piecewise exponential, defined by (12) with $T = 5$ s, $P_{insp} = -2.4$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$.
<figure><img src="image_8.png" /><figcaption>FIGURE 7. Time evolution of the alveolar pressure $P_{alv,CO_2}$ for both scenarios: piecewise constant applied pressure, defined by (11) with $T = 5$ s, $P_{insp} = -2$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$ (dot line) and piecewise exponential, defined by (12) with $T = 5$ s, $P_{insp} = -2.4$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$ (plain line).</figcaption></figure>

<table><caption>TABLE 1. Mean outputs over a breathing cycle provided by the resolution of the 1D problem for both scenarios: piecewise constant applied pressure, defined by (11) with $T = 5$ s, $P_{insp} = -2$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$ and piecewise exponential, defined by (12) with $T = 5$ s, $P_{insp} = -2.4$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$.</caption><thead><tr><th>1D model</th><th>Timeslot</th><th>Exponential</th></tr></thead><tbody><tr><td>V̄<sub>O₂</sub> (mL/min)</td><td>258</td><td>256</td></tr><tr><td>V̄<sub>CO₂</sub> (mL/min)</td><td>-200</td><td>-192</td></tr><tr><td>RQ</td><td>0.78</td><td>0.75</td></tr><tr><td>P̄<sub>a,O₂</sub> (mmHg)</td><td>100.1</td><td>98.3</td></tr><tr><td>P̄<sub>a,CO₂</sub> (mmHg)</td><td>40.3</td><td>40.5</td></tr><tr><td>V<sub>T</sub> (L)</td><td>0.49</td><td>0.49</td></tr></tbody></table>

<table><caption>TABLE 2. Mean outputs over a breathing cycle provided by the resolution of the 0D problem for both scenarios: piecewise constant applied pressure, defined by (11) with $T = 5$ s, $P_{insp} = -2$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$ and piecewise exponential, defined by (12) with $T = 5$ s, $P_{insp} = -2.4$ cmH$_2$O, $P_{exp} = 0$ cmH$_2$O, $IE = 0.4$.</caption><thead><tr><th>0D model</th><th>Timeslot</th><th>Exponential</th></tr></thead><tbody><tr><td>V̄<sub>O₂</sub> (mL/min)</td><td>250</td><td>247</td></tr><tr><td>V̄<sub>CO₂</sub> (mL/min)</td><td>-179</td><td>-172</td></tr><tr><td>RQ</td><td>0.72</td><td>0.70</td></tr><tr><td>P̄<sub>a,O₂</sub> (mmHg)</td><td>92.9</td><td>90.9</td></tr><tr><td>P̄<sub>a,CO₂</sub> (mmHg)</td><td>40.9</td><td>41.1</td></tr><tr><td>V<sub>T</sub> (L)</td><td>0.49</td><td>0.49</td></tr></tbody></table>

to reach the acinus through the physiological dead space. However, it is expected that the 0D and 1D models provide similar results at least in the healthy case, under normal breathing conditions.

Indeed, for both considered scenarios, the two models behave in a really similar manner as stated by the results summarized in Table 1 for the 1D model and Table 2 for the 0D model.
<table><caption>TABLE 3. 1D vs 0D models in hyperventilation: piecewise constant applied pressure defined by (11) with $T = 2$ s, $P_{\text{insp}} = -4.9$ cmH₂O, $P_{\text{exp}} = 0$ cmH₂O, $\text{IE} = 0.1$.</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td>V&#773;<sub>O₂</sub> (mL/min)</td><td>272</td><td>262</td></tr><tr><td>V&#773;<sub>CO₂</sub> (mL/min)</td><td>-286</td><td>-201</td></tr><tr><td>RQ</td><td>1.05</td><td>0.77</td></tr><tr><td>P&#773;<sub>a,O₂</sub> (mmHg)</td><td>114.2</td><td>103.7</td></tr><tr><td>P&#773;<sub>a,CO₂</sub> (mmHg)</td><td>37.3</td><td>40.3</td></tr><tr><td>V&#773;<sub>T</sub> (L)</td><td>0.40</td><td>0.40</td></tr></tbody></table>

As observed in Table 1, and Table 2 the main observable quantities of interest are only slightly affected by the choice of the pressure model in the normal regime with healthy parameters.

Note that, contrary to the 1D model, the 0D model does not take into account the transit time through the bronchial tree which shall induce a time delay due to the transport and the mixing of gases in the dead space region. However, note that for both chosen normal breathing scenarios, the 1D model predicts higher gases exchanges. One reason may be that the time delay is beneficial to maintain higher exchanges in average over a time period and could be due to the fact that, at expiration, the gases stay longer in the acinus region enabling higher exchanges. Indeed, as already stated the advection-diffusion 1D model naturally takes into account time delay and mixing of gases.

It may also have an impact on the gas distribution in particular in regimes that involve small or large characteristic times (e.g. hyperventilation) or pathological parameters. It may also influence the sensitivity of the models with respect to the breathing scenarios or the optimality of the reference breathing scenario.

### 3.2.1. Dependency of the 1D and 0D models with respect to the inspiration fraction.
To explore the effect of the transit time in the bronchial tree, let us consider breathing scenarios with different inspiratory ratios.

In order to investigate the differences between both models, let us first consider three hyperventilation piecewise constant scenarios. A hyperventilation scenario is characterized by a breathing pattern such that the minute volume defined by $60V_T/T$ is greater than in the reference case (which corresponds to a value of 6 in the considered reference case). We consider only a piecewise constant breathing pattern in order to compare the models with no time delay in the applied pressure.

In Table 3, corresponding to a small period ($T = 2$ s), a large inspiratory pressure ($P_{\text{insp}} = -4.9$ cmH₂O), and a small inspiratory ratio ($\text{IE} = 0.1$) the output values for the 1D model are similar yet slightly higher as for the reference case (see Table 1), whereas for the 0D model the exchanges are improved compared to the reference case (see Table 2). Here, as the inspiration ratio is small compared to the transit time in the bronchial tree, the exchanges can not be improved for the 1D model as they are with the 0D model where the time delay is not taken into account.

Table 4 shows that for an other hyperventilation scenario in which the period is small but the inspiratory ratio is large, namely $P_{\text{insp}} = -2.0$ cmH₂O, $P_{\text{exp}} = 0$ cmH₂O, $\text{IE} = 0.7$, the 1D model predicts better exchanges than its 0D counterpart. One reason may be that the equilibrium between the alveolar pressures and the blood pressures is reached more rapidly in the 0D problem than in the 1D one.

Finally, in Table 5, both models behave similarly for $\text{IE} = 0.5$. This tends to show that to optimize the exchange for the 0D model a small inspiration fraction is recommended, which is not the case for the 1D model, for which bigger inspiration fraction is recommended. Note that in every chosen hyperventilation scenario the exchanges are higher than in the reference case given in Table 1 and Table 2. To investigate deeper the dependency with respect to the inspiration fraction, we focus on the blood partial pressure of carbon dioxide for various breathing scenarios with $\text{IE}$ varying
<table><caption>TABLE 4. 1D vs 0D models in hyperventilation: piecewise constant applied pressure defined by (11) with $T = 2$ s, $P_{\text{insp}} = -2.0$ cmH₂O, $P_{\text{exp}} = 0$ cmH₂O, $IE = 0.7$.</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td>V&#775;<sub>O₂</sub> (mL/min)</td><td>264</td><td>274</td></tr><tr><td>V&#775;<sub>CO₂</sub> (mL/min)</td><td>-231</td><td>-295</td></tr><tr><td>RQ</td><td>0.87</td><td>1.08</td></tr><tr><td>P&#773;<sub>a,O₂</sub> (mmHg)</td><td>104.6</td><td>118.1</td></tr><tr><td>P&#773;<sub>a,CO₂</sub> (mmHg)</td><td>39.1</td><td>37.0</td></tr><tr><td>V<sub>T</sub> (L)</td><td>0.40</td><td>0.40</td></tr></tbody></table>

<table><caption>TABLE 5. 1D vs 0D models in hyperventilation: piecewise constant applied pressure defined by (11) with $T = 2$ s, $P_{\text{insp}} = -2.0$ cmH₂O, $P_{\text{exp}} = 0$ cmH₂O, $IE = 0.5$.</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td>V&#775;<sub>O₂</sub> (mL/min)</td><td>271</td><td>271</td></tr><tr><td>V&#775;<sub>CO₂</sub> (mL/min)</td><td>-272</td><td>-267</td></tr><tr><td>RQ</td><td>1.00</td><td>0.99</td></tr><tr><td>P&#773;<sub>a,O₂</sub> (mmHg)</td><td>111.9</td><td>114.5</td></tr><tr><td>P&#773;<sub>a,CO₂</sub> (mmHg)</td><td>37.7</td><td>38</td></tr><tr><td>V<sub>T</sub> (L)</td><td>0.4</td><td>0.4</td></tr></tbody></table>

<figure><img src="image_9.png" /><figcaption>FIGURE 8. $\bar{P}_{a,CO₂}$ vs IE for the 1D and 0D models with a piecewise constant applied pressure defined by (11) with $T = 5$ s, $P_{\text{insp}} = -2.0$ cmH₂O, $P_{\text{exp}} = 0$ cmH₂O.</figcaption></figure>

and the period and amplitude of the applied pressure fixed, both for the piecewise constant and piecewise exponential patterns. We consider a normal breathing scenario and study the differences of predicted exchanges for the 0D model and the 1D model with respect to the inspiration fraction.

In Fig. 8, below a value $IE_0 \approx 0.25$, the blood partial pressure of carbon dioxide predicted by the 0D model is lower than the one predicted by the 1D model, showing that, on average, the exchange are higher for the 1D model. In particular, if the targeted blood partial pressure of carbon dioxide is between 40 and 41 mmHg, then the range of inspiratory fractions is larger for the 1D model than for the 0D one. Moreover, no scenario gives partial pressure below 40mmHg for the 0D model, to the contrary of the 1D model for which the exchanges are significantly higher for $IE > IE_0$.
<figure><img src="image_10.png" /><figcaption>FIGURE 9. $\bar{P}_{a,CO_2}$ vs IE for the 1D and 0D models and a piecewise exponential applied pressure defined by (12) with $T = 5$ s, $P_{\text{insp}} = -2.5$ cmH$_2$O, $P_{\text{exp}} = 0$ cmH$_2$O.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIGURE 10. Time evolution of $\dot{V}_{CO_2}$ over a breathing period for three different values of the inspiratory ratio: $IE = 0.1, 0.25, 0.5$. Comparison of the 1D model (plain line) and the 0D model (dash - dot line).</figcaption></figure>

Considering a piecewise exponential pattern, Fig. 9 gives the same behavior, except that the minimal blood partial pressure of carbon dioxide is reached for lower values of IE than in the piecewise constant case. This is probably due to the maximum airflow which is much larger in the piecewise constant case, which allows the carbon dioxide to be washed out more easily.

Thus, for bigger inspiratory fraction, the 1D model predicts better exchanges, whereas for small ones it is the 0D model. This may be explained by the transit time in the bronchial tree: for small inspiratory fraction, fresh air fails to reach the acini and for large inspiratory fraction the 0D model reaches more rapidly equilibrium, minimizing the exchanges.

In order to understand whether the difference in the behavior of both models with respect to the inspiratory fraction appears only in average, we consider the time evolution of the carbon dioxide flux over a respiration cycle. We see in Fig. 10 that, when the inspiration fraction is below $IE_0$, then the fluxes are greater for the 0D model on the whole period, whereas it is the contrary for $IE > IE_0$. For $IE = IE_0$ then the flux $\dot{V}_{CO_2}$ is greater at the beginning of inspiration for the 0D model then the exchanges are better for in the 1D case and both models behave similarly at expiration.

### 3.2.2. Sensitivity analysis with respect to the breathing parameters.
In this subsection, we investigate the sensitivity of the models with respect to the breathing scenario by making the time period $T$, the
inspiratory pressure $P_{\text{insp}}$ and the inspiratory fraction $\text{IE}$ vary for a healthy case. First, we consider a piecewise constant pattern with $T \in [2.8, 6]$ s, $P_{\text{insp}} \in [-3, -1]$ cmH$_2$O and $\text{IE} \in [0.2, 0.8]$, then a piecewise exponential pattern with $T \in [2.8, 6]$ s, $P_{\text{insp}} \in [-3.5, -1.5]$ cmH$_2$O and $\text{IE} \in [0.2, 0.8]$.

The output variables are the average blood partial pressures of each gas over a time period, namely $\bar{P}_{a,\text{O}_2}$ and $\bar{P}_{a,\text{CO}_2}$ given by (15).

In the series of Figures 11-14, 15-18, 19-22 and 23-26 each dot of coordinates $(\bar{P}_{a,\text{CO}_2}, \bar{P}_{a,\text{O}_2})$ corresponds to a simulation for the 0D model and 1D model respectively, with piecewise constant or piecewise exponential patterns. These dots are colored with respect to the input parameters, namely $T$, $|P_{\text{insp}}|$, $\text{IE}$ as well as the work of the applied pressure $W$ (defined by (16)). Moreover, for each simulation, we zoom around $\bar{P}_{a,\text{CO}_2} \approx 40$; the dots are then represented with correction by the quadratic regression curve.

First, we note that for a given $\bar{P}_{a,\text{CO}_2}$, the blood partial pressure of oxygen $\bar{P}_{a,\text{O}_2}$ is fixed in the 0D case. For the 1D case, there is more dispersion but $\bar{P}_{a,\text{O}_2}$ has small variations around a given value. Thus, both models provide an operating point: for instance for $\bar{P}_{a,\text{CO}_2} = 40$ mmHg, $\bar{P}_{a,\text{O}_2} \approx 100$ mmHg, which are the standard values for healthy patients. Yet, there is no unique scenario that gives this point. As a consequence fixing $\bar{P}_{a,\text{CO}_2}$ or $\bar{P}_{a,\text{O}_2}$ is not sufficient to discriminate between scenarios and other criteria shall be investigated.

<figure><img src="image_12.png" /><figcaption>FIGURE 11. $\bar{P}_{a,\text{O}_2}$ vs $\bar{P}_{a,\text{CO}_2}$ computed with the 0D model for a timeslot applied pressure. Each dot corresponds to a given scenario and is colored with respect to the time period $T$.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>FIGURE 12. $\bar{P}_{a,\text{O}_2}$ vs $\bar{P}_{a,\text{CO}_2}$ computed with the 0D model for a timeslot applied pressure. Each dot is colored with respect to inspiratory pressure $P_{\text{insp}}$.</figcaption></figure>

**0D model, piecewise constant applied pressure.** For the 0D model, there is no particular pattern with respect to the inspiratory fraction (see Fig. 13) whereas some emerge for the period and are appearing for the inspiratory pressure as well as for the pressure work. In Fig. 11 for decreasing carbon dioxide blood pressure (and increasing oxygen blood pressure), the average period decreases: for high gas exchange the 0D model favors short periods. The $\bar{P}_{a,\text{O}_2}-\bar{P}_{a,\text{CO}_2}$ curve is even more clearly structured with respect to the inspiratory pressure (see Fig. 12) or the pressure work (see Fig. 14). The general tendency is: the lower the applied pressure, the lower the pressure work and exchange rate.

**1D model, piecewise constant applied pressure.** Let us now consider the same cases for the 1D model. No structure emerges with respect to $T$, see Fig. 15, nor $\text{IE}$, see Fig. 17. No tendency appears when $\bar{P}_{a,\text{CO}_2}$ increases and $\bar{P}_{a,\text{O}_2}$ decreases. Yet for a given carbon dioxide blood partial pressure $\bar{P}_{a,\text{CO}_2}$, the oxygen blood partial pressure $\bar{P}_{a,\text{O}_2}$ increases as the period decreases. Moreover, when
<figure><img src="image_14.png" /><figcaption>FIGURE 13. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 0D model for a timeslot applied pressure. Each dot is colored with respect to inspiratory fraction IE.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>FIGURE 14. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 0D model for a timeslot applied pressure. Each dot is colored with respect to the pressure work $W$.</figcaption></figure>

following the curve from left to right in Fig. 16 and Fig. 18, in average, $|P_{\text{insp}}|$ and $W$ decrease respectively. Furthermore, for a given carbon dioxide blood pressure $\bar{P}_{a,CO_2}$, the oxygen blood pressure $\bar{P}_{a,O_2}$ increases as $|P_{\text{insp}}|$ increases and it comes at the cost of an increase of the pressure work $W$.

<figure><img src="image_16.png" /><figcaption>FIGURE 15. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D model for a timeslot applied pressure. Each dot is colored with respect to the time period $T$.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>FIGURE 16. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D model for a timeslot applied pressure. Each dot is colored with respect to the inspiratory pressure $|P_{\text{insp}}|$.</figcaption></figure>

**0D model, piecewise exponential applied pressure.** Let us consider the outputs of the 0D model for piecewise exponential scenario. In Figures 19-22, we observe the same behaviors as the ones observed for the piecewise constant applied pressure, except that for the amplitude of the applied pressure, less structures emerge (compare Fig. 20 with Fig. 12). In Fig. 22 we note that a given couple of ($\bar{P}_{a,CO_2}$, $\bar{P}_{a,O_2}$) corresponds to a fixed pressure work $W$. In average, the larger the work, the larger the gas exchanges.

**1D model, piecewise exponential applied pressure.** Finally, let us consider the outputs of the 1D model for piecewise exponential scenario. The observations are quite similar to the one
<figure><img src="image_18.png" /><figcaption>FIGURE 17. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D model for a timeslot applied pressure. Each dot is colored with respect to the inspiratory IE.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>FIGURE 18. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D model for a timeslot applied pressure. Each dot is colored with respect to the pressure work $W$.</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>FIGURE 19. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 0D Model for a piecewise exponential applied pressure. Each dot is colored with respect to the period $T$.</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>FIGURE 20. $P_{a,O_2}$ vs $P_{a,CO_2}$ computed with the 0D Model for a piecewise exponential applied pressure. Each dot is colored with respect to the inspiratory pressure $|P_{\text{insp}}|$.</figcaption></figure>

made previously: there is no pattern with respect to the inspiratory fraction IE (see Fig. 25), for a fixed carbon dioxide blood partial pressure $\bar{P}_{a,CO_2}$ any period will give an oxygen blood partial pressure $\bar{P}_{a,O_2}$ slightly varying around a given value, $\bar{P}_{a,O_2}$ increasing when the period $T$ decreases (see Fig. 23). Regarding the amplitude of the inspiratory pressure, the general tendency is, that the larger $|P_{\text{insp}}|$ is, the larger the gas exchanges are, as shown in Fig. 24. As before and as observed in Fig. 26 the curve remains structured with respect to the pressure work.

From these observations, in particular in the 1D case, if we have a target $\bar{P}_{a,CO_2}$ then the $\bar{P}_{a,O_2}$ is almost fixed and if we want to minimize the applied pressure work $W$ then a period $T$ is picked up but several values of the inspiratory fraction IE and amplitude pressure $|P_{\text{insp}}|$ may be admissible.

Thus, from these curves, it is unclear why the observed stereotypical scenario corresponds to an inspiratory fraction of 0.4. The aim of the next section is to explore which cost function could select the particular scenario observed in practice.
<figure><img src="image_22.png" /><figcaption>FIGURE 21. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 0D Model for a piecewise exponential applied pressure. Each dot is colored with respect to the inspiratory fraction $IE$.</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>FIGURE 22. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 0D Model for a piecewise exponential applied pressure. Each dot is colored with respect to the pressure work $W$.</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>FIGURE 23. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D Model for a piecewise exponential applied pressure. Each dot is colored with respect to the period $T$.</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>FIGURE 24. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D Model for a piecewise exponential applied pressure. Each dot is colored with respect to $|P_{\text{insp}}|$.</figcaption></figure>

3.2.3. *Optimal breathing scenario.* In this subsection, we investigate which cost function the observed standard breathing scenario is minimizing. As we already mentioned, there exist chemical sensors in the brain that measure the carbon dioxide blood partial pressure. Thus, the respiratory breathing pattern may adapt to maintain this partial pressure at a given level or in a given interval. In the case of our models, a given carbon dioxide blood partial pressure corresponds to a fixed oxygen blood partial pressure. Yet, as noted in the previous subsection, neither in the 0D case, nor in the 1D case, any particular scenario is selected, since there is no uniqueness of the scenario that provides a given output. Furthermore, minimizing the applied pressure work is not sufficient to discriminate between several scenarios. Consequently, as there are mechanical captors that are sensible to the lung elongation during inspiration [17], we consider an additional criterion which measures the difference between the maximal inhaled volume and the standard functional residual capacity, namely $D_V = (V_{\text{max}} - V_{\text{FRC}})^2$. As in the previous subsection, we consider four cases: the
<figure><img src="image_26.png" /><figcaption>FIGURE 25. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D Model for a piecewise exponential applied pressure. Each dot is colored with respect to IE.</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>FIGURE 26. $\bar{P}_{a,O_2}$ vs $\bar{P}_{a,CO_2}$ computed with the 1D Model for a piecewise exponential applied pressure. Each dot is colored with respect to $W$.</figcaption></figure>

0D model and the 1D model both for piecewise constant and piecewise exponential applied pressure. In the obtained series of Figures 27-29, 30-32, 33-35 and 36-38, the point of coordinates $(D_V, W)$ corresponds to one model and a fixed breathing scenario and each dot is colored with respect to the period $T$, the inspiratory pressure amplitude $|P_{\text{insp}}|$ or the inspiratory fraction IE associated to the specific scenario. We keep only the simulations for which the carbon dioxide blood partial pressure $\bar{P}_{a,CO_2}$ belongs to the range [39 - 41] mmHg.

**0D model, piecewise constant applied pressure.** In Fig. 27, 28, 29 the inputs that minimize both $D_V$ and $W$ among the simulations performed correspond to small periods (around 3 s), small amplitudes and rather small inspiratory fractions. These values do not correspond to the standard breathing scenario. Note that for the performed test cases as $|P_{\text{insp}}|$ increases, $D_V$ and $W$ increase as well.

<figure><img src="image_28.png" /><figcaption>FIGURE 27. $W$ vs $D_V$ computed with the 0D model for a timeslot applied pressure. Each dot is colored with respect to the time period $T$.</figcaption></figure>

<figure><img src="image_29.png" /><figcaption>FIGURE 28. $W$ vs $D_V$ computed with the 0D model for a timeslot applied pressure. Each dot is colored with respect to the inspiratory pressure $P_{\text{insp}}$.</figcaption></figure>

**1D model, piecewise constant applied pressure.** In Fig. 30, 31, 32 we observe that, among the performed test cases, the scenario which gives the smallest values of $W$ and $D_V$ is not the
<figure><img src="image_30.png" /><figcaption>FIGURE 29. $\mathcal{W}$ vs $\mathcal{D}_V$ computed with the 0D model for a timeslot applied pressure. Each dot is colored with respect to the inspiratory fraction $IE$.</figcaption></figure>

same as for 0D case. The selected inspiratory fraction and period are bigger than in the 0D case, and in particular it seems to select a value of $IE$ around 0.7, which is large compared to the standard acknowledge value of 0.4. As for the 0D model, we observe that $\mathcal{W}$ and $\mathcal{D}_V$ are correlated to $P_{\text{insp}}$. The differences in the 0D case and in the 1D case between the scenario for which $\bar{P}_{a,\text{CO}_2} \in [39, 41]$ mmHg and for $\mathcal{W}$ and $\mathcal{D}_V$ are minimal may be explained by the transit time induced by the bronchial tree. The 1D model can not ensure enough gas exchange when the period, the inspiratory pressure and inspiratory ratio are small, whereas the 0D model allows the exchanges to take place without time delay.

<figure><img src="image_31.png" /><figcaption>FIGURE 30. $\mathcal{W}$ vs $\mathcal{D}_V$ computed with the 1D model for a timeslot applied pressure. Each dot are colored with respect to the time period $T$.</figcaption></figure>

<figure><img src="image_32.png" /><figcaption>FIGURE 31. $\mathcal{W}$ vs $\mathcal{D}_V$ computed with the 1D model for a timeslot applied pressure. Each dot are colored with respect to the inspiratory pressure $P_{\text{insp}}$.</figcaption></figure>

**0D model, piecewise exponential applied pressure.** By considering the 0D model for various piecewise exponential scenario we observe in Fig. 34, Fig. 35 that, once again, to have the smallest values of $\mathcal{W}$ and $\mathcal{D}_V$ one has to select small inspiratory fraction $IE$ (yet larger than in the piecewise constant case), and small inspiratory pressure $|P_{\text{insp}}|$. Note that this inspiratory pressure is bigger than the one in the piecewise constant scenario. This can be explained by the pressure pattern and the need to sufficiently ventilate the lung to obtain a carbon dioxide blood partial pressure $\bar{P}_{a,\text{CO}_2}$ in the range [39 – 41] mmHg. Note moreover that the period for small $\mathcal{W}$ and $\mathcal{D}_V$ in Fig. 33 is bigger than the one for the piecewise constant scenario described in Fig. 27.
<figure><img src="image_33.png" /></figure>

FIGURE 32. $W$ vs $D_V$ computed with the 1D model for a timeslot applied pressure. Each dot are colored with respect to the inspiratory fraction $IE$.

<figure><img src="image_34.png" /></figure>

FIGURE 33. 0D Model $W$ vs $D_V$ for a piecewise exponential pressure, colored with respect to the period $T$.

FIGURE 34. 0D Model $W$ vs $D_V$ for a piecewise exponential pressure, colored with respect to the inspiratory pressure $|P_{\text{insp}}|$.

<figure><img src="image_35.png" /></figure>

FIGURE 35. 0D Model $W$ vs $D_V$ for a piecewise exponential pressure, colored with respect to the inspiratory fraction $IE$.

**1D model, piecewise exponential applied pressure.** Finally in Fig. 36, 37, 38 that consider the 1D model for various exponential piecewise scenarios, we observe that, to obtain the smallest values of $W$ and $D_V$, the period is around 5 s and the inspiratory fraction is higher than in the
0D case and around the desired value of 0.4. It was not the case of the 1D model with a piecewise constant scenario which favors larger inspiratory fraction.

<figure><img src="image_36.png" /><figcaption>FIGURE 36. 1D Model $W$ vs $D_V$ for a piecewise exponential pressure, colored with respect to the period $T$.</figcaption></figure>

<figure><img src="image_37.png" /><figcaption>FIGURE 37. 1D Model $W$ vs $D_V$ for a piecewise exponential pressure, colored with respect to the inspiratory pressure $|P_{\text{insp}}|$.</figcaption></figure>

As a general remark, it seems that the piecewise exponential pressure provide scenario for which the pressure work is smaller than when considering piecewise constant pressure, while maintaining proper level of gas exchanges (remember that the simulations correspond to a carbon dioxide blood partial pressure $\bar{P}_{a,CO_2}$ in the range [39 – 41] mmHg which ensures also a good level of $\bar{P}_{a,O_2}$.

**Optimal breathing pattern.** Next, in order to focus on the scenarios that enable to optimize the applied pressure work $W$ without stretching the lung too much, we select the tenth better scenarios in the following sense: we choose the following cost function, $\alpha D_V^2 + \beta W^2$ where both components are normalized in order to obtain the same order of magnitude. Tables 6-9 gather together the average and standard deviation of the input data $T$, $P_{\text{insp}}$, IE and output data $\bar{P}_{a,CO_2}$, $\bar{P}_{a,O_2}$, $V_T$ of the tenth better scenario with different target for the carbon dioxide blood partial pressure $\bar{P}_{a,CO_2}$. In Table 6, we consider only the $\bar{P}_{a,CO_2}$ between 39 and 41 mmHg for a piecewise constant applied pressure. The main differences between the 0D and 1D case are that the period and the inspiratory fraction are smaller in the 0D case. Neither models provide excepted values in terms of period, inspiratory fraction and tidal volume even if they provide standard values for the blood partial pressures. Next, considering the scenarios for which the cost function is minimal (among

<figure><img src="image_38.png" /><figcaption>FIGURE 38. 1D Model $W$ vs $D_V$ for a piecewise exponential pressure, colored with respect to the inspiratory fraction $\text{IE}$.</figcaption></figure>
<table><caption>TABLE 6. Timeslot pressure for a $\bar{P}_{a,CO_2}$ between 39 and 41 mmHg</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td><i>T</i> (s)</td><td>3.10 ± 0.22</td><td>3.64 ± 0.21</td></tr><tr><td><i>P</i><sub>insp</sub> (cmH<sub>2</sub>O)</td><td>-1.64 ± 0.08</td><td>-1.6 ± 0.0</td></tr><tr><td>IE</td><td>0.38 ± 0.08</td><td>0.59 ± 0.05</td></tr><tr><td><i>P</i>&#773;<sub>a,CO<sub>2</sub></sub> (mmHg)</td><td>40.84 ± 0.14</td><td>40.89 ± 0.06</td></tr><tr><td><i>P</i>&#773;<sub>a,O<sub>2</sub></sub> (mmHg)</td><td>93.42 ± 1.10</td><td>96.15 ± 0.63</td></tr><tr><td><i>V</i><sub>T</sub> (L)</td><td>0.36 ± 0.01</td><td>0.37 ± 0.01</td></tr><tr><td><i>W</i> (cmH<sub>2</sub>O · L · s<sup>-1</sup>)</td><td>0.19 ± 0.01</td><td>0.16 ± 0.01</td></tr><tr><td><i>D</i><sub>V</sub> (L<sup>2</sup>)</td><td>0.13 ± 0.01</td><td>0.16 ± 0.01</td></tr></tbody></table>

<table><caption>TABLE 7. Timeslot pressure for a $\bar{P}_{a,CO_2}$ under 40 mmHg</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td><i>T</i> (s)</td><td>3.02 ± 0.18</td><td>3.82 ± 0.15</td></tr><tr><td><i>P</i><sub>insp</sub> (cmH<sub>2</sub>O)</td><td>-1.82 ± 0.06</td><td>-1.8 ± 0.0</td></tr><tr><td>IE</td><td>0.41 ± 0.08</td><td>0.60 ± 0.07</td></tr><tr><td><i>P</i>&#773;<sub>a,CO<sub>2</sub></sub> (mmHg)</td><td>39.82 ± 0.14</td><td>39.84 ± 0.11</td></tr><tr><td><i>P</i>&#773;<sub>a,O<sub>2</sub></sub> (mmHg)</td><td>100.77 ± 0.89</td><td>103.03 ± 0.62</td></tr><tr><td><i>V</i><sub>T</sub> (L)</td><td>0.40 ± 0.01</td><td>0.42 ± 0.01</td></tr><tr><td><i>W</i> (cmH<sub>2</sub>O · L · s<sup>-1</sup>)</td><td>0.24 ± 0.01</td><td>0.20 ± 0.01</td></tr><tr><td><i>D</i><sub>V</sub> (L<sup>2</sup>)</td><td>0.17 ± 0.02</td><td>0.20 ± 0.01</td></tr></tbody></table>

<table><caption>TABLE 8. Piecewise exponential pressure for a $\bar{P}_{a,CO_2}$ between 39 and 41 mmHg</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td><i>T</i> (s)</td><td>5.08 ± 0.41</td><td>5.06 ± 0.67</td></tr><tr><td><i>P</i><sub>insp</sub> (cmH<sub>2</sub>O)</td><td>-2.66 ± 0.10</td><td>-2.30 ± 0.0</td></tr><tr><td>IE</td><td>0.34 ± 0.02</td><td>0.41 ± 0.02</td></tr><tr><td><i>P</i>&#773;<sub>a,CO<sub>2</sub></sub> (mmHg)</td><td>40.91 ± 0.03</td><td>40.92 ± 0.05</td></tr><tr><td><i>P</i>&#773;<sub>a,O<sub>2</sub></sub> (mmHg)</td><td>92.76 ± 0.28</td><td>95.37 ± 0.53</td></tr><tr><td><i>V</i><sub>T</sub> (L)</td><td>0.51 ± 0.03</td><td>0.47 ± 0.03</td></tr><tr><td><i>W</i> (cmH<sub>2</sub>O · L · s<sup>-1</sup>)</td><td>0.11 ± 0.01</td><td>0.09 ± 0.01</td></tr><tr><td><i>D</i><sub>V</sub> (L<sup>2</sup>)</td><td>0.27 ± 0.02</td><td>0.25 ± 0.02</td></tr></tbody></table>

the performed simulations) and for which $\bar{P}_{a,CO_2}$ is smaller than 40 mmHg, we observe in Table 7 that we have the same type of scenarios but with an average value of $\bar{P}_{a,O_2}$ bigger than previously. This is not surprising since we select simulations for which $\bar{P}_{a,CO_2} \le 40$ and thus with, possibly, better exchanges than in the case $\bar{P}_{a,CO_2} \in [39, 41]$.

Next we present the same results in the case of a piecewise exponential applied pressure. We observe in Tables 8 and 9 that the selected scenario in the 1D case corresponds almost to the acknowledged standard values, both for the inputs and the outputs. It is not the case for the 0D case for which IE is still below the expected value in the case where the $\bar{P}_{a,CO_2}$ is between 39 and 41 mmHg.

Let us now change the cost function by adding a constraint on $\bar{P}_{a,CO_2}$ in order to minimize its distance to 40 mmHg. The new considered cost function is the same as before but we add $(\bar{P}_{a,CO_2} - 40)^2$. Tables 10 and 11 show, again, that, for a piecewise constant pressure, the selected IE is lower than the standard value for the 0D model and bigger for the 1D model. Note moreover
<table><caption>TABLE 9. Piecewise exponential pressure for a $\bar{P}_{a,CO_2}$ under 40 mmHg</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td><i>T</i> (s)</td><td>4.96 ± 0.48</td><td>5.08 ± 0.64</td></tr><tr><td><i>P</i><sub>insp</sub> (cmH2O)</td><td>-2.90 ± 0.14</td><td>-2.50 ± 0.0</td></tr><tr><td>IE</td><td>0.41 ± 0.06</td><td>0.46 ± 0.02</td></tr><tr><td><i>P̅</i><sub>a,CO<sub>2</sub></sub> (mmHg)</td><td>39.91 ± 0.05</td><td>39.88 ± 0.08</td></tr><tr><td><i>P̅</i><sub>a,O<sub>2</sub></sub> (mmHg)</td><td>100.14 ± 0.36</td><td>102.57 ± 0.53</td></tr><tr><td><i>V</i><sub>T</sub> (L)</td><td>0.58 ± 0.03</td><td>0.52 ± 0.04</td></tr><tr><td><i>W</i> (cmH<sub>2</sub>O · L · s<sup>-1</sup>)</td><td>0.14 ± 0.01</td><td>0.11 ± 0.01</td></tr><tr><td><i>D</i><sub>V</sub> (L<sup>2</sup>)</td><td>0.37 ± 0.04</td><td>0.32 ± 0.02</td></tr></tbody></table>

<table><caption>TABLE 10. Timeslot pressure</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td><i>T</i> (s)</td><td>3.10 ± 0.25</td><td>3.94 ± 0.27</td></tr><tr><td><i>P</i><sub>insp</sub> (cmH2O)</td><td>-1.84 ± 0.08</td><td>-1.8 ± 0.0</td></tr><tr><td>IE</td><td>0.36 ± 0.08</td><td>0.61 ± 0.09</td></tr><tr><td><i>P̅</i><sub>a,CO<sub>2</sub></sub> (mmHg)</td><td>40.03 ± 0.12</td><td>40.00 ± 0.07</td></tr><tr><td><i>P̅</i><sub>a,O<sub>2</sub></sub> (mmHg)</td><td>99.41 ± 0.80</td><td>102.08 ± 0.56</td></tr><tr><td><i>V</i><sub>T</sub> (L)</td><td>0.40 ± 0.02</td><td>0.42 ± 0.01</td></tr><tr><td><i>W</i> (cmH<sub>2</sub>O · L · s<sup>-1</sup>)</td><td>0.23 ± 0.02</td><td>0.19 ± 0.01</td></tr><tr><td><i>D</i><sub>V</sub> (L<sup>2</sup>)</td><td>0.16 ± 0.02</td><td>0.20 ± 0.01</td></tr></tbody></table>

<table><caption>TABLE 11. Piecewise exponential pressure</caption><thead><tr><th></th><th>0D Model</th><th>1D Model</th></tr></thead><tbody><tr><td><i>T</i> (s)</td><td>5.00 ± 0.52</td><td>5.10 ± 0.61</td></tr><tr><td><i>P</i><sub>insp</sub> (cmH2O)</td><td>-2.82 ± 0.06</td><td>-2.50 ± 0.0</td></tr><tr><td>IE</td><td>0.39 ± 0.02</td><td>0.42 ± 0.03</td></tr><tr><td><i>P̅</i><sub>a,CO<sub>2</sub></sub> (mmHg)</td><td>40.09 ± 0.10</td><td>40.10 ± 0.12</td></tr><tr><td><i>P̅</i><sub>a,O<sub>2</sub></sub> (mmHg)</td><td>98.93 ± 0.69</td><td>101.21 ± 0.78</td></tr><tr><td><i>V</i><sub>T</sub> (L)</td><td>0.57 ± 0.04</td><td>0.52 ± 0.04</td></tr><tr><td><i>W</i> (cmH<sub>2</sub>O · L · s<sup>-1</sup>)</td><td>0.13 ± 0.01</td><td>0.10 ± 0.01</td></tr><tr><td><i>D</i><sub>V</sub> (L<sup>2</sup>)</td><td>0.35 ± 0.03</td><td>0.30 ± 0.03</td></tr></tbody></table>

that the tidal volume is lower than the standard value. In the case of a piecewise exponential applied pressure we obtain a strong good fit with the standard values for the 1D model, which is not the case of the 0D model as the tidal volume is bigger and thus the associated pressure work is bigger than the one obtained in the 1D case.

## 4. CONCLUSION AND PERSPECTIVES

We have presented a new 1D-0D coupled system developed to describe the transport and diffusion of oxygen and carbon dioxide in the bronchial tree through the blood. This model is a coupled non-linear model that enables to recover acknowledged standard values of gas fluxes and blood partial pressures without any fitting of the input physiological parameters. This model takes naturally into account the time delay induced by the air transport along the bronchial tree.

We have compared it to its 0D counterpart. For standard breathing scenarios it appears that the time delay favors the gas exchanges. It can be explained by the fact that for the 0D model
<table><caption>TABLE 12. Parameters</caption><thead><tr><td>Physiological model parameters</td><td>Normal breathing values</td><td>References</td></tr></thead><tbody><tr><td>Elastance (E)</td><td>4 cmH<sub>2</sub>O · L<sup>-1</sup></td><td>[15, 11]</td></tr><tr><td>Resistance (R)</td><td>2 cmH<sub>2</sub>O · s · L<sup>-1</sup></td><td>[15, 35]</td></tr><tr><td>Functional Residual Capacity (V<sub>FRC</sub>)</td><td>2.5 L</td><td>[15, 17]</td></tr><tr><td>Reduced atmospheric pressure (P<sub>atm</sub>)</td><td>713 mmHg</td><td>[11]</td></tr><tr><td>Diffusion coefficient of O<sub>2</sub> in air (D<sub>O<sub>2</sub></sub>)</td><td>0.2 cm<sup>2</sup> · s<sup>-1</sup></td><td>[19]</td></tr><tr><td>Diffusion coefficient of CO<sub>2</sub> in air (D<sub>CO<sub>2</sub></sub>)</td><td>0.14 cm<sup>2</sup> · s<sup>-1</sup></td><td>[19]</td></tr><tr><td>Mole fraction in fresh air for O<sub>2</sub> (χ&#773;<sub>O<sub>2</sub></sub>)</td><td>21%</td><td>[11]</td></tr><tr><td>Mole fraction in fresh air for CO<sub>2</sub> (χ&#773;<sub>CO<sub>2</sub></sub>)</td><td>0.04%</td><td>[11]</td></tr><tr><td>Membrane diffusing capacity of O<sub>2</sub> (D<sub>m,O<sub>2</sub></sub>)</td><td>21 mL · min<sup>-1</sup> · mmHg<sup>-1</sup></td><td>[11]</td></tr><tr><td>Membrane diffusing capacity of CO<sub>2</sub> (D<sub>m,CO<sub>2</sub></sub>)</td><td>420 mL · min<sup>-1</sup> · mmHg<sup>-1</sup></td><td>[11]</td></tr><tr><td>Capillary volume (V<sub>c</sub>)</td><td>70 mL</td><td>[11, 35]</td></tr><tr><td>Transit time of red blood cells in V<sub>c</sub> (τ<sub>b</sub>)</td><td>0.75 s</td><td>[35]</td></tr><tr><td>Pressure of O<sub>2</sub> in deoxygenated blood (P<sub>v,O<sub>2</sub></sub>)</td><td>40 mmHg</td><td>[11]</td></tr><tr><td>Pressure of CO<sub>2</sub> in deoxygenated blood (P<sub>v,CO<sub>2</sub></sub>)</td><td>45 mmHg</td><td>[11]</td></tr></tbody></table>

the equilibrium between alveolar pressures and blood pressures is reached more rapidly whereas the transport in the bronchial tree maintains a gap between the alveolar and blood partial pressures.

We have also investigated the sensitivity of both models with respect to breathing scenarios by considering two types of pleural applied pressure: a piecewise constant and a piecewise exponential (modeling the action of the respiratory muscles). As no time delay is modelled in the 0D model the physiological gas exchanges can be reached for small breathing period, small inspiratory ratio and small inspiratory pressure amplitude. It is not the case for the 1D model for which the period and the inspiratory ratio have to be larger. Furthermore for both models we have observed that a given carbon dioxide arterial partial pressure is associated to a unique or almost unique oxygen arterial partial pressure. However, numerous scenarios provide these pairs. This non uniqueness makes the breathing process quite robust with respect to the breathing pattern yet it leads to a natural question: which criteria enables to select a particular breathing pattern? To answer this question, we have explored several cost functions: the pressure work, the total volume variation of the lung to which we further add the deviation from a targeted carbon dioxide arterial partial pressure. In the present study, the stereotypical scenario appears to be a piecewise exponential one, minimizing weighted linear combination of the three mentioned quantities, underlying the role of carbon dioxide captors, lung distension together with low effort in the breathing scenario selection. The proposed 1D system could then further be adapted to model pathologies such as asthma, emphysema or fibrosis and understand their impact on breathing patterns, output data and possible minimizing criteria. It could also be useful to explore the effect of mechanical ventilation and optimize the delivered applied pressure. In particular, in the provided test cases, piecewise constant applied pressure (as delivered by mechanical ventilators) seems not to be optimal and exponential pleural applied pressure seems to favors gas exchanges with lower lung extension and applied pressure work and, thus, shall induce lower lung injuries.

## APPENDIX A. NUMERICAL DATA

### REFERENCES

[1] A. Albanese, L. Cheng, M. Ursino, and N. W. Chbat. An integrated mathematical model of the human car-diopulmonary system: model development. *American Journal of Physiology-Heart and Circulatory Physiology*, 310(7):H899-H921, 2016. PMID: 26683899.
[2] L. Baffico, C. Grandmont, and B. Maury. Multiscale modeling if the respiratory tract. *Mathematical Models and Methods in Applied Sciences*, 20(01):59–93, 2010.

[3] J. H. T. Bates. *Lung Mechanics: An Inverse Modeling Approach*. Cambridge University Press, Cambridge, 2009.

[4] A. Ben-Tal. Simplified models for gas exchange in the human lungs. *Journal of Theoretical Biology*, 238(2):474–495, 2006.

[5] L. Berger, R. Bordas, K. Burrowes, V. Grau, S. Tavener, and D. Kay. A poroelastic model coupled to a fluid network with applications in lung modelling. *International journal for numerical methods in biomedical engineering*, 32(1), 2016.

[6] L. Boudin, C. Grandmont, B. Grec, and S. Martin. A coupled model for the dynamics of gas exchanges in the human lung with Haldane and Bohr's effects. *Journal of Theoretical Biology*, 573:111590, 2023.

[7] L. Cheng, O. Ivanova, H.-H. Fan, and M. C. K. Khoo. An integrative model of respiratory and cardiovascular control in sleep-disordered breathing. *Respiratory Physiology & Neurobiology*, 174(1):4–28, Nov. 2010.

[8] L. de Rochefort, L. Vial, R. Fodil, X. Maître, B. Louis, D. Isabey, G. Caillibotte, M. Thiriet, J. Bittoun, E. Durand, and G. Sbirlea-Apiou. In vitro validation of computational fluid dynamic simulation in human proximal airways with hyperpolarized 3He magnetic resonance phase-contrast velocimetry. *Journal of Applied Physiology*, 102(5):2012–2023, 2007. PMID: 17289906.

[9] T. Gemci, V. Ponyavin, Y. Chen, H. Chen, and R. Collins. Computational model of airflow in upper 17 generations of human respiratory tract. *Journal of Biomechanics*, 41(9):2047–2054, 2008.

[10] F. S. Grodins, J. Buell, and A. J. Bart. Mathematical analysis and digital simulation of the respiratory control system. *Journal of Applied Physiology*, 22(2):260–276, Feb. 1967.

[11] A. C. Guyton and J. E. Hall. *Textbook of Medical Physiology*. Elsevier Saunders, Philadelphia, 11th edition, 2006.

[12] R. P. Hämäläinen and A. A. Viljanen. Modelling the respiratory airflow pattern by optimization criteria. *Biological Cybernetics*, 29(3):143–149, 1978.

[13] M. Ismail, A. Comerford, and W. Wall. Coupled and reduced dimensional modeling of respiratory mechanics during spontaneous breathing. *International Journal for Numerical Methods in Biomedical Engineering*, 29(11):1285–1305, 2013.

[14] F. Jiang, T. Hirano, C. Liang, G. Zhang, K. Matsunaga, and X. Chen. Multi-scale simulations of pulmonary airflow based on a coupled 3d-1d-0d model. *Computers in Biology and Medicine*, 171:108150, 2024.

[15] A. T. Johnson. *Biomechanics and Exercise Physiology: Quantitative Modeling*. CRC Press, Boca Raton, 2007.

[16] C. Karamaoun, A. Van Muylem, and B. Haut. Modeling of the nitric oxide transport in the human lungs. *Frontiers in Physiology*, 7:255, 2016.

[17] B. M. Koeppen and B. A. Stanton, editors. *Berne & Levy Physiology*. Elsevier, Philadelphia, 7th edition, 2018.

[18] S. Martin and B. Maury. Modeling of the oxygen transfer in the respiratory process. *ESAIM: Mathematical Modelling and Numerical Analysis*, 47(4):935–960, 2013.

[19] B. Mauroy. *Hydrodynamique dans le poumon, relations entre flux et géométries*. PhD thesis, École normale supérieure de Cachan, Paris, 2004.

[20] B. Maury. *The Respiratory System in Equations*. Number 7 in MS&A - Modeling, Simulations & Applications. Springer-Verlag, 2013.

[21] J. Mead. Control of respiratory frequency. *Journal of Applied Physiology*, 15(3):325–336, 1960.

[22] C. Ngo, S. Dahlmanns, T. Vollmer, B. Misgeld, and S. Leonhardt. An object-oriented computational model to study cardiopulmonary hemodynamic interactions in humans. *Computer Methods and Programs in Biomedicine*, 159:167–183, 2018.

[23] F. Noël and B. Mauroy. Interplay between optimal ventilation and gas transport in a model of the human lung. *Frontiers in Physiology*, 10:488, 2019.

[24] A. B. Otis, W. O. Fenn, and H. Rahn. Mechanics of breathing in man. *Journal of Applied Physiology*, 2(11):592–607, 1950.

[25] M. Paiva and L. A. Engel. Theoretical studies of gas mixing and ventilation distribution in the lung. *Physiological Reviews*, 67(3):750–796, 1987.

[26] C. Patte, M. Genet, and D. Chapelle. A quasi-static poromechanical model of the lungs. *Biomechanics and Modeling in Mechanobiology*, 21(2):527–551, 2022.

[27] N. Pozin, S. Montesantos, I. Katz, M. Pichelin, I. Vignon-Clementel, and C. Grandmont. A tree-parenchyma coupled model for lung ventilation simulation. *International Journal for Numerical Methods in Biomedical Engineering*, 33(11):e2873, 2017. e2873 CNM-Sep-16-0190.R1.

[28] B. Sapoval and M. Filoche. Role of diffusion screening in pulmonary diseases. *Advances in Experimental Medicine and Biology*, 605:173–178, 2008.

[29] B. Sapoval, M. Filoche, and E. R. Weibel. Smaller is better – but not too small: A physical scale for the design of the mammalian pulmonary acinus. *Proceedings of the National Academy of Sciences of the United States of America*, 99:10411–10416, 2002.
[30] K. B. Saunders, H. N. Bali, and E. R. Carson. A breathing model of the respiratory system: The controlled system. *Journal of Theoretical Biology*, 84(1):135-161, May 1980.

[31] A. J. Swan, A. R. Clark, and M. H. Tawhai. A computational model of the topographic distribution of ventilation in healthy human lungs. *Journal of Theoretical Biology*, 300:222-231, 2012.

[32] A. J. Swan and M. H. Tawhai. Evidence for minimal oxygen heterogeneity in the healthy human pulmonary acinus. *Journal of Applied Physiology*, 110(2):528-537, 2011.

[33] A. Van Muylem, C. Noël, and M. Paiva. Modeling of impact of gas molecular diffusion on nitric oxide expired profile. *Journal of Applied Physiology (Bethesda, Md.: 1985)*, 94(1):119-127, Jan. 2003.

[34] E. R. Weibel. *Morphometry of the Human Lung*. Academic Press, New York, 1963.

[35] J. B. West. *Respiratory Physiology: The Essentials*. Lippincott Williams and Wilkins, Philadelphia, 9th edition, 2011.

S.M. benefited from the support of the France 2030 program through the IdEx Université Paris Cité (ANR 18 IDEX 0001).

<footnote n="1">¹ LJLL, SORBONNE UNIVERSITÉ, CNRS, & INRIA PARIS, FRANCE, & ULB, BELGIUM
*Email address: celine.grandmont@inria.fr*</footnote>

<footnote n="2">²TIPS, UNIVERSITÉ LIBRE DE BRUXELLES, BRUSSELS, BELGIUM
*Email address: cyril.karamaoun@ulb.be*</footnote>

<footnote n="3">³MAP5, UNIVERSITÉ PARIS CITÉ, CNRS, F-75006 PARIS - FRANCE
*Email address: sebastien.martin@u-paris.fr*</footnote>

<footnote n="4">⁴LJLL, SORBONNE UNIVERSITÉ, CNRS, & INRIA PARIS
*Email address: frederique.noel@sorbonne-universite.fr*</footnote>