Numerical investigation of the dynamics for a normalized time-fractional diffusion equation

1.   Introduction

In recent decades, time-fractional diffusion equations have received significant attention across various fields. Time-fractional derivatives extend classical models by incorporating memory effects, which better capture real-world processes where the rate of change is not constant ^[1]. They are particularly useful for modeling phenomena with irregular or non-standard diffusion behaviors, such as in heterogeneous media ^[2], complex fluids ^[3], or biological systems ^[4]. These time-fractional derivatives allow for a more accurate representation of systems with long-range dependencies, thus providing deeper insights and improved predictions in various scientific and engineering fields. Some specific applications of time-fractional derivatives are as follows. The time-fractional Burgers equation, applying a fractional differential method to the classical Burgers equation ^[5,6], was used to model a range of physical processes such as turbulence, shock waves, and traffic flow, incorporating memory effects over time. In biology, the time-fractional method was applied to the reaction-diffusion equation, specifically Fisher's equation, which models species propagation ^[7,8]. In the field of finance, European and double-barrier options were evaluated using the time-fractional Black–Scholes equation ^[9,10,11]. In particular, for double barrier options, the closer the underlying asset price is to the lower barrier, the more the Black–Scholes model tends to overestimate the option's value. Furthermore, a smaller $\alpha$ exacerbates this price bias. The reaction-diffusion equation, known as the Allen–Cahn (AC) equation, finds applications in many fields such as physics, materials science, and biology. By applying the time-fractional method to the AC equation, these approaches accurately model complex dynamical systems by incorporating memory effects and history dependence associated with time ^[12,13].

There have been many numerical methods for the time-fractional diffusion equations ^[14,15,16]. For simplicity of exposition, let us consider the following one-dimensional conventional time-fractional diffusion equation:

where $u(x, t)$ is the concentration at $x$ and $t$, and $u_0(x)$ is the initial condition,

where $\Gamma(z) = \int_0^\infty \tau^{z-1}e^{-\tau}d \tau$ is the gamma function. It is noted that when $\alpha = 1$, Eq (1) reduces to the conventional diffusion equation ^[17]. Let us define a weight function $w_\alpha^t(s)$ as follows:

Then, Eq (4) can be rewritten as

For different values of $\alpha = 0.1, 0.5,$ and $0.9$, the weight functions $w_\alpha^t(s)$ at $t = 1$ are illustrated in Figure 1(a). It can be seen that the weight functions $w_\alpha^t(s)$ remain flat for small $\alpha$ and show a sharp transition near time $t$ for large $\alpha$. Figure 1(b) shows the weight functions for different time values of $t = 0.5, \; 1$, and $2$ with $\alpha = 0.5$. As $t$ increases, the functions simply translate to the right direction.

We note that $\int_0^t w_\alpha^t(s)ds$ approaches infinity as $t$ increases, for any values of $0 < \alpha < 1$. That is,

which approaches infinity as $t$ increases, for any values of $0 < \alpha < 1$. From Eq (7), we can see there are scaling differences associated with values of $\alpha$ when comparing the effects of $\alpha$ on the dynamics of the time-fractional diffusion equations because $W_\alpha(t)$ depends on both $\alpha$ and time $t$. This is physically sound and can be inferred intuitively from previous studies ^[18,19] related to wave propagation and diffusion problems.

To resolve these scaling differences associated with $\alpha$ values, we propose a normalized time-fractional diffusion equation and conduct numerical investigations of the dynamics of the proposed equation. In this study, we propose the following normalized time-fractional diffusion equation:

where

where $(1-\beta)/t^{1-\beta}$ is the normalizing factor, which makes the right-hand side term in Eq (11) $\partial u / \partial x$ when $\partial u / \partial x$ is constant. That is

Let us define a weight function $w_\beta^t(s)$ as follows:

Then, from Eq (7), we have

which is independent of the fractional order $\beta$ and time $t$, unlike that of the conventional time-fractional derivative, $W_\alpha(t) = t^{1-\alpha} / \Gamma(2-\alpha)$. To the authors' knowledge, this is the first time that the normalized time-fractional diffusion equation is proposed, where the total integration of the weight function is always one for all time-fractional orders and times. For $\beta$ values of $0.1$, $0.5$, and $0.9$, the weight functions $w_\beta^t(s)$ at $t = 1$ are as shown in Figure 2(a). We can observe that the weight functions $w_\beta^t(s)$ are flat when $\beta$ is small and exhibit a sharp transition near time $t$ when $\beta$ is large. Figure 2(b) shows the weight functions for different times $t = 0.5, \; 1$, and $2$ with $\beta = 0.5$.

Figure 3 shows the temporal evolutions of $W_\alpha(t) = \int_0^t w_\alpha^t(s)ds = t^{1-\alpha}/\Gamma(2-\alpha)$ for $\alpha = 0.1, \; 0.5$, and $0.9$. $W_\alpha(t)$ is an increasing function with respect to time $t$ for a fixed fractional order $\alpha$. At early times, $W_\alpha(t)$ increases with respect to the fractional order $\alpha$, whereas at later times, $W_\alpha(t)$ decreases with respect to the fractional order $\alpha$ for a fixed time $t$. However, $W_\beta(t)$ is independent of the fractional order $\beta$ and time $t$.

