Many problems in nature have multiscale characteristics. For example, material weathering, material fracture due to atomistic defects such as impurities, chemical reaction of many substances mixed with fast and slow reaction rates, or some biological problems with a slow growth. In these processes large-scale or long-term phenomena are often influenced by important small-scale or short-term effects, and the variation of long-term properties will also effect the short-term response. A method that fully resolve short-term effects may not be effective for catching a long-term property. Therefore we must effectively deal with these two-way coupled problems.
In the past two decades, people are interested in developing efficient multiscale methods for different temporal/spatial multiscale problems. For spatial multiscale problem, people mainly develop computable models such as atomistic-to-continuum coupling (AtC) models methods to simulate material behaviors involved with the atomistic scale, where the simulation of materials at the full atomistic scale is not possible due to its computatonal cost (See, e.g. [Tadmor1996, Lin2003, Weinan2007, Shenoy1999, Weinan2006, Shapeev2011, Ortner2012, Abdulle2012, Luskin2013] and most recently [Wang2018] for AtC of both simple and complex lattice structures). For temporal multiscale problem, one of the most prominent class of techniques is the Heterogeneous Multiscale Method (HMM) [Engquist2005, E2003, Weinan2003, Abdulle2012(2)] that aims at efficient decoupling of a macro-scale and micro-scale problem, where the latter one enters the macro-scale problem in terms of temporal averages. In the context of ordinary differential equations, complexity plays a smaller role. The authors of [Ariel2009, Ariel2009(2)] numerically realize analytical averaging techniques for problems that are also characterized by resonance effects. Frei and Richter addressed the Navier-Stokes flow problem with multiscale characteristics in time [Frei2020, Richter2021(1)]. There are other applications and developments related to temporal multiscale problems (see, e.g. [E2011, Richter2021(2), Trucu2013]). To our knowledge, none of those has considered fractional differential equations.
In this paper, we consider a temporal multiscale problem with fractional differentiations. The HMM [Weinan2003, E2003, E2011, Abdulle2012, Engquist2005] is a general framework for designing multiscale algorithms, and in particular, temporal multiscale problems. Here, we consider a nonlinear coupled system of fractional ODEs, which is simplified from describing the problem of atherosclerosis [Frei2020, Yang2016]. In [Frei2020], a multiscale approach has been rigirously developed for a linear form of ODEs as a simplification of the evolution of plaques in blood vessels, where a periodic force is applied to the flow due to the periodic heart pulse nature. Considering that the plaque accumulation or evolution may have memory effects in time and that a blood flow may be highly nonlinear, we thus aim to develop the method for a fractional order multiscale nonlinear ODE system with a periodic applied force. This development will biuld up a fundation to deal with many types of fractional order multiscale ODEs/PDEs in the future.
The study of fractional calculus has attracted the attention of many scholars in the last 20 years, and it has been used in the modeling of problems in science, such as fluid mechanics [Song2016], biomathematics [Kosztolowicz2020], anomalous diffusion [Zhaoyang2021, Sandev2019, Kilbas2006], etc. On the contrary to the integer one, the fractional derivative has memory effects and non-local characteristics, which naturally appears in various problems in real world and leads to an increase in computational cost [Gao2011]. In many cases, we are interested in the behavior of solution in a long period of if time. For a temporal multiscale problem with macroscopic and microscopic features in time, we simply increase the time step size for a long time computation, the accuracy can not be guaranteed and the temporal microscopic behaviour is not visible under a large time step. Our basic idea is to construct a local periodic problem to approximate the microscopic bahaviour of the original problem. The advantage is that the fast variable with microscopic feature in each macro-step can be solved effectively and independently. We prove that the constructed local periodic problem can effectively approximate the original problem by analyzing the error of the coupled micro-macro multiscale method. Numerical tests are implemented to show the correctness of our algorithm. By using numerical schemes proposed, we find a significant acceleration effect compared with the fully resolved simulation at least when the time scale separation is larger.
The rest of this manuscript is organized as follows. In Section 2, we introduce a temporal multiscale problem with fractional derivatives considered and prove several auxiliary lemmas to obtain effective approximate periodic equations. In section 3, we conduct the error analysis at the continuous level. In Section 4, we present the temporal discretisation of the multiscale scheme for solving the proposed problem. Numerical tests are presented in Section 5 to illustrate the accuracy and computational efficiency of our method. Some conclusions and remarks are provided in Section 6.
2 Model problem and preliminaries
We consider a nonlinear coupled system of fractional ordinary differential equations.
Problem 1. On I=[0.T], let , ,
Here, and represent the fast variable and the slow variable, respectively. is the Caputo fractional derivative of order [Podlubny1999]. is a small parameter that controls the change of . The force term is a local periodic function with period 1 without loss of generality, that is, , is the integer ceiling function.
In addition to describing the chemical kinetics system [Singh2017], when , problem 1 is a simplified atherosclerosis model with a certain type of memory effect included in the vessel wall or plaque growth. Macrophages in the vessel wall take up low density lipoproteins (LDL), which carry cholesterol and triglycerides to the tissues, and are finally transformed into foam cells, which are engorged with lipids [Hahn2009]. The flow dynamics in the vessels are on the scale of milliseconds to seconds. However, the growth of plaque is usually described in terms of months, which is a slow process. For the formulation of the plaque growth equation, we consider the existence of a variety of cells in the vessel wall such as endothelial cells and vascular smooth muscle cells [Yang2016, Libby2002], and thus the diffusion of macrophages and foam cells is most likely to be anomalous. It may thus be more suitable to use a fractional operator representation.
Remark 1.1. Problem 1 is a simplification of the problem of atherosclerosis described by a coupled PDEs. The goal of this paper is to present an idea to deal with the atherosclerosis problem based on this nonlinear fractional ODE model. The idea may be applied to the PDE problem, though analysis and implementation would be much more complicated and considered in future (Preliminary results have been obtained).
To apply the multiscale method and ensure the existence of the solution, it is necessary to have some assumptions about and .
Assumption 2. Let , . and its derivative are bounded
and Lipschitz continuous with respect to slow and fast variables
Assumption 3. Let , . is bounded
and is differentiable with a bounded derivative
Remark 3.1. From Eq. (5), we can obtain the Lipschitz condition with respect to slow and fast variables
Remark 3.2. For the linear form given in [Frei2020], i.e., , assumption 3 is reduced to . (5) ensures that problem 1 has a local periodicity. Furthermore, we would like to point out that assumption 3 is reasonable if this model is used to describe the simplified atherosclerosis problem [Frei2020]. Its derivative is bounded as long as the growth of plaque does not change the topology of the computing domain.
We introduce a new variable
where is the Riemann-Liouville fractional integral operator and defined by
By inserting , we have
LEMMA 4. Let , , and let assumption 2 be satisfied. It holds that
with a constant .
Applying the fractional integral operator on both sides of the second equation of (1), we obtain
Similarly, we can also get an estimate of the new variable
By Lipschitz continuity of , we have
We can thus approximate the averaged evolution equation for by
The discretization of Eq. (14) in a macro-time step is not effective because it is involved with the dynamic evolution of to on the fast scale. We need deal with and approximate the fast variable solution in a more effective way. We thus construct an auxiliary problem on a period to approximate the original problem 1. Next, we give two lemmas to derive an equivalent local periodic problem, which may provide a fast scale solution without initial value.
LEMMA 5. Let , and . For the first equation of (1), we have the following estimate
LEMMA 6. Let , and . For a fixed , if or , the periodic condition to the equation given below holds
Let , we have
By the mean value theorem, we obtain the following
with the solution
Hence, we obtain . ∎
Based on the above lemma, for this temporal multiscale problem, we construct the following auxiliary problem.
Problem 7. On I=[0.T], let , ,
Remark 7.1. If the assumption of is periodic, i.e., , the periodic condition of Eq. (24) should be [Frei2020].
Remark 7.2. Lemma 6 implies that , will show local periodic properties. Since we focus on long-term simulation, the initial non-periodic part has little effect on it, which is confirmed in numerical experiments.
3 Error analysis of the fractional multiscale problem
In this section, we will prove that the difference between variable in the auxiliary problem 24 and variable in original problem 1 is small. The error estimate between and is also obtained. Therefore, the solution of problem 24 can be used as an effective approximation of the solution of problem 1.
LEMMA 8. Let be fixed and . For the local periodic problem
The following lemma will give the error estimate between the local periodic solution and the fast variable . For a fixed , we consider the family of periodic solutions
with a constant .
For , we have
the function fulfils the following equation
Let , it holds
For , we use the fractional differential operator on both sides of the second equation of (1), we have
By Gronwall inequality, we can easily get the following estimates
Lemma 9 proves that the initial value problem 1 can be transformed into a local periodic problem 24 with fixed slow variable . Next we will prove the last lemma to show the final error estimate between the solution of the original problem 1 and the solution of the multiscale approximation problem 24.
Let , we have
Lipschitz continuity of in assumption 2 gives
By using lemma 4 and lemma 8, we get the following estimates
where . Applying Gronwall’s inequality technique it is easy to get
which satisfies for . Finally, we insert and use lemma 9 to obtain
We have completed the theoretical analysis of an approximate multiscale problem, which is the fundation for the numerical calculation of fractional temporal multiscale problems.
4 Numerical algorithms
For Eq. (1), the L1 approximation [Gao2012] for the Caputo fractional derivative of order is given by
where , .
Then, the governing equation (1) is discretized as follows
Here, for the calculation of the fast variable , we will use the first-order explicit/implicit Euler scheme for discretization. The calculation process is shown as follows
For the calculation of Eq. (24), we first consider the identification of periodic solutions in a period. For a given initial value , iterate to , and then update the initial value to and enter the loop until the given error limit is reached. This scheme is more general than the averaging acceleration scheme applied in [Richter2021(1)]. It is suitable for calculating fractional periodic problems but the calculation time is longer. Below we give the theoretical basis for this shooting method.
Let and be the trial solution and the periodic solution respectively. is the initial value for trial, for fixed , it is not difficult to obtain the following estimates
We can see that the longer the time, the better the identification efficiency of periodic solutions.
From lemma 6, we can see that the decay of the nonstationary solution to periodic solution may be slow, which depends on . If necessary, we can ignore the case of nonstationary solutions and focus on the acceleration of periodic problems.
Based on Algorithm 2, we give the following multiscale method.
5 Numerical experiments
In this section, we carry out several numerical experiments to test the correctness and effectiveness of the proposed numerical algorithm. Specifically, we construct exact solutions to verify the correctness of fully resolved simulation and multiscale method, and then compare the calculation time and accuracy.
5.1 Local test
Example 1. The coupled differential equations are as follow
The approximate local periodic equation is
Here, , .
This problem has the following exact solution
In the simulation, we set the time step , and . The first-order explicit Euler scheme is used to discretize the fast variable equation. It can be seen from Fig. 1 that the numerical computational results of Algorithm 1 and Algorithm 2 are in good agreement with the exact solution. In Algorithm 2, we take . We would like to point out that Algorithm 2 is also suitable for solving other types of boundary value problems (see, e.g., [Ahmad2011, Ahmad2021]).
5.2 Global test
We further consider a long term numerical simulation to test the effectiveness of the multiscale method.
Example 2. The linear coupled system of fractional differential equations are given by
where , .
The corresponding approximate equations are as follow
We has an exact solution for (53)
For the discretization of fast variable, we use the first-order implicit Euler scheme. For Algorithm 1, we set the time step , , and for long-term numerical simulation. For Algorithm 3, we set the micro-scale time step and the macro-scale time step to carry out the simulation. Fig. 2 presents the comparison of the numerical results with analytical solution for a long time. Fig. 3 shows the calculation results of the multiscale method and the comparison with the exact solution. It is not difficult to see that the multiscale method has good accuracy.
In order to better test the effectiveness of the multiscale method, we set the micro-scale time step as . Quantitative error measurements and computational time are shown in Table 1. We observe that the error decreases as the macro-scale time step decreases, and the computational cost increases accordingly. The convergence rates show that our numerical method is first order accurate. It should be noted that due to the memory characteristic of the fractional operator, the macro-scale time step we choose should be much smaller compared to the computation of the integer time multiscale problem.
Next, we test the impact of the micro-scale time step on the performance of the multiscale method. As shown in Table 2, the micro-scale time step has little effect on the performance of the proposed method. The reason may be that the equation for fast variable is of integer order.
For a fully resolved simulation, the error is 0.0223, but the time is up to 19 hours. If we need to calculate for a longer time, as shown in the second equation of (48), the calculation of fractional operator will greatly increase the calculation cost, and the multiscale method presents a huge advantage. According to this example, the calculation time is saved by more than 1000 times, and the maximum error does not exceed 4% (). When the minimum macro-scale time step is , we find that the acceleration is nearly 300 times and the maximum error is only 2%. It is worth pointing out that
the discretization scheme of Eq. (53) is of low-order accuracy, which may cause large error results.
Example 3. Next we test two Riccati equations coupled with each other
where . We set and the scale separation parameter .
The approximate equations are as follow