On transient waves in linear viscoelasticity

TL;DR

This paper reviews the wave-front expansion method, including the Buchen-Mainardi algorithm, and discusses its applications to fundamental models of linear viscoelasticity, both ordinary and fractional.

Contribution

It provides a comprehensive overview of the wave-front expansion technique and its application to viscoelastic models, highlighting recent developments and specific use cases.

Findings

Effective application of wave-front expansion to viscoelastic models

Extension of the Buchen-Mainardi algorithm to fractional models

Detailed discussion of practical applications in viscoelasticity

Abstract

The aim of this paper is to present a comprehensive review of method of the wave-front expansion, also known in the literature as the Buchen-Mainardi algorithm. In particular, many applications of this technique to the fundamental models of both ordinary and fractional linear viscoelasticity are thoroughly presented and discussed.

Figures

Equations

$$c=\frac{1}{\sqrt{\rho J_0}},$$

$$\Psi (t)=\frac{1}{J_0}\frac{dJ(t)}{dt},\text{mbox}ast>0,$$

$$\frac{\partial^{2}r}{\partial x^{2}}-\frac{1}{c^{2}}\,\Big{[}1+\Psi(t)\ast\,\Big{]}\frac{\partial^{2}r}{\partial t^{2}}=0$$

$$\left[\frac{{\partial }^2}{\partial x^2}-{\mu }^2(s)\right]\tilde{r}(s,x)=0,$$

$$\mu (s)\equiv s{\left[\rho s\tilde{J}(s)\right]}^{1/2}.$$

$$\tilde{r}(s,x)={\tilde{r}}_0(s)\exp \left[-x\mu (s)\right].$$

$$\mu (s)\stackrel{s\to \infty }{\sim }\sum _{k=0}^{\infty }{b}_k{s}^{1-{\beta }_k},0={\beta }_0<{\beta }_1<\cdots$$

$${\mu }_{+}(s)\equiv \sum _{k=0}^m{b}_k{s}^{1-{\beta }_k},{\beta }_k\le 1,k=0,1,\dots ,m.$$

$${\mu }_{-}(s)\equiv \sum _{k=0}^{\infty }{b}_k{s}^{1-{\beta }_k}-{\mu }_{+}(s).$$

$$\tilde{r}(s,x)\stackrel{s\to \infty }{\sim }{\tilde{r}}_0(s)\exp \left[-x{\mu }_{+}(s)\right]\tilde{R}(s,x),$$

$$\tilde{R}(s,x)\equiv \exp \left[-x{\mu }_{-}(s)\right].$$

$$J(t)\sim J_0+O(t^{\alpha }),\text{mbox}ast\to 0^{+}(\iff s\to \infty ),$$