The contents of this paper are as follows. In Section 2, numerical solution algorithms for the conventional and normalized time-fractional diffusion equations are presented. In Section 3, numerical experiments are provided. Finally, Section 4 presents conclusion and potential progress for future study.

2.   Numerical solutions

2.1.   Conventional time-fractional diffusion equation

Let $\Omega = (L_x, \; R_x)$ be the computational domain, which is discretized as follows: $\Omega_h = \{ x_i | x_i = L_x+(i-1)h, \; i = 1, \ldots, N_x \}$, where $h = (R_x-L_x)/(N_x-1)$ for some positive integer $N_x$, see Figure 4.

Let $u_i^n = u \left(x_i, \; t_n \right)$ and $t_n = (n-1) \Delta t$, where $\Delta t$ is the time step. Equation (4) can be approximated by the following numerical quadrature formula:

where we have used the identity $(1-\alpha)\Gamma(1-\alpha) = \Gamma(2-\alpha)$ and approximated $\partial u (x_i, s) / \partial s$ over the interval $[t_p, \; t_{p+1}]$ using the finite difference $(u_i^{p+1}-u_i^p)/\Delta t$ in Eq (15). Therefore, we have the following finite difference discretization of Eq (1) using Eq (16):

where

Here, we use the zero Dirichlet boundary condition: $u_0^{n+1} = 0$ and $u_{N_x}^{n+1} = 0.$ Then, Eq (17) can be rewritten as follows:

which can be rearranged as follows:

Note that an implicit temporal discretization is used for stability. Although fully explicit schemes are generally sufficient in terms of stability and efficiency for second-order partial differential equations, as indicated in ^[20,21], an implicit scheme is necessary due to the presence of the source term, $-\sum_{p = 1}^{n-1} w_p^n (u_i^{p+1}-u_i^p)/\Delta t$, in Eq (20).

Equation (20) is a tridiagonal system with a zero Dirichlet boundary condition and we can use the Thomas algorithm ^[22] to efficiently solve this system. Thus, the solution vector $\mathbf{u}^{n+1}$ can be found by solving the tridiagonal system using the Thomas algorithm:

where $A$ is a tridiagonal matrix with zero Dirichlet at $i = 1$ and $i = N_x$. The detailed procedure is provided below.

where $F = \sum_{p = 1}^{n-1}w_p^n(u_i^{p+1}-u_i^p)/\Delta t$.

2.2.   Normalized time-fractional diffusion equation

Equation (11) can be approximated by the following numerical quadrature formula:

Therefore, we have the following finite difference discretization of Eq (8) using Eq (21):

where

We note that the weight parameter $w_p^n$ satisfies the following condition for any value of $n$:

Then, Eq (22) can be rewritten as follows:

Equation (25) is a tridiagonal system with a zero Dirichlet boundary condition and we can use the Thomas algorithm to efficiently solve this system. Thus, the solution vector $\mathbf{u}^{n+1}$ can be found by solving the tridiagonal system using the Thomas algorithm:

where $A$ is a tridiagonal matrix with zero Dirichlet at $i = 1$ and $i = N_x$.

3.   Numerical experiments

In this section, we present several numerical experiments in a finite domain $\Omega \times (0, T)$, where $\Omega = (0, 1)$ and $T$ is a final time, to investigate the effects of $\alpha$ and $\beta$ on the evolution dynamics of the conventional and normalized time-fractional diffusion equations. As the first numerical test, we consider the following low-frequency initial condition:

Figure 5(a) and (b) show $u(x, T)$ for different values of $\alpha = 0.1, 0.5$, and $0.9$; and $\beta = 0.1, 0.5$, and $0.9$, respectively. Here, $T = 0.001$ is used. As shown in Figure 5(a), for the conventional time-fractional diffusion equation, we observe that the temporal evolution is faster when the value of $\alpha$ is smaller. This phenomenon is attributed to the different weight function associated with different $\alpha$ values, as shown in Figure 1(a). When $\alpha$ is smaller, the total sum of the weight function at early times is smaller as shown in Figure 3, which results in an effectively larger diffusion process and faster temporal evolution. This can be understood as follows: if we divide both sides of Eq (1) by the smaller weight value, we obtain an effectively larger diffusion coefficient, which leads to an increased diffusion process. However, in the case of the normalized time-fractional diffusion equation, there is little variation with respect to different values of $\beta$ compared to the conventional time-fractional diffusion equation, as seen in Figure 5(b).

In the second numerical experiment, let us consider the following high-frequency initial condition:

