I'm trying to model a diffusion process characterized by subsequent logistic diffusion processes. To give a better idea of what I'm trying to model, you can take a look at the attached figure. In such a figure, I represented a double-logistic diffusion process, where the function $f(t)$ remains stable, then has a steep increase, then remains stable, then a steep increase, and finally remains stable.
A single logistic diffusion process is usually governed by: $$ f(t) = \frac{L}{ 1 + \exp(-\kappa(t-t_{0}))},$$
where $L$ is the curve's maximal value, $\kappa$ is the steepness of the curve, and $t_{0}$ is the t-value of the sigmoid's midpoint.
How can I model a process where logistic diffusion processes repeat for $k$ times? Any ideas or references will be much appreciated!
A possible way of dealing with the problem is to view it as a dynamical system and to use phase portraits. Let me rename $f(t) = x(t)$ and rescale things such that $x$ obeys the logistic equation $dx/dt \equiv \dot{x} = - kx(1-x)$. The phase portrait is then a plot of $\dot{x}$ versus $x$. For $k > 0$ this is an inverted parabola:
As one can see, the zero at $x=1$ corresponds to a stable equilibrium, often called the carrying capacity in a population dynamics context. The solution of the logistic equation is the sigmoid given in the question,
$$ x(t) = \frac{1}{1 + \exp(-k(t - t_0))} \; . \quad (1) $$
From this one can see that it takes an infinite amount of time to reach equilibrium, $x=1$. Alternatively, this can be inferred from the integral (or by solving for $t = t(x)$),
$$ t - t_0 = \frac{1}{k} \int_{1/2}^x \frac{dy}{y(1-y)} = \frac{1}{k} \, \ln \frac{x}{1-x} \; , \quad (2) $$
which diverges for $x \to 1$.
To obtain the behaviour shown in the question (a smoothed staircase function), the phase portrait must be similar to that of a rotating plane pendulum. Qualitatively, the phase portrait needs to look as follows:
This is similar to quasi-periodically replicating parabolae, but in such a way that there are no intermediate zeros (attractors). The closer the minima of $\dot{x}(x)$ are to zero, the flatter will be the plateaux in the graph of $x = x(t)$. The physical analogue is a plane pendulum starting out from an inverted position (unstable equilibrium) with just enough energy to perform a few revolutions until the energy drops adiabatically, whereupon a stable equilibrium is reached.
This can indeed be modelled with a "pendulum equation" of the type
$$ \dot{x} = k \omega_0 \sqrt{1 - (1/k^2) \sin^2 x} \; , \quad (3) $$
which describes the 'wiggly' behaviour of the plot above. The monotonous initial and final pieces of the graph have been attached by replacing the square root by its quadratic Taylor expansions at $x \le - 10$ and $x \ge 10$.
[For a pendulum of length $\ell$, $x(t)$ is the total angle covered in time $t$, $\omega_0^2 = g/\ell$ is the unperturbed (harmonic) angular frequency and $k$ the elliptic modulus depending on the amplitude of the oscillation. The energy of the pendulum is then $2mg\ell k^2$ with $mg$ the usual gravitational force.] It is important to have $k = 1 + \epsilon$, $0 < \epsilon \ll 1$, [hence $E > 2mg\ell$] to avoid zeros in $\dot{x}$ corresponding to intermediate equilibria. After all, the dynamics should not get stuck before the final plateau in $x(t)$ is reached.
Unfortunately, as is well known from the discussion of the plane pendulum, the differential equation $(3)$ yields an elliptic integral, namely
$$ k \omega_0 t = \int_0^x \frac{dy}{\sqrt{1 - (1/k^2) \sin^2 y}} \; , \quad (4) $$
which replaces (2). Solving for $x$ one finds
$$ x(t) = \text{am} (k\omega_0 t, 1/k) \; , \quad (5) $$
with $\text{am}$ denoting the Jacobi amplitude. Maple knows this as JacobiAM, and hence one can plot it straightforwardly (say for k = 1.01):
Adding the monotonous bits of the second plot to (3) by hand should yield horizontal asymptotes.
I have tried to model the generic behaviour displayed in Plot 2 by various trigonometric functions. The integrals can often be done, but then one has to worry about correctly continuing inverse trig (arc) functions beyond their fundamental domains, which I found somewhat tedious.
The pendulum analogue works nicely and seems appropriate for a physics site. I reiterate, though, that all one needs in principle is the generic phase portrait shown in Plot 2.