Renewal age and size structured numerical methods

Fix a finite horizon $T>0$. Let mortality and fertility be bounded Lipschitz functions $\mu ,\beta :[0,\mathrm{\infty })\to [0,\mathrm{\infty })$, and let $n_0\in L^1(0,\mathrm{\infty })$ be nonnegative. The linear age-structured renewal transport problem is

$n_t(t,a)+n_a(t,a)=-\mu (a)n(t,a),a>0,0<t\le T,$

with initial density $n(0,a)=n_0(a)$ and inflow condition

$n(t,0)=B(t),B(t):={\int }_0^{\mathrm{\infty }}\beta (a)n(t,a)da\mathrm{.}$

Here $a$ is age, $B$ is the birth rate, and the characteristic age path solves the first-order initial value problem $a^{\mathrm{\prime }}(s)=1$ in the sense of First-order systems, initial value problems, and solutions on intervals, hence $a(s)=a(0)+s$. Define survival from age $a$ over duration $s\ge 0$ by

$S(a,s):=\exp \text{\hspace{0.17em}⁣}(-{\int }_a^{a+s}\mu (r)dr)\mathrm{.}$

A weak solution on $[0,T]$ is a nonnegative $n\in C([0,T];L^1(0,\mathrm{\infty }))$ with continuous $B:[0,T]\to [0,\mathrm{\infty })$ satisfying $B(t)={\int }_0^{\mathrm{\infty }}\beta (a)n(t,a)da$ and, for every $\varphi \in C_c^1([0,T]\times [0,\mathrm{\infty }))$ with $\varphi (T,a)=0$,

$-{\int }_0^T\text{\hspace{0.17em}⁣}{\int }_0^{\mathrm{\infty }}n(t,a)({\varphi }_t+{\varphi }_a-\mu (a)\varphi )dadt={\int }_0^{\mathrm{\infty }}{n}_0(a)\varphi (0,a)da+{\int }_0^{T}B(t)\varphi (t,0)dt\mathrm{.}$

A mild solution is a weak solution that also satisfies, for almost every $a\ge 0$,