Figures 6(a) and (b) show $u(x, T)$ for different values of $\alpha = 0.1, 0.5$, and $0.9$; and $\beta = 0.1, 0.5$, and $0.9$, respectively. We observe that the temporal evolutions are faster than those of the low-frequency initial condition for both the conventional and normalized time-fractional diffusion equations. As illustrated in Figure 6(a), the temporal evolution for the conventional time-fractional diffusion equation accelerates more than that of the low-frequency initial condition as the value of $\alpha$ decreases. A smaller $\alpha$ results in a smaller total sum of the weight function, which leads to an effectively larger diffusion process and consequently faster temporal evolution. However, for the normalized time-fractional diffusion equation, there is a little variation concerning different $\beta$ values, unlike the conventional time-fractional diffusion equation, as shown in Figure 6(b).

For the third numerical test, let us consider the following random initial condition:

where $rand(x)$ is a random number between $-1$ and $1$. Figures 7(a) and (b) show $u(x, T)$ for different values of $\alpha = 0.1, 0.5$, and $0.9$; and $\beta = 0.1, 0.5$, and $0.9$, respectively. The random initial condition can be considered as a combination of multiple frequency modes. As expected from the previous computational results, we observe large variation with respect to $\alpha$ values and small variation with respect to $\beta$ values.

In the fourth numerical test, we consider the following initial condition:

This initial condition contains multiple modes with gradually changing profiles. Figure 8(a) shows the numerical solution $u(x, T)$ for a relatively short time, $T = 0.0005$, and Figure 8(b) and 8(c) display $u(x, T)$ for a relatively long time, $T = 0.05$. Figures 8(a) and 8(b) present results for different values of $\alpha = 0.1, 0.5$, and $0.9$, and Fig. 8(c) shows results for different values of $\beta = 0.1, 0.5$, and $0.9$. Because the initial condition is a combination of multiple frequency modes, we observe significant evolutions at high frequencies and minor changes at low frequencies.

In the final numerical test, we examine the evolution process of the conventional time-fractional method and the normalized time-fractional method for different fractional orders $\alpha = 0.1, 0.9$ and $\beta = 0.1, 0.9$. We consider the following initial condition:

In this test, we used a final time of $T = 5$. Figures 9(a) and (b) represent the evolution of the conventional time-fractional diffusion equation for the different $\alpha = 0.1$ and $0.9$, respectively. When $\alpha$ is small, the temporal evolution is fast at early times but suddenly slows down at later times. This phenomenon, as shown in Figure 3, can be attributed to the small weight function values at early times and the large weight function values at later times. When $\alpha$ is large, the opposite behavior can be observed. Next, Figures 9(c) and (d) illustrate the evolution of the normalized time-fractional diffusion equation for the different values of $\beta = 0.1$ and $0.9$, respectively. In this case, the temporal evolution depends on the $\beta$ values but it is relatively independent of time.

4.   Conclusions

In summary, the evolution dynamics of the conventional and normalized time-fractional diffusion equations are different. The total integration of the weight function for the conventional time-fractional diffusion equation is not normalized, continues to increase to infinity as time progresses; it depends on the values of $\alpha$ and time $t$. Due to these dependencies, it is difficult to study the effects of $\alpha$ on the evolutionary dynamics. As an alternative, we proposed a normalized time-fractional diffusion equation that ensures that the total integration of the weight function remains one, regardless of the values of $\beta$ and time $t$. The proposed weight functions provide an intuitive understanding of the fractional time derivative, which represents a weighted average of the temporal history of the time derivative. Numerical experiments indicate that the normalized time-fractional diffusion equation shows slower dynamics for simple initial conditions as the fractional order parameter $\beta$ decreases. The diffusion equation is a fundamental model for various phenomena such as heat conduction, fluid flow, chemical dispersion, cellular diffusion processes, and the evolution of option prices and market risks in finance. As future work, we will apply the proposed model to the Black–Scholes equation ^[23,24], the AC equation ^[25,26,27], and the nonlocal Cahn–Hilliard equation ^[28], and investigate the development of a fast algorithm to reduce memory requirements ^[29].

Author contributions

Chaeyoung Lee: Investigation, Writing–original draft, Validation, Writing–Review & Editing, Software; Yunjae Nam: Formal analysis, Visualization, Writing–original draft, Writing–Review & Editing; Minjoon Bang: Writing–original draft, Writing–Review & Editing, Software, Visualization; Seokjun Ham: Validation, Writing–original draft, Writing–Review & Editing, Visualization; Junseok Kim: Conceptualization, Investigation, Formal analysis, Supervision, Writing–original draft, Writing–Review & Editing. All authors have read and agreed to the published version of the manuscript.

Acknowledgments

The first author (C. Lee) was supported by the National Research Foundation of Korea (NRF) grant funded by the Korea government (MSIT) (No. 2022R1C1C2003896). The corresponding author (J.S. Kim) was supported by the Brain Korea 21 FOUR from the Ministry of Education of Korea. The authors are grateful to the reviewers for their insightful and constructive feedback, which has significantly enhanced the quality of this paper.

Conflict of interest

The authors declare that there is no conflict of interest.

Appendix

The following listings 1 and 2 are Python codes for the conventional and normalized time-fractional diffusion equations, respectively.

Listing 1. Python code for a conventional time-fractional heat equation

Listing 2. Python code for a normalized time-fractional heat equation