$$\{\begin{array}{rl}& J_0=0\Rightarrow b_0=0,{\beta }_1=\alpha /2;\\ & J_0>0\Rightarrow b_0=1/c,{\beta }_1=\alpha .\end{array}$$

$${\mu }_{-}(s)=\sum _{k=m+1}^{\infty }{b}_k{s}^{1-{\beta }_k},\text{mbox}with{\beta }_k>1,k=m+1,m+2,\dots$$

$$\tilde{R}(s,x)\stackrel{s\to \infty }{\sim }\sum _{k=0}^{\infty }{v}_k(x)s^{-{\lambda }_k},0\le {\lambda }_0<{\lambda }_1<\cdots ,$$

$$\mathcal{O}\tilde{R}(s,x)=\left\{\frac{d^2}{d{x}^2}-2{\mu }_{+}(s)\frac{d}{dx}-[{\mu }^2(s)-{\mu }_{+}^2(s)]\right\}\tilde{R}(s,x)=0,$$

$$\tilde{R}(s,0)=1.$$

$$\mathbb{L}\stackrel{\epsilon \to 0}{\sim }\sum _{i=0}^{\infty }{\epsilon }^{{\nu }_i}{\mathbb{L}}_i,0={\nu }_0<{\nu }_1<\cdots ,$$

$$\mathbb{L}f=0,$$

$$f\stackrel{\epsilon \to 0}{\sim }\sum _{k=0}^{\infty }{\epsilon }^{{\lambda }_k}{v}_k,{\lambda }_0<{\lambda }_1<\cdots .$$

$${\mathbb{L}}_0{v}_0=0,$$

$${\mathbb{L}}_0{v}_k=-\sum _{ij}{\mathbb{L}}_i{v}_j,k=1,2,3\dots ,$$

$${\mathbb{L}}_0{v}_0=0,$$

$${\mathbb{L}}_0{v}_k=-\sum _i{\mathbb{L}}_i{v}_{j(i,k)},k=1,2,3\dots ,$$

$${\mathbb{L}}_0=\frac{\partial }{\partial x},$$

$${\mathbb{L}}_i=p_i\frac{{\partial }^2}{\partial x^2}+q_i\frac{\partial }{\partial x}+r_i,i=1,2,3\dots ,$$

$${\lambda }_0=0,v_k(0)={\delta }_{k0},$$

$${\lambda }_k=\sum _{i=1}^N{m}_i{\nu }_i,$$

$$\frac{\partial v_0}{\partial x}=0,\frac{\partial v_k}{\partial x}=-\sum _i\left(p_i\frac{{\partial }^2}{\partial x^2}+q_i\frac{\partial }{\partial x}+r_i\right)v_{j(i,k)},k=1,2,3\dots ,$$

$$v_k(x)=\sum _{\ell =0}^k{A}_{k,\ell }\frac{x^{\ell }}{\ell !}.$$

$$\{\begin{array}{l}{A}_{k,\ell }={\delta }_{k\ell }\ell =0,\\ A_{k,\ell }=-\sum _i\left[\sum _J\left(p_i{A}_{J,\ell +1}+q_i{A}_{J,\ell }+r_i{A}_{J,\ell -1}\right){\delta }_{J,j(i,k)}\right]1\le \ell \le k,\\ A_{k,\ell }=0\ell >k,\end{array}$$

$$\tilde{r}(s,x)$$

Full text

Paper published in Wave Motion 74C (2017) pp. 191–212, DOI: 10.1016/j.wavemoti.2017.07.008.

1. 1 Introduction

2. 2 Buchen-Mainardi wave-front expansion

3. 3 The (Ordinary) Maxwell Model

4. 3.1 Wave-front Expansion

5. 3.2 Long time asymptotic expansion

6. 3.3 Numerical Results

7. 4 The Fractional Maxwell Model of order $\alpha$

8. 4.1 Wave-front expansion for the general case

9. 4.2 The Fractional Maxwell Model of order 1/2

10. 4.2.1 Wave-front Expansion

11. 4.2.2 Long time asymptotic expansion

12. 4.2.3 Numerical Results

13. 4.3 The Fractional Maxwell Model of order 3/4

14. 4.3.1 Wave-front Expansion

15. 4.3.2 Long time asymptotic expansion

16. 4.3.3 Numerical Results

17. 5 The (Ordinary) Voigt Model

18. 5.1 Wave-front Expansion

19. 5.2 Long time asymptotic expansion

20. 5.3 Numerical Results

21. 6 The Fractional Voigt Model of order $\alpha$

22. 6.1 Wave-front expansion for the general case

23. 6.2 The Fractional Voigt Model of order 1/2

24. 6.2.1 Wave-front Expansion

25. 6.2.2 Long time asymptotic expansion

26. 7 Conclusions

1. Introduction

In the last three decades, Linear Viscoelasticity [2; 10; 13; 14; 22] has presented itself as a fundamental tool for modelling various physical situations, from seismology to biophysics (see e.g. [5; 8; 23]). In this framework, models involving fractional derivatives in the constitutive equations play an important role, allowing for a natural formal background for describing systems with naturally occurring memory effects (see e.g. [6; 13; 22; 24; 29; 30; 31]). There are also a few models that show a continuous transition from a purely fractional nature to an ordinary one, which appear to be of particular interest for their applications in hemodynamics (see [3; 8]).

An important phenomenon that tend to emerge in various branches of physics is the appearance of transient effects, i.e. processes that lead to a burst of energy in a system, caused by a sudden change of state, that last for a very short time. These processes emerge in the description of mechanical system as transient waves. It is therefore very important to understand how these short-lived waves propagate within certain media.

The first attempt to formally describe the behaviour of the wave-front of a transient wave, as it propagates in a viscoelastic material, was proposed in the seventies by Buchen and Mainardi, with the seminal paper [1]. Recently, the so called formalism of the Buchen-Mainardi wave-front expansion, was applied to the Bessel models of linear viscoelasticity (see [3]). It is worth remarking that the result of this analysis (see [5]) has shown a consistency between the model proposed in [8] by Giusti and Mainardi, and the expected behaviour of a pressure wave as it propagates within a large artery.

It is worth remarking that, in recent years, many authors have tackled the problem of the wave-front expansion in viscoelastic media. Among the many who have worked on this matter (see e.g. [20; 21; 25; 27]), we feel appropriate to highlight the remarkable results obtained by Hanyga and Seredynska, for further details see e.g. [14; 15; 17].

The aim of this paper is then to review and generalize the Buchen-Mainardi algorithm for the wave-front expansion in viscoelastic media. Specifically, after a thorough review of the algorithm, we show an explicit implementation of the formalism for all the main fundamental model of both ordinary and fractional linear viscoelasticity. The article is therefore organized as follows:

In Section 2, we thoroughly revisit and generalize the Buchen-Mainardi algorithm.

In Section 3 and 4 we discuss the formalism for the Maxwell model and to the fractional Maxwell model of general order $\alpha$. Then, we explicitly solve the fractional case for two specific realizations of the model, precisely $\alpha=1/2,\,3/4$.

In Section 5 and 6 we then lay the groundwork for the application of the formalism to the Voigt model and to the fractional Voigt model of general order $\alpha$. Then, for sake of brevity, we explicitly present the discussion for the most interesting fractional realization of the model, i.e. $\alpha=1/2$.

Finally, we conclude the paper presenting some final considerations.

2. Buchen-Mainardi wave-front expansion

In this section we wish to clarify and slightly generalize the renown algorithm for performing asymptotic wave-front expansions for singular viscoelastic models, first presented by Buchen and Mainardi in [1]. Let us consider an homogeneous and isotropic viscoelastic medium of density $\rho$. This body is assumed to be semi-infinite in extent (i.e. $x\geq 0$) and unperturbed for $t<0$.

For $t\geq 0$, the accessible portion of the body (i.e. $x=0$) is subjected to a perturbation (an input signal) denoted by $r_{0}(t)$. The problem now is to compute the response $r(t,x)$ of the material for $x>0$. This response function can be either the stress or the strain, or some other quantities that depend on the specific implementation of the model (for example, in [1] the response function could also represent the displacement $u(t,x)$ or the particle velocity $\dot{u}(t,x)$).

Let $J(t)$ be the uniaxial creep compliance corresponding to our viscoelastic body, for which we set $J_{0}\equiv J(0^{+})$. Moreover, defining the wave-front velocity as

[TABLE]

and denoting with $\Psi(t)$ the creep memory function, defined as

[TABLE]

thus we have that the equation of motion reads (see e.g. [22])

[TABLE]

where the $\ast$ denotes the convolution in the Laplace sense.

Now, taking the Laplace transform of Eq. (2.1), switching from the time domain to the Laplace domain, we get the following second order differential equation,

[TABLE]

where

[TABLE]

Taking into account the boundary conditions, one can easily find a formal solution for Eq. (2.2) in the Laplace domain, indeed

[TABLE]

As discussed in [1], the aim of this procedure is to compute the asymptotic expansion of $r(t,x)$ in the neighborhood of the pulse onset. This is done by inverting, term by term, the expansion of (2.4) as $s\to\infty$ that, taking profit of Watson’s lemma or by means of the Tauberian theorems (see e.g. [26]), will correspond to an expansion for $t\to(x/c)^{+}$.

Let us first discuss the behavior of $\mu(s)$ as $s\to\infty$. For a singular model it is known (see e.g. [22]) that $\mu(s)$ will have an expansion in terms of decreasing powers of $s$. Thus, in a general fashion, one could write this expansion as

[TABLE]

Now, it is important to define two quantities whose relevance will appear clear later in this discussion. Then, let us denote $\mu_{+}(s)$ the sum of the first $m+1$ ($m\in\mathbb{N}_{0}$) terms of Eq. (2.5) that depend on positive powers of $s$, i.e.

[TABLE]

Moreover, we denote with $\mu_{-}(s)$ the remainder of the series (2.5), i.e.

[TABLE]

Therefore, taking profit of these definitions, we can rewrite the solution in Eq. (2.4), for $s\to\infty$, as follows

[TABLE]

with

[TABLE]

According to [1; 22], we shall refer to $\mu_{+}(s)$ as the principal part of the expansion of $\mu(s)$, which is directly connected to the asymptotic expansion of the creep compliance, either in the Laplace domain or in the time domain. Taking profit of the general asymptotic expansion of the creep compliance for a singular viscoelastic model, i.e.

[TABLE]

with $0<\alpha\leq 1$, one can easily show that

[TABLE]

It is worth remarking that this restriction on the values of $\alpha$ allows for a finite propagation speed of the wave-front, as discussed in [15].

Given that

[TABLE]

we can expand $\widetilde{R}(s,x)$ in terms of decreasing powers of $s$ as follows,

[TABLE]

for some functions $v_{k}(x)$ that will encode the dependence on the spatial direction.

Then, plugging the expression in Eq. (2.8) into Eq. (2.2) one finds that $\widetilde{R}(s,x)$ should also satisfy the following differential equation

[TABLE]

together with the condition

[TABLE]

In order to perform such an expansion for $\widetilde{R}(s,x)$ we should take profit of the following theorem by Friedlander and Keller, that states

Theorem 1** ([7]).**

Let $\mathbb{L}$ be a linear differential operator such that its asymptotic expansion with respect to a parameter $\varepsilon$, as $\varepsilon\to 0$, is given by

[TABLE]

where $\mathbb{L}_{i}$ are some suitable linear differential operators.

Let also assume $f$ to be a solution of

[TABLE]

with an asymptotic expansion as $\varepsilon\to 0$

[TABLE]

If the asymptotic expansion of $\mathbb{L}$ is given by termwise application of the $\mathbb{L}_{i}$ to the $v_{j}$ and if $\mathbb{L}_{0}v_{k}\neq 0$ for $k>0$, then the coefficients $v_{k}$ satisfy the following recursive system of equations

[TABLE]

where the summation is taken over values of $i,\,j$ such that $\nu_{i}+\lambda_{j}=\lambda_{k}$. The value of $\lambda_{0}$ is arbitrary, but $\lambda_{k}$ for $k>0$ is the $(k+1)$-st number in the increasing sequence formed from the set $\lambda_{0}+\sum_{i=1}^{\infty}m_{i}\,\nu_{i}$, where $m_{i}$ are non-negative integers.

Remark*.*

Notice that $j$ is a dummy index, indeed from the condition $\nu_{i}+\lambda_{j}=\lambda_{k}$ one can easily infer that $j=j(i,k)$. Therefore, one can rewrite the previous system as

[TABLE]

for some function $j$ of $i$ and $k$ defined by the relation between the exponents of the asymptotic expansions.

If we now expand $\mu(s)$ and $\mu_{+}(s)$ in powers of $s$, as $s\to\infty$, in Eq. (2.13) and then divide the whole equation for the term containing the highest power of $s$, we get a differential operator of the form required by the previous theorem. In particular, for these specific application we will have that $\varepsilon=1/s$ (which actually goes to zero as $s\to\infty$) and the operators $\mathbb{L}_{i}$ will be some functions of the derivative with respect to $x$. The functional dependence on $\partial/\partial x$ of these operators will depend on the expansion of $\mu(s)$.

If we implement the previous discussion, we obtain quite generally

[TABLE]

where the coefficients $p_{i}$, $q_{i}$ and $r_{i}$, together with the exponents $\nu_{i}$, are determined by the asymptotic expansion of $\mu(s)$ and, therefore, from the asymptotic behavior of $s\,\widetilde{J}(s)$ as $s\to\infty$.

Now, if we impose the condition $\widetilde{R}(s,0)=0$ on the asymptotic expansion of $\widetilde{R}(s,x)$ in Eq. (2.12) we infer that

[TABLE]

where $\delta_{kh}$ is the usual Kronecker delta function.

Hence, from the theorem above we get that

[TABLE]

where $N$ is the label of the last non-vanishing sub-operators $\mathbb{L}_{i}$ in the asymptotic expansion of $\mathbb{L}$, and that the coefficients $v_{k}(x)$ will be given by the recursive system of equations

[TABLE]

where the function $j(i,k)$ is deduced by the condition $\nu_{i}+\lambda_{j}=\lambda_{k}$.

It is easily understood that the solution for this system of equations is given by functions $v_{k}(x)$ such that $v_{k}(x)\in\mathbb{R}_{k}[x]$, i.e.

[TABLE]

Plugging these polynomials into Eq. (2.17) we get a recursive system for the coefficients $A_{k,\ell}$,

[TABLE]

where we decide to put a “redundant” sum over a dummy index $J$ in order to lighten the notation, avoiding cumbersome subscripts like $A_{j(i,k),\,\ell}$.

Now, plugging the asymptotic expansion in Eq. (2.12) into Eq. (2.8), we have that

[TABLE]

Then, defining

[TABLE]

with the convention $1/c=0$ when $J_{0}=0$, we get

[TABLE]

For a mathematical discussion of the functions $\widetilde{\Phi}_{k}(s,x)$ we invite the interested reader to refer to [1; 22].

Now, inverting the result in Eq. (2.23) back to the time domain (see [1; 22]) we get the asymptotic expansion for $r(t,x)$ as $t\to(x/c)^{+}$, i.e.

[TABLE]

where the functions $v_{k}(x)$ are computed by means of Eq. (2.18) and (2.19) and the functions $\Phi_{k}(t,x)$ are the Laplace inverse of the functions $\widetilde{\Phi}_{k}(s,x)$ in Eq. (2.22).

3. The (Ordinary) Maxwell Model

The ordinary Maxwell model of linear viscoelasticity is defined in terms of its constitutive equation

[TABLE]

where $a_{1}$ and $b_{1}$ are strictly positive constants.

Applying the Laplace transform to both sides we get

[TABLE]

thus,

[TABLE]

from which we infer that

[TABLE]

It is interesting to remark that the creep compliance for the Maxwell model in the time domain is easily obtained by inverting the previous expression, indeed

[TABLE]

where $J_{0}\equiv a_{1}/b_{1}$ and $J_{1}\equiv 1/b_{1}$.

3.1. Wave-front Expansion

For the (ordinary) Maxwell model we have that

[TABLE]

where we have taken profit of the relation $c=1/\sqrt{\rho\,J_{0}}$.

Now, as $s\to\infty$ we have that

[TABLE]

From this expression one can easily read off $\mu_{+}(s)$ and $\mu_{-}(s)$, i.e.

[TABLE]

Then,

[TABLE]

from which we get

[TABLE]

Therefore, the differential operator $\mathscr{O}$ is given by

[TABLE]

Then, the rescaled operator is easily computed, i.e.

[TABLE]

Hence, $\mathbb{L}$ can be expressed in terms of decreasing powers of $s$ as follows,

[TABLE]

where

[TABLE]

from which we can also deduce that $\nu_{0}=0$, $\nu_{1}=1$ and $N=1$.

Now, the coefficients $\lambda_{k}$ are easily computed from the previous results, indeed

[TABLE]

thus, $\lambda_{k}=k$ for $k\in\mathbb{N}$. Moreover, from the condition

[TABLE]

we immediately infer that $j=k-\nu_{i}$.

Then, taking profit of the system in Eq. (2.19), we can compute the values of the coefficients $A_{k,\ell}$ for the Maxwell model. In particular, given the previous discussion we have that the system in Eq. (2.19) can be rewritten as

[TABLE]

In particular, the second line can be rewritten, dropping the sum and the delta, as

[TABLE]

Thus, accordingly with the discussion in Section 2, we can immediately read off the values of all $p_{i}$’s, $q_{i}$’s and $r_{i}$’s from the operator $\mathbb{L}_{1}$, indeed

[TABLE]

whereas $p_{i}=q_{i}=r_{i}=0\,,\,\,\forall i\geq 2$.

If we focus on the coefficients such that $1\leq l\leq k$ we get

[TABLE]

Thus,

[TABLE]

Considering an initial step input (i.e. $r_{0}(s)=1/s$), then the function $\widetilde{\Phi}_{k}(s,x)$ is given by

[TABLE]

that, inverting it back to the time domain, gives

[TABLE]

Given the previous results, we are now capable to write down the explicit wave-front expansion for the Maxwell model, i.e.

[TABLE]

where the coefficients $A_{k,\ell}$ are determined as in Eq. (3.23).

3.2. Long time asymptotic expansion

From Eq. (3.6) one can easily deduce that

[TABLE]

thus,

[TABLE]

Now, inverting back to the time domain we get

[TABLE]

where erfc stands for the complementary error function.

3.3. Numerical Results

In the following we show separately the plots of both the asymptotic expansions for the ordinary Maxwell model, at fixed values of $x$. Then, we highlight an explicit matching between the wave-front expansion and the long time expansion, setting $\rho=J_{0}=J_{1}=1$.

4. The Fractional Maxwell Model of order $\alpha$

Similarly to the ordinary case, the fractional Maxwell model of linear viscoelasticity is defined by means of its constitutive equation

[TABLE]

where the constants $a_{1}$ and $b_{1}$ are required to be strictly positive.

Here, the time derivatives are replaced by the well known fractional Caputo derivative of order $\alpha\in\mathbb{R}^{+}$ (see [22]), defined, for a certain function of time $f(t)$, as

[TABLE]

that, assuming appropriate initial condition of physical interest for viscoelastic materials, can be replaced as follows when the Laplace transform is applied to both sides of a constitutive equation

[TABLE]

for further details on the initial conditions we invite the reader to refer to [9; 13; 22].

Then we can finally infer that

[TABLE]

where $\tau^{\alpha}=a_{1}$.

4.1. Wave-front expansion for the general case

For the fractional Maxwell case we find

[TABLE]

Again, as $s\to\infty$ one can easily find that

[TABLE]

where we have taken profit of the relations $c=\sqrt{b_{1}/(\rho\,a_{1})}$ and $\tau^{\alpha}=a_{1}$.

At this point, since $\mu_{+}(s)$ has to contain only positive powers of $s$, it is recommended to distinguish between different cases. As $1-n\alpha\geq 0$, then we have to satisfy the condition $n\leq[1/\alpha]$. So, if $\alpha=1/2$, we find that $n=1,2$, if $\alpha>1/2$, $n=1$ and if $\alpha<1/2$ there are several values of $n$ (e.g. $\alpha=1/4$ implies $n=1,2,3,4$).

4.2. The Fractional Maxwell Model of order 1/2

If we fix $\alpha=1/2$, we have that the (4.4) becomes

[TABLE]

4.2.1. Wave-front Expansion

As a consequence, we can perform the wave-front expansion, and we have

[TABLE]

From the previous general case, as $s\to\infty$ we find that

[TABLE]

and, from the latter, we can easily read

[TABLE]

Then,

[TABLE]

from which we get

[TABLE]

Therefore, the differential operator $\mathscr{O}$, defined in (2.13), is given by

[TABLE]

At this point, the rescaled operator is easily computed, i.e.

[TABLE]

Hence, $\mathbb{L}$ can be expressed in terms of decreasing powers of $s$ as follows,

[TABLE]

where

[TABLE]

and, from which we can also deduce that $\nu_{0}=0$, $\nu_{1}=1/2$, $\nu_{2}=1$ and $N=2$.

Now, the coefficients $\lambda_{k}$ are easily computed taking into account the previous results, indeed

[TABLE]

thus, $\lambda_{k}=k/2$ for $k\in\mathbb{N}$. Moreover, from the condition

[TABLE]

we immediately infer that $j=k-2\,\nu_{i}$.

Thus, following the algorithm in Section 2, we can immediately read off the values of all $p_{i}$’s, $q_{i}$’s and $r_{i}$’s from the expressions for $\mathbb{L}_{1}$ and $\mathbb{L}_{2}$, indeed

[TABLE]

while $p_{i}=q_{i}=r_{i}=0\,,\,\,\forall i\geq 3$.

Finally, if we focus on the coefficients such that $1\leq l\leq k$ we obtain the relation

[TABLE]

Considering an initial step input (i.e. $r_{0}(s)=1/s$), then the function $\widetilde{\Phi}_{k}(s,x)$ is given by

[TABLE]

that, inverting it back to the time domain, gives

[TABLE]

Hence, given the previous results, the wave-front expansion for the fractional Maxwell model of order $1/2$ is given by

[TABLE]

as $t\to\left(x/c\right)^{+}$.

4.2.2. Long time asymptotic expansion

From Eq. (4.8) one can easily deduce that

[TABLE]

thus,

[TABLE]

Thereafter, inverting the latter expression back to the time domain we get [28]

[TABLE]

where $W_{\gamma,\delta}$ represents the Wright function

[TABLE]

which satisfies

[TABLE]

For further information about Wright functions see [12].

4.2.3. Numerical Results

In the following we show separately the plots of both the asymptotic expansions for the fractional Maxwell model of order $1/2$. Similarly to the ordinary Maxwell model, we consider plots at fixed position $x$. Then, we show an explicit matching between the wave-front expansion and the long time expansion, setting $a_{1}=b_{1}=\rho=1$ so that $c=1$.

4.3. The Fractional Maxwell Model of order 3/4

Now we set $\alpha=3/4$, so that the creep compliance in the Laplace domain reads

[TABLE]

and than we repeat the procedure discussed above.

4.3.1. Wave-front Expansion

First write

[TABLE]

Secondly, we compute the limit for $s\to\infty$ of (4.34), namely

[TABLE]

Then, from the previous expression, we immediately read off

[TABLE]

also inferring that

[TABLE]

from which we get

[TABLE]

Hence, the differential operator $\mathscr{O}$ is given by

[TABLE]

that can be rescaled as

[TABLE]

At this point, $\mathbb{L}$ can be expressed in terms of decreasing powers of $s$ as follows,

[TABLE]

where

[TABLE]

and from these expansions we immediately deduce that $\nu_{0}=0$, $\nu_{1}=1/2$, $\nu_{2}=3/4$, $\nu_{3}=1$ and $N=3$.

Given that, the coefficients $\lambda_{k}$ are easily computed from the previous results, so that

[TABLE]

accordingly, $\lambda_{k}=k/4$ for $k\in\mathbb{N}$. Moreover, considering the condition

[TABLE]

we conclude that $j=k-4\,\nu_{i}$.

Thus, following the procedure described by the algorithm, we can directly realise the values of all the coefficients $p_{i}$’s, $q_{i}$’s and $r_{i}$’s from the expressions for every $\mathbb{L}_{i}$, indeed

[TABLE]

while $p_{i}=q_{i}=r_{i}=0\,,\,\,\forall i\geq 4$.

Focusing our attention on the coefficients such that $1\leq l\leq k$, we get

[TABLE]

Now, we consider an initial step input (i.e. $r_{0}(s)=1/s$), so that the function $\widetilde{\Phi}_{k}(s,x)$ is given by

[TABLE]

that, inverting back to the time domain, gives

[TABLE]

Then, from Eq. (2.24) we have that the wave-front expansion for the fractional Maxwell model of order $3/4$ is given by

[TABLE]

as one approaches the wave-front. It is important to recall that $W_{\gamma,\delta}$ represents the Wright function introduced above for the fractional Maxwell of order $1/2$.

4.3.2. Long time asymptotic expansion

From Eq. (4.34) one can easily deduce that

[TABLE]

so that,

[TABLE]

Now, inverting the latter back to the time domain, we finally get

[TABLE]

4.3.3. Numerical Results

In the following we show separately the plots of both the asymptotic expansions for the fractional Maxwell model of order $3/4$. Then, we show an explicit matching between the wave-front expansion and the long time expansion, setting $a_{1}=b_{1}=\rho=1$ so that $c=1$.

5. The (Ordinary) Voigt Model

The ordinary Voigt model is defined in terms of its constitutive equation that reads

[TABLE]

where $m\geq 0$ and $b_{1}$ is a strictly positive constant.

In the Laplace domain the latter is rewritten as

[TABLE]

from which we can immediately read off the Laplace transform of the creep compliance,

[TABLE]

where $J_{1}=1/m$ and $\tau_{\varepsilon}=b_{1}/m$.

For sake of completeness, the creep compliance defined in Eq. (5.3) can be inverted back to the time domain and it gives

[TABLE]

from which one can infer that the wave-front velocity is infinite (i.e. $J_{0}=0$).

5.1. Wave-front Expansion

Here we have that

[TABLE]

and from asymptotic behaviour for $s\to\infty$ we have

[TABLE]

Hence, it is straightforward to see that

[TABLE]

Moreover, we have that

[TABLE]

Consequently, the last expressions tell us that

[TABLE]

Now, the differential operator $\mathscr{O}$ for the ordinary Voigt model is given by

[TABLE]

If we then rescale $\mathscr{O}$ as follows

[TABLE]

we have that both $\mathscr{O}$ and $\mathbb{L}$ enjoy the same homogeneous solutions111Which are the physically significant one, as discussed in Section 2., provided that we keep ourself away form the singularities of $\mathscr{O}$ in the Laplace domain. Clearly, this rescaled operator is expressed in terms of negative powers of $s$ as follows

[TABLE]

where

[TABLE]

Furthermore, from Eq. (5.13) we can also infer that $\nu_{1}=1/2$, $\nu_{2}=1$, $\nu_{3}=3/2$ and that $N=3$. Thus from the theorem by Friedlander and Keller we deduce that

[TABLE]

Therefore, if we proceed following the general procedure of the Buchen-Mainardi Algorithm, we can then compute the coefficients $A_{k,\ell}$. In particular, for the coefficients for which $1\leq l\leq k$ we find that

[TABLE]

Considering an initial step input (i.e. $r_{0}(s)=1/s$), then the function $\widetilde{\Phi}_{k}(s,x)$ reads

[TABLE]

that, inverting back to the time domain, gives

[TABLE]

At this point, we are able to write the wave-front expansion for the Voigt model, namely

[TABLE]

where the coefficients $A_{k,\ell}$ are determined as follows

[TABLE]

5.2. Long time asymptotic expansion

From Eq. (5.5) we can immediatly deduce that

[TABLE]

The inversion of this function to the time domain, as $s\to 0$ or equivalently as $t\to\infty$, can be performed by means of the saddle-point approximation method. Indeed, as discussed in [19], for the Voigt model we get the following long time asymptotic expansion,

[TABLE]

5.3. Numerical Results

Here are shown the plots of both the asymptotic expansions for the ordinary Voigt model, computed at different fixed values of $x$. Then, we underline a numerical matching between the wave-front expansion and the long time expansion. In particular, in order to roughly show the matching, we perform a polynomial interpolation between the two approximations.

6. The Fractional Voigt Model of order $\alpha$

The constitutive equation fractional Voigt model reads

[TABLE]

where $m\geq 0$ and $b_{1}$ is a strictly positive constant, as for the ordinary case. As for the fractional Maxwell model, we have replaced the ordinary time derivatives with the Caputo’s fractional derivative, as defined in (4.2).

In the Laplace domain we can rewrite the constitutive equation as follows

[TABLE]

from which we immediately read off the Laplace transform of the creep compliance,

[TABLE]

where $J_{1}=1/m$ and $\tau_{\varepsilon}=b_{1}/m$.

6.1. Wave-front expansion for the general case

Following the same procedure extensively discussed above, we have

[TABLE]

Now, as $s\to\infty$ we find that

[TABLE]

Again, at this point, we distinguish different cases as shown for the fractional Maxwell model. As $1-\alpha(n+1/2)\geq 0$, then we have to satisfy the condition $n\leq[1/\alpha-1/2]$.

In this paper, we only consider the fractional Voigt model of order $\alpha=1/2$, so that $n\leq 2$ in our expansion.

6.2. The Fractional Voigt Model of order 1/2

As said before, fixing $\alpha=1/2$, Eq. (6.3) turns into

[TABLE]

6.2.1. Wave-front Expansion

Then,

[TABLE]

that, as $s\to\infty$, reads

[TABLE]

From (6.8) it is straightforward to evaluate the expression for $\mu_{+}(s)$, namely

[TABLE]

and, consequently,

[TABLE]

thus, we get

[TABLE]

Therefore, the differential operator $\mathscr{O}$ is given by

[TABLE]

Then, rescaling $\mathscr{O}$, we find $\mathbb{L}$, i.e.

[TABLE]

where

[TABLE]

Moreover, from Eq. (6.14) we can also infer that $\nu_{1}=1/2$, $\nu_{2}=3/4$, $\nu_{3}=1$, $\nu_{4}=5/4$ and that $N=4$. Again, thanks to the theorem by Friedlander and Keller we deduce that

[TABLE]

Hence, the coefficients $A_{k,\ell}$ are easily computed. Specifically, the coefficients for which $1\leq l\leq k$ are given by

[TABLE]

Considering an initial step input (i.e. $r_{0}(s)=1/s$), then the function $\widetilde{\Phi}_{k}(s,x)$ is given by

[TABLE]

that it is not easy to invert back, at least analytically, to the time domain.

6.2.2. Long time asymptotic expansion

From Eq. (6.7) one can easily infer that

[TABLE]

The inversion of this function to the time domain, as $s\to 0$ or equivalently as $t\to\infty$, appears to be particularly involved. Nevertheless, if we consider only the first term of the expansion, we get the following long time asymptotic behaviour,

[TABLE]

where $\Theta$ is the Heaviside step function. A numerical inversion of Eq. (6.23) has shown negligible deviations from Eq. (6.24) once further terms of the expansion are taken into account.

7. Conclusions

In this article we revisited and generalized the wave-front formalism for viscoelatic transient phenomena, first proposed by Buchen and Mainardi in [1], also providing a wide variety of applications to the most relevant linear viscoelastic models.

Specifically, in Section 2 we outlined a way to set the problem of the production and propagation of transient waves in (linear) viscoelastic media. Then, we precisely derived the Buchen-Mainardi algorithm in a very general fashion, also laying out the procedure for generic (linear) media. In Section 3 we then performed the wave-front expansion for an ordinary Maxwell solid considering, as source of the perturbations of the material, an instantaneous step input. Subsequently, we computed the long time response of the material and we compared the plots of the two approximations. The result of this analysis is a nice matching of the two asymptotic expansions.

In Section 4 we then tackled the same problem in the framework of the fractional Maxwell model. The section begins with a discussion of the wave-front expansion for a generic fractional Maxwell model of order $\alpha\in\mathbb{R}^{+}$. Then, in order to be able to perform some explicit computations, we decided to inspect two realizations of the model, respectively $\alpha=1/2$ and $\alpha=3/4$. Although this analysis shows that we still have a quite nice matching between the short time and the long time approximations, this seems to emerge on a time scale which is way shorter than the one that we found for the ordinary model. Nevertheless, the matching still appears within the range of validity of the wave-front approximation (which is obtained from a truncation of the series in Eq. (2.24) to the 30th term).

In Section 5 we turned our attention to the ordinary Voigt model. As above, after computing both the wave-front expansion and the long time approximation (the latter was obtained by Jeffreys in [19]), we graphically compared the two asymptotic behaviors. Unfortunately, the wave-front expansion breaks way before showing any explicit matching, and this is independent from the number of terms considered for the series in Eq. (2.24) (aside from minor corrections).

Finally, in Section 6 we performed the same procedure for the fractional Voigt model of order $1/2$. We decided not to show plots for the two approximations because, as already seen in the Maxwell models, the wave-front expansion breaks even earlier than the one of the ordinary model. Moreover, the long time approximation shows very little deviations from the step function.

As a final remark, it is important to stress that the long time expansions presented in this paper have been obtained following a very simplified approach. The reason for this is that the main focus of the paper is on the (short time) wave-front expansion, and the comparison between the two asymptotic approximation is used to qualitatively infer to which extent the latter can be trusted, in some very simple cases. Clearly, if we wanted to provide a more precise description of the long time expansion we should recur to some more sophisticated such as the one described in [16; 27]. Furthermore, it is also important to remark, for sake of completeness, that in [16; 27] the authors address an important result by Kolsky [20] concerning a peculiar universality of viscoelastic pulse shapes.

Acknowledgments

The authors acknowledge the two anonymous reviewers for the constructive comments and suggestions which have helped to improve the manuscript significantly.

The work of the authors has been carried out in the framework of the activities of the National Group of Mathematical Physics (GNFM, INdAM).

Moreover, the work of A.G. has been partially supported by GNFM/INdAM Young Researchers Project 2017 “Analysis of Complex Biological Systems”.