$n(t,a)=\{\begin{array}{ll}{n}_0(a-t)S(a-t,t),& a\ge t,\\ B(t-a)S(0,a),& 0\le a<t,\end{array}$

The formula follows the age characteristics $a^{\mathrm{\prime }}(s)=1$ and their mortality survival. No delay equation or finite-age cutoff is identified with this model by definition.

Fix a finite horizon $T>0$ and a nondegenerate fixed size interval $[0,X]$ with $X>0$. Let $g,\mu ,\gamma :[0,X]\to [0,\mathrm{\infty })$ be bounded Lipschitz functions, with $g(x)\ge g_{\min }>0$. Here $g$ is size growth, $\mu$ is ordinary mortality, and $\gamma$ is the rate at which a dividing parent is removed. For each parent size $y\in [0,X]$, let $K_y$ be a finite nonnegative Borel measure of daughter sizes, supported on $[0,y]$. Assume $y\mapsto K_y(E)$ is measurable for every Borel set $E\subseteq [0,X]$, and that the number moment $k_0(y):=K_y([0,y])$ is uniformly bounded; the size moment $k_1(y):={\int }_{[0,y]}x{K}_y(dx)$ is then finite.

Because $K_y$ may have atoms, represent the population by a finite nonnegative Borel measure ${\nu }_t$ on $[0,X)$ rather than requiring a density. Given an initial measure ${\nu }_0$, a weak solution is a curve $t\mapsto {\nu }_t$ that is continuous when tested against continuous functions on $[0,X]$ vanishing at $X$ and satisfies, for every $\varphi \in C^1([0,T]\times [0,X])$ with $\varphi (t,X)=0$,

${\int }_{[0,X)}\varphi (T,x){\nu }_T(dx)-{\int }_{[0,X)}\varphi (0,x){\nu }_0(dx)={\int }_0^T\text{\hspace{0.17em}⁣}{\int }_{[0,X)}({\varphi }_t+g{\varphi }_x-(\mu +\gamma )\varphi )(t,x){\nu }_t(dx)dt+{\int }_0^T\text{\hspace{0.17em}⁣}{\int }_{[0,X)}\text{\hspace{0.17em}⁣}({\int }_{[0,y]}\varphi (t,x)K_y(dx)){\nu }_t(dy)dt\mathrm{.}$

This is the measure-valued transport equation with daughter production ${\int }_{[0,X)}{K}_y(dx){\nu }_t(dy)$, no incoming source at $x=0$, and absorbing outflow at $x=X$. The test condition $\varphi (t,X)=0$ encodes removal when a characteristic reaches the upper boundary.

Let ${\mathrm{\Pi }}_{[0,X]}(x):=\min \{X,\max \{0,x\}\}$ and extend the growth field to $\tilde{g}(x):=g({\mathrm{\Pi }}_{[0,X]}(x))$ on $\mathbb{R}$. This is a bounded globally Lipschitz vector field. By Size characteristics, positivity, and discrete number and size-moment balances ↗, its IVP has a unique solution ${\chi }_y$ on $[0,T]$ from every initial size $y\in [0,X]$. A model characteristic is ${\chi }_y$ restricted to $[0,\min \{T,{\tau }_y\}]$, where ${\tau }_y$ is its first hitting time of $X$ (taking ${\tau }_y=+\mathrm{\infty }$ if there is no hit on $[0,T]$, and ${\tau }_y=0$ if $y=X$). On this restricted interval, $\tilde{g}=g$. For $y<X$ the restricted interval is nondegenerate, and the characteristic is a first-order IVP solution in the sense of First-order systems, initial value problems, and solutions on intervals. For $y=X$ the restriction is only the stopped path on $\{0\}$; no derivative on that singleton is asserted, while the extended path ${\chi }_X$ remains an IVP solution on $[0,T]$. The interval is fixed; moving boundaries, zero growth and discontinuous coefficients are outside this definition.

For Linear age-structured renewal transport problem, put $M(a)={\int }_0^a\mu (r)dr$, $K(a)=\beta (a)e^{-M(a)}$, and

$G(t)={\int }_0^{\mathrm{\infty }}\beta (x+t)\exp \text{\hspace{0.17em}⁣}(-{\int }_x^{x+t}\mu (r)dr)n_0(x)dx\mathrm{.}$

There is a unique mild solution. Its birth rate is the unique continuous solution of

$B(t)=G(t)+{\int }_0^{t}K(a)B(t-a)da,$

and its density, for almost every $a\ge 0$, is

$n(t,a)=\{\begin{array}{ll}{n}_0(a-t)\exp \text{\hspace{0.17em}⁣}(-{\int }_{a-t}^a\mu (r)dr),& a\ge t,\\ B(t-a)e^{-M(a)},& 0\le a<t\mathrm{.}\end{array}$

Fix $0<T<A$, bounded nonnegative Lipschitz rates $\mu ,\beta :[0,A]\to [0,\mathrm{\infty })$, and initial density $n_0^A\in L_{+}^1(0,A)$. The finite-age renewal model with absorbing outflow is the age transport problem on $0<a<A$ with the same renewal inflow $B_A(t)={\int }_0^A\beta (a)n^A(t,a)da$ as in Linear age-structured renewal transport problem. A cohort follows age speed $1$, retains its mortality survival $S(a,s)$ until it reaches $A$, and is then removed; no cohort enters from ages above $A$. The initial density is $n_0^A$ on $[0,A]$. A mild solution is the characteristic formula on $0<a<A$, with an initial-cohort branch for $a\ge t$, a newborn branch for $a<t$, and zero extension to $a\ge A$. Since $A>T$, no cohort born during $[0,T]$ reaches the outflow boundary.

Input a partition $0=x_0<x_1<\cdots <x_J=X$ with $J\ge 1$, a time step $\mathrm{\Delta }t>0$, a nonnegative current mass vector $M^k=(M_j^k{)}_{j=1}^J$, and the coefficients and kernel from Size transport with a nonlocal division kernel. Partition $[0,X]$ into disjoint cells $C_j=[x_{j-1},x_j)$ for $j<J$ and $C_J=[x_{J-1},X]$, with widths $h_j=x_j-x_{j-1}>0$ and centers $x_j^c=(x_{j-1}+x_j)\mathrm{/}2$. Set $u_j^k=M_j^k\mathrm{/}{h}_j$, the incoming flux $F_0^k=0$, the internal and outgoing upwind fluxes $F_j^k=g(x_j)u_j^k$ for $1\le j\le J$, and $d_j=\mu (x_j^c)+\gamma (x_j^c)$.

For each parent cell $C_{\mathrm{\ell }}$, form nonnegative daughter coefficients

$A_{j\mathrm{\ell }}=\frac{1}{h_{\mathrm{\ell }}}{Q}_M(y\mapsto K_y(C_j);x_{\mathrm{\ell }-1},x_{\mathrm{\ell }},1),$

using Composite midpoint rule on a uniform grid on the parent-size interval. Since the one-panel midpoint rule is $Q_M(f;x_{\mathrm{\ell }-1},x_{\mathrm{\ell }},1)=h_{\mathrm{\ell }}f(x_{\mathrm{\ell }}^c)$, this gives $A_{j\mathrm{\ell }}=K_{x_{\mathrm{\ell }}^c}(C_j)\ge 0$. If $\mathrm{\Delta }t(g(x_j)\mathrm{/}{h}_j+d_j)>1$ for any cell, return invalid-cfl without an update. Otherwise the finite-volume update is

$M_j^{k+1}=M_j^k-\mathrm{\Delta }t(F_j^k-F_{j-1}^k)-\mathrm{\Delta }t{d}_j{M}_j^k+\mathrm{\Delta }t{\sum }_{\mathrm{\ell }=1}^J{A}_{j\mathrm{\ell }}{M}_{\mathrm{\ell }}^k\mathrm{.}$

The flux $F_J^k$ is absorbing outflow. The admissible variable-width CFL is $\mathrm{\Delta }t(g(x_j)\mathrm{/}{h}_j+d_j)\le 1$ for every cell. If a grid or any required coefficient/kernel evaluation fails, return failure without a partial update. Otherwise return $M^{k+1}$ together with the fluxes and daughter coefficients used in the step.

Input a finite horizon $T>0$, a time step $\mathrm{\Delta }t>0$, an age cutoff $A>0$, bounded nonnegative rates $\mu ,\beta$ on $[0,A]$, and a finite positive initial quadrature $(a_i,m_i^0{)}_{i=1}^I$ supported in $[0,A)$ for the initial age measure on that interval. Set $N=\lfloor T\mathrm{/}\mathrm{\Delta }t\rfloor$; the procedure returns levels at $t_k=k\mathrm{\Delta }t$ for $0\le k\le N$, with final time $N\mathrm{\Delta }t\le T$. One option, when a nonnegative real-valued initial density $n_0$ is evaluable at every point of $[0,A]$, is to choose any finite partition $0=b_0<\cdots <b_I=A$ and apply Composite midpoint rule on a uniform grid with $n=1$ separately on each cell $[b_{i-1},b_i]$: set $\mathrm{\Delta }{a}_i=b_i-b_{i-1}$, $a_i=(b_{i-1}+b_i)\mathrm{/}2$, and $m_i^0=Q_M(n_0;b_{i-1},b_i,1)=\mathrm{\Delta }{a}_i{n}_0(a_i)$ for $1\le i\le I$. The cell widths need not be equal. A cohort at age $a_i^k$ and mass $m_i^k$ follows the exact age characteristic

$a_i^{k+1}=a_i^k+\mathrm{\Delta }t$

and updates its survival mass by forward Euler,

$m_i^{k+1}=(1-\mathrm{\Delta }t\mu (a_i^k))m_i^k\mathrm{.}$

If $\mathrm{\Delta }t\mathrm{\parallel }\mu {\mathrm{\parallel }}_{\mathrm{\infty }}>1$, return status invalid-step with no updates. Otherwise initialize the birth and outflow records as empty and set ${\nu }_h^0={\sum }_i{m}_i^0{\delta }_{a_i^0}$. For each $k=0,\dots ,N-1$, first calculate $B_h^k={\sum }_i\beta (a_i^k)m_i^k$. Advance every existing cohort by $a_i^{k+1}=a_i^k+\mathrm{\Delta }t$ and update its survival mass by $m_i^{k+1}=(1-\mathrm{\Delta }t\mu (a_i^k))m_i^k$. If an updated age is at least $A$, record that updated mass as discrete outflow and discard the cohort; otherwise retain it. Then insert a newborn cohort of age $0$ and mass $\mathrm{\Delta }t{B}_h^k$ at time $t_{k+1}$. Record $B_h^k$ and the step's outflow, and form ${\nu }_h^{k+1}$ from the retained cohorts and this newborn. A failed rate evaluation returns failure with no partial update. The procedure returns the finite sequence $({\nu }_h^k{)}_{k=0}^N$ and its birth/outflow records. The bounded-Lipschitz quadrature defect of the initial measure and the Euler survival error are separate inputs to any consistency claim.

For Linear age-structured renewal transport problem, let ${\beta }_{\mathrm{\infty }}=\mathrm{\parallel }\beta {\mathrm{\parallel }}_{\mathrm{\infty }}$. The solution map $S_{\mu ,\beta }(t):n_0\mapsto n(t,\cdot )$ is linear and positive, and

$\mathrm{\parallel }{S}_{\mu ,\beta }(t)n_0{\mathrm{\parallel }}_1\le e^{{\beta }_{\mathrm{\infty }}t}\mathrm{\parallel }{n}_0{\mathrm{\parallel }}_1\mathrm{.}$

If $n,\tilde{n}$ use possibly different nonnegative bounded coefficients $(\mu ,\beta )$ and $(\tilde{\mu },\tilde{\beta })$ and initial data $n_0,{\tilde{n}}_0$, put $\delta =\mathrm{\parallel }\mu -\tilde{\mu }{\mathrm{\parallel }}_{\mathrm{\infty }}+\mathrm{\parallel }\beta -\tilde{\beta }{\mathrm{\parallel }}_{\mathrm{\infty }}$ and ${\tilde{\beta }}_{\mathrm{\infty }}=\mathrm{\parallel }\tilde{\beta }{\mathrm{\parallel }}_{\mathrm{\infty }}$. Then

$\mathrm{\parallel }n(t)-\tilde{n}(t){\mathrm{\parallel }}_1\le e^{{\beta }_{\mathrm{\infty }}t}(\mathrm{\parallel }{n}_0-{\tilde{n}}_0{\mathrm{\parallel }}_1+\delta t{e}^{{\tilde{\beta }}_{\mathrm{\infty }}t}\mathrm{\parallel }{\tilde{n}}_0{\mathrm{\parallel }}_1)\mathrm{.}$

Under Size transport with a nonlocal division kernel, for every initial size $y\in [0,X]$ the extended IVP specified there has a unique solution ${\chi }_y$ on $[0,T]$. The solution is strictly increasing and agrees with the original growth field until it first reaches $X$, so its restriction through that first hit is the model characteristic.

For Finite-volume size transport with positive division quadrature, assume nonnegative initial masses, nonnegative quadrature coefficients $A_{j\mathrm{\ell }}$, and the stated CFL. Then all $M_j^k$ remain nonnegative. Define the cell number and center-weighted size moment by $N_h^k={\sum }_j{M}_j^k$ and $W_h^k={\sum }_j{x}_j^c{M}_j^k$. With $k_{0,\mathrm{\ell }}^h={\sum }_j{A}_{j\mathrm{\ell }}$ and $k_{1,\mathrm{\ell }}^h={\sum }_j{x}_j^c{A}_{j\mathrm{\ell }}$, the exact balances are

$N_h^{k+1}-N_h^k=\mathrm{\Delta }t(-F_J^k-{\sum }_j{d}_j{M}_j^k+{\sum }_{\mathrm{\ell }}{k}_{0,\mathrm{\ell }}^h{M}_{\mathrm{\ell }}^k),$

$W_h^{k+1}-W_h^k=\mathrm{\Delta }t({\sum }_{j=1}^{J-1}(x_{j+1}^c-x_j^c)F_j^k-x_J^c{F}_J^k-{\sum }_j{x}_j^c{d}_j{M}_j^k+{\sum }_{\mathrm{\ell }}{k}_{1,\mathrm{\ell }}^h{M}_{\mathrm{\ell }}^k)\mathrm{.}$

Thus a division kernel preserves the discrete size moment only when its first-moment production cancels the parent-removal term; this identity is not automatic for an arbitrary kernel.

For Finite-age renewal model with absorbing outflow, define $K(a)=\beta (a)S(0,a)$ and

$G_A(t)={\int }_0^{A-t}\beta (x+t)S(x,t)n_0^A(x)dx,0\le t\le T\mathrm{.}$

There is a unique continuous nonnegative birth rate satisfying

$B_A(t)=G_A(t)+{\int }_0^{t}K(a)B_A(t-a)da,$

and the finite-age density is, for almost every $0<a<A$,

$n^A(t,a)=\{\begin{array}{ll}{n}_0^A(a-t)S(a-t,t),& t\le a<A,\\ B_A(t-a)S(0,a),& 0\le a<t\mathrm{.}\end{array}$

It is the unique mild solution, has $B_A(t)={\int }_0^A\beta (a)n^A(t,a)da$, and obeys

$\mathrm{\parallel }{n}^A(t,\cdot ){\mathrm{\parallel }}_1\le e^{{\beta }_{\mathrm{\infty }}t}\mathrm{\parallel }{n}_0^A{\mathrm{\parallel }}_1,{\beta }_{\mathrm{\infty }}=\mathrm{\parallel }\beta {\mathrm{\parallel }}_{\mathrm{\infty }}\mathrm{.}$

For $t\le T<A$, its exact cumulative number outflow is

$O_A(t)={\int }_{A-t}^A{n}_0^A(x)S(x,A-x)dx\mathrm{.}$

Input a uniform age partition $0=a_0<a_1<\cdots <a_J=A$ with $\mathrm{\Delta }a=a_j-a_{j-1}$, nonnegative cell masses $M_j$, a bounded nonnegative fertility function $\beta :[0,A]\to [0,{\beta }_{\mathrm{\infty }}]$, and one rule setting for each cell. A setting is either midpoint with an integer panel count $q_j\ge 1$, or capped adaptive with finite depth $D_j\in \mathbb{N}$, positive integer panel-call budget $P_j\ge 1$, and indicator tolerance ${\tau }_j>0$. If no settings are supplied, use midpoint with $q_j=1$ in every cell. In each cell $C_j=[a_{j-1},a_j]$, let $u_j=M_j\mathrm{/}\mathrm{\Delta }a$ and approximate ${\int }_{C_j}\beta (a)u_{j}da$ by its selected positive rule:

1. Midpoint mode: use $u_j{Q}_M(\beta ;a_{j-1},a_j,q_j)$ from Composite midpoint rule on a uniform grid and record status midpoint.
2. Capped adaptive mode: call Finite-depth adaptive Simpson with an indicator on $\beta$ over $C_j$ with inputs $(D_j,P_j,{\tau }_j)$, then multiply its estimate by $u_j$ and retain its status and indicator. The fertility evaluator returns a finite value at every requested node or reports failure.

If an evaluation fails, return failure and no birth rate. Otherwise return

$B_h={\sum }_{j=1}^J{u}_j{Q}_j(\beta ),$

the per-cell statuses, and the sum of any reported adaptive indicators. Both rules have nonnegative quadrature weights summing to $\mathrm{\Delta }a$, so $0\le B_h\le {\beta }_{\mathrm{\infty }}{\sum }_j{M}_j$. An adaptive accepted status means only that its indicators meet the requested budget; neither that indicator nor a cap-reached estimate bounds the true quadrature error. Any convergence theorem using this procedure must separately assume or prove a true quadrature-defect bound.

Let $\mu ,\beta :[0,\mathrm{\infty })\to [0,\mathrm{\infty })$ be bounded $C^1$ functions with bounded derivatives, let $n_0\in C^1([0,\mathrm{\infty }))\cap L^1(0,\mathrm{\infty })$ be nonnegative, and fix $T>0$, $\mathrm{\Delta }t>0$, and $A>T$. Put ${\widehat{n}}_0=n_0{\mathbf{1}}_{[0,A-T)}$ and ${\tau }_A={\int }_{A-T}^{\mathrm{\infty }}{n}_0(a)da$. Choose the initial positive cohort quadrature for ${\widehat{n}}_0(a)da$ with support in $[0,A-T)$ and bounded-Lipschitz defect ${ϵ}_0$ on test functions of norm at most $1$. For the procedure Age-cohort quadrature along exact characteristics, set $N=\lfloor T\mathrm{/}\mathrm{\Delta }t\rfloor$, consider its levels $0\le k\le N$, and write ${\beta }_{\mathrm{\infty }}:=\mathrm{\parallel }\beta {\mathrm{\parallel }}_{\mathrm{\infty }}$. Assume $\mathrm{\Delta }t\mathrm{\parallel }\mu {\mathrm{\parallel }}_{\mathrm{\infty }}\le 1$. Define

$\mathrm{\parallel }\eta {\mathrm{\parallel }}_{B{L}^{\ast }}:={\sup }_{\mathrm{\parallel }\varphi {\mathrm{\parallel }}_{\mathrm{\infty }}+\mathrm{Lip}(\varphi )\le 1}\mid \int \varphi d\eta \mid$

for a finite signed measure $\eta$. Then, for fixed $T$, there is $C_T$ independent of $\mathrm{\Delta }t$, initial cohort count and cutoff such that

${\max }_{k\mathrm{\Delta }t\le T}\mathrm{\parallel }{\nu }_h^k-n(t_k,a)da{\mathrm{\parallel }}_{B{L}^{\ast }}\le C_T({ϵ}_0+\mathrm{\Delta }t+{\tau }_A)\mathrm{.}$