Zener Model with General Fractional Calculus: Thermodynamical Restrictions

Abstract

We studied a Zener-type model of a viscoelastic body within the context of general fractional calculus and derived restrictions on coefficients that follow from the dissipation inequality, which is the entropy inequality under isothermal conditions. We showed, for a stress relaxation and a wave propagation, that the restriction that follows from the entropy inequality is sufficient to guarantee the existence and uniqueness of the solution. We presented numerical data related to the solution of a wave equation for several values of parameters.

1. Introduction

Concerning various applications of the “classical” Riemann–Liouville and Caputo fractional derivatives in various problems of physics and mechanics, there exists a very large amount of literature (see for example [1,2,3,4,5] and references therein). On the other hand, the concept of general fractional calculus (GFC) has recently gained an increasing interest from a theoretical point of view. We mention the review article of V. E. Tarasov [6], where a development of ideas concerning the generalization of the Riemann–Liouville and Caputo fractional integrals and derivatives are presented. The idea of GFC may be traced to the work of A. N. Kochubei [7]. Further developments may be found in the papers of Y. Luchko [8,9], as well as in [1,2].

The main idea of the GFC concept is to describe dynamical systems with the nonlocality in time and space.

The following definitions of GFC are introduced; see [6,8,10] and references given therein:

$$\begin{array}{ccc}\hfill I_{\left(M\right)}^t\left[\tau \right]f\left(\tau \right)& =& \left(M\ast f\right)\left(t\right)={\int }_0^{t}M\left(t-\tau \right)f\left(\tau \right)d\tau ,\hfill \\ \hfill D_{\left(K\right)}^t\left[\tau \right]f\left[\tau \right]& =& \frac{d}{dt}\left(K\ast f\right)\left(t\right)=\frac{d}{dt}{\int }_0^{t}K\left(t-\tau \right)f\left(\tau \right)d\tau ,\hfill \\ \hfill {}^C{D}_{\left(K\right)}^t\left[\tau \right]f\left[\tau \right]& =& \left(K\ast f^{\left(1\right)}\right)\left(t\right)={\int }_0^{t}K\left(t-\tau \right)f^{\left(1\right)}\left(\tau \right)d\tau ,t\ge 0,\hfill \end{array}$$

(1)

where $f^{\left(1\right)}\left(t\right)=\frac{df}{dt}$ and $f\ast g$ denotes convolution. Equation (1)${}_1$ defines the integral in GFC, whereas (1)${}_2$ and (1)${}_3$ define the generalized fractional derivative in the Riemann–Liouville sense, and generalized fractional derivative in the Caputo sense, respectively. Kernels M and K in (1) belong to the class of Sonin kernels, defined as follows. Let $a<b$. Then,

$$C_{\left(a,b\right)}\left(0,\infty \right)=\left\{f:f\left(t\right)=t^{p}Y\left(t\right),a<p<b,Y\left(t\right)\in C[0,\infty )\right\},$$

where $C[0,\infty )$ is the space of continuous functions on $[0,\infty )$. Then, the set of Sonin kernels $S_{-1}$ is defined as a set of pairs of functions $\left(M\left(t\right),K\left(t\right)\right)$, such that $M\left(t\right),K\left(t\right)\in C_{\left(-1,0\right)}\left(0,\infty \right),$ which satisfies the so-called Sonin condition

$$\left(M\ast K\right)\left(t\right)={\int }_0^{t}M\left(t-\tau \right)K\left(\tau \right)d\tau =1,t\ge 0;$$

(2)

$\left(M\left(t\right),K\left(t\right)\right)$ is called a Sonin pair. By [11], if $M\left(t\right)$ and $K\left(t\right),t\in (0,\infty )$, we have locally integrable derivatives (that belong to $L_{loc}^1\left({\mathbb{R}}_{+}\right)$) that satisfy

$$\underset{t\to 0}{lim}tM\left(t\right)=0,\underset{t\to 0}{lim}tK\left(t\right)=0.$$

Then,

$${\int }_0^t{M}^{\left(1\right)}\left(\tau \right)\left[K\left(t\right)-K\left(t-\tau \right)\right]d\tau =K\left(t\right)M\left(t\right),$$

(3)

for all $t>0$, where $M^{\left(1\right)}$ and $K^{\left(1\right)}$ exist. Moreover, (2) implies

$$\begin{array}{ccc}\hfill D_{\left(K\right)}^t{I}_{\left(M\right)}^t\left[\tau \right]f\left(\tau \right)& =& f\left(t\right),\hfill \\ \hfill I_{\left(M\right)}^t{D}_{\left(K\right)}^t\left[\tau \right]f\left[\tau \right]& =& f\left(t\right)-f\left(0\right),t>0;\hfill \end{array}$$

(4)

see [7,12]. Then, (2) implies

$$s\tilde{M}\left(s\right)\tilde{K}\left(s\right)=1,\Re s>0,$$

(5)

where $\widetilde{(·)}$ denotes the Laplace transform. This is used in [7] for the construction of M when K is given.

We note that, in the frame of the Riemann–Liouville fractional calculus, one has

$$M\left(t\right)=\frac{t^{\alpha -1}}{\Gamma \left(\alpha \right)},K\left(t\right)=\frac{t^{-\alpha }}{\Gamma \left(1-\alpha \right)},t>0,$$

so that (2), (3) and (5) hold true.

Our aim in this work is to analyze the Zener model of a viscoelastic body with specific K and M satisfying (2). Moreover, we will derive the thermodynamical admissibility conditions that guarantee the dissipation.

We consider kernels M and K in two forms, called cases. The first one is proposed by Hanyga [12] as

Case H

$$\begin{array}{ccc}\hfill M_1\left(t\right)& =& \frac{t^{-\beta }}{\Gamma \left(1-\beta \right)}+\frac{t^{\alpha -\beta }}{\Gamma \left(\alpha -\beta +1\right)},t>0,\hfill \\ \hfill K_1\left(t\right)& =& t^{\beta -1}{E}_{\alpha ,\beta }\left(-t^{\alpha }\right),\beta \ge 0,0<\alpha \le \beta ,t>0,\hfill \end{array}$$

(6)

where $E_{\alpha ,\beta }\left(t\right)={\displaystyle \sum _{k=0}^{\infty }}\frac{t^k}{\Gamma \left(k\alpha +\beta \right)},t\ge 0,\alpha >0,\beta \in \mathbb{C},$ is a two-parameter Mittag–Lefler function [13]. Conditions $\beta \ge 0$ and $0<\alpha \le \beta$ imply that K is a singular, locally integrable completely monotone function; see [14], p.144. It was shown in [12], Theorem 4.1, that any singular, unbounded in a neighborhood of zero, locally integrable, completely monotone function is a Sonin kernel, i.e., satisfies (5). Thus, functions (6) make a pair that belongs to $S_{-1}$. The Laplace transforms of $M_1$ and $K_1$ are

$${\tilde{M}}_1\left(s\right)=\frac{s^{\alpha }+1}{s^{1+\left(\alpha -\beta \right)}},{\tilde{K}}_1\left(s\right)=\frac{s^{\alpha -\beta }}{s^{\alpha }+1},\Re s>0.$$

(7)

The second Sonin pair, proposed by Zacher [15], is

Case Z

$$\begin{array}{ccc}\hfill M_2\left(t\right)& =& \left[\frac{t^{\alpha -1}}{\Gamma \left(\alpha \right)}exp\left(-\mu t\right)+\mu {\int }_0^t\frac{{\tau }^{\alpha -1}}{\Gamma \left(\alpha \right)}exp\left(-\mu \tau \right)d\tau \right]\hfill \\ \hfill K_2\left(t\right)& =& \frac{t^{-\alpha }}{\Gamma \left(1-\alpha \right)}exp\left(-\mu t\right),,t\ge 0\hfill \end{array}$$

(8)

where $\alpha \in \left(0,1\right),\mu \ge 0$. The Laplace transforms of $M_2$ and $K_2$ are

$${\tilde{M}}_2\left(s\right)=\frac{1}{{\left(s+\mu \right)}^{\alpha }}\left(1+\frac{\mu }{s}\right),{\tilde{K}}_2\left(s\right)=\frac{1}{{\left(s+\mu \right)}^{1-\alpha }},\Re s>0.$$

(9)

Our goal is to investigate mechanical models with constitutive equations in the framework of GFC.

More precisely, we consider the one-dimensional generalized Zener model [16] given as

$$a{}^C{D}_{\left(K_i\right)}^t\sigma \left(x,t\right)+\sigma \left(x,t\right)=b{}^C{D}_{\left(K_i\right)}^t\epsilon \left(x,t\right)+\epsilon \left(x,t\right),t\in \left(0,\infty \right),x\in \mathbb{R},$$

(10)

where $i=1,2,$ and $\sigma \left(x,t\right)$ denotes the stress $\epsilon \left(x,t\right)$, which is the strain in a body at time t and position x, and $a>0,b>0$ are constants. The case $a=b$ is trivial; the results are well known, so we do not consider it.

As in our earlier papers, Refs [17,18,19], we consider a generalized Zener constitutive equation for a viscoelastic body. We point out that our intention in this paper is to present that fractional derivatives with kernels $M_1$ and $M_2$ given above (Cases H and Z) require a completely different implementation in solving the wave equation in relation to the method of solving the same equation with the Riemann–Louvile fractional derivatives in the constitutive equation.

Concerning our previous paper [19], we have two remarks. The first one is that the obtained sufficient conditions for the thermodynamical admissibility of the constitutive equation (10) for both cases were obtained by our original approach given in the quoted paper, with the aim of analyzing the dissipations inequality. The second remark is that one can assume, as in [19], that the body force f and the initial conditions $u_0,v_0$ are random since they may incorporate trough epistemic randomness or errors in the measuring devices. Such an analysis will be considered in our future work since it involves additional extensive analysis. Moreover, the environmental noise in random fluctuations in transient dynamics of interdisciplinary physical models is an important issue that deserves further investigations; see [20,21,22,23,24] for the “classical” fractional derivatives in various applications. In several papers, the authors of [21,25]—see also references therein—use a memristor (elements of electric circuits) in the analysis of fluctuations of various nonlinear models of nanoelectronics. Their stochastic approach gives a new insight concerning physical models through the prediction of memristor behaviour.

We comment on the possible applicability of our constitutive equation (10). In the case of the Caputo derivative, in [26], the description of experimental data for certain dental materials was successful. Since we have one additional parameter $\mu$ in (8) and (10), the modeling of experimental results within the new framework is expected to be even more precise.

Let us briefly present the content of the paper. After a short introduction concerning the generalized functions framework, we analyze in Section 2 the thermodynamical restrictions for the Zener constitutive equation in cases H and Z. Starting with Section 3, up to the end of the paper, we will consider case Z with kernels (8). Section 3 is devoted to the stress relaxation and the spatially one-dimensional wave equation in a viscoelastic material when the coefficients satisfy the thermo-dynamical restriction; that is, condition $b>a.$ This condition implies the existence and the unicity of a solution of this equation. Moreover, in Section 3, we analyze the regularity of a solution and give an example confirming the numerical evidence of the properties of a solution.

Notation and Notions

We use the usual notation; for example, $\mathbb{R}=(-\infty ,\infty )$, ${\mathbb{R}}_{+}=(0,\infty ),$$C^k\left(\mathbb{R}\right)$ is the space of functions with continuous k derivatives. We refer to the classical distribution theory (cf. [27] or any other one about distributions) for the mathematical preliminaries related to the space of tempered distributions ${\mathcal{S}}^{\prime }\left(\mathbb{R}\right)$, the Fourier and the Laplace transforms. Recall that the space of smooth functions, where all derivatives rapidly decrease as $\left|x\right|\to \infty$—that is, are bounded by any power ${(1+|x\left|\right)}^{-m},m>0$—is denoted by $\mathcal{S}\left(\mathbb{R}\right)=\mathcal{S}$. Its dual is the space of tempered distributions ${\mathcal{S}}^{\prime }$. Elements $f\in {\mathcal{S}}^{\prime }$ supported by $[0,\infty )$ ($f\in {\mathcal{S}}_{+}^{\prime }\left({\mathbb{R}}^d\right))$ have the form $f=D^{k}F,D^k={\partial }^k/\partial t^k$, where F is a continuous polynomially bounded function supported by $[0,\infty )$ (equal to zero in $(0,\infty )$). The Laplace transform will be considered for the so-called exponentially bounded distributions, the linear combination of the distributions of the form $f\left(t\right)=e^{at}g\left(t\right)$, where $g\in {\mathcal{S}}_{+}^{\prime }\left(\mathbb{R}\right)$. Recall that, for the exponentially bounded function f ($f\left(t\right)\le M{e}^{at},t>0,a\in \mathbb{R}$),

$$\mathcal{L}\left(f\right)\left(s\right)=\tilde{f}\left(s\right)={\int }_0^{\infty }{e}^{-st}f\left(t\right)dt,\Re s>a.$$

The Fourier transform is an isomorphism of $\mathcal{S}$ onto $\mathcal{S}$ and of ${\mathcal{S}}^{\prime }$ onto ${\mathcal{S}}^{\prime }$. It is given by

$$\mathcal{F}\left(\varphi \right)\left(\xi \right)=\widehat{\varphi }\left(\xi \right)={\int }_{-\infty }^{\infty }{e}^{-ix\xi }f\left(x\right)dx,\xi \in \mathbb{R},\varphi \in \mathcal{S}.$$

Recall that Heaviside’s function H is the characteristic function of $[0,\infty )$. Its derivative over $\mathbb{R}$ is the delta distribution $\delta \left(x\right)$. Also in ${\mathcal{S}}_{+}\left({\mathbb{R}}_t\right)$, its derivative is $\delta \left(t\right).$ Let us recall that $\mathcal{L}\left(t_{+}\right)=1/s^2,\mathcal{L}\left(H\left(t\right)\right)=1/s,\Re s>0$, where $t_{+}=H\left(t\right)t,t\in \mathbb{R}$.

We denote the dual pairing of a test function and a distribution through the integral sign, explaining such an integral when it exists in the sense of classical functions or in the sense of distribution pairing.

In order to simplify the exposition, we will assume in the main theorem of Section 3 that a function $f(x,t),x\in \mathbb{R},t>0$, the initial data $u_0\left(x\right),v_0\left(x\right),x\in \mathbb{R}$, and their derivatives up to an imposed order belong to the space of integrable continuous functions with a suitable decay property, in order to have the existence of a convolution in x of these functions with a fundamental solution P, for which, we show that it is a distribution.

2. Thermodynamical Restrictions

We consider the constitutive Equation (10) with kernel K of both cases. Our intention is to derive the restrictions on the coefficients in (10) that follow from the second law of thermodynamics under isothermal conditions. We assume that $\epsilon (x,t)=0$ and $\sigma (x,t)=0$ for each $x\in \mathbb{R}$, and $t<0$. Further on, we assume that $\epsilon \in C^1\left([0,\infty )\right)$. The second law of thermodynamics, under appropriate isothermal conditions, requires that, for any cycle of duration $T>0$, where cycle here means $\epsilon \left(x,0\right)=\epsilon \left(x,T\right)=0$, there exists $D>0$ such that the dissipation inequality

$$D\left(x\right)={\int }_0^T\sigma (x,t)\frac{\partial \epsilon \left(x,t\right)}{\partial t}\left(t,x\right)dt\ge 0,$$

(11)

holds for every $x\in (-\infty ,\infty )$. Because of that, in the analysis that follows, we shall write (11) without x. Inequality (11) is used for $\epsilon$, which does not satisfy the conditions of a cycle. For example, in [28], the use of (11) for any sufficiently smooth $\epsilon$ is proposed, which also satisfies $\epsilon (x,t)=0,t\in (-\infty ,0]$. We follow [28] since it does not require a definition of a cycle. Note that the cycles are differently defined in various papers. In [29], it is required that the entropy inequality holds (11) is just a special case of it) for a specially defined D-cyclic process. We refer to [17,30,31,32] for a more detailed analysis of a dissipativity condition.

Applying the Fourier transform now on a function (or a distribution) depending on $t\ge 0$ (those supported by $[0,\infty )$ with respect to variable t) to (10)${}_2$, we obtain the Fourier transform of ${}^C{D}_{\left(K\right)}^{t}f\left(t\right)$ as

$$\begin{array}{ccc}\hfill \mathcal{F}\left({}^C{D}_{\left(K\right)}^{t}f\left(t\right)\right)\left(\omega \right)& =& \mathcal{F}\left({\int }_0^{t}K\left(t-\tau \right)f^{\left(1\right)}\left(\tau \right)d\tau \right)\hfill \\ & =& \mathcal{F}\left(K\left(t\right)\right)\left(\omega \right)\mathcal{F}\left(f^{\left(1\right)}\left(t\right)\right)\left(\omega \right)\hfill \\ & =& \mathcal{F}\left(\frac{t^{-\alpha }}{\Gamma \left(1-\alpha \right)}H\left(t\right)\right)\left(\omega \right)\left(i\omega \right)\mathcal{F}\left(f\left(t\right)\right)\left(\omega \right)\hfill \\ & =& \frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}\mathcal{F}\left(f\left(t\right)\right)\left(\omega \right),\omega \in \mathbb{R}.\hfill \end{array}$$

2.1. Restrictions on the Coefficients for Case Z

In this subsection, we use notation $K_2=K$ and $H_2=H$. Applying the Fourier transform (on a function or distribution depending on t), we obtain

$$\begin{array}{ccc}\hfill \mathcal{F}\left({}^C{D}_{\left(K\right)}^{t}f\left(t\right)\right)\left(\omega \right)& =& \mathcal{F}\left({\int }_0^{t}K\left(t-\tau \right)f^{\left(1\right)}\left(\tau \right)d\tau \right)\hfill \\ & =& \mathcal{F}\left(K\left(t\right)\right)\left(\omega \right)\mathcal{F}\left(f^{\left(1\right)}\left(t\right)\right)\left(\omega \right)\hfill \\ & =& \mathcal{F}\left(\frac{t^{-\alpha }}{\Gamma \left(1-\alpha \right)}H\left(t\right)\right)\left(\omega \right)\left(i\omega \right)\mathcal{F}\left(f\left(t\right)\right)\left(\omega \right)\hfill \\ & =& \frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}\mathcal{F}\left(f\left(t\right)\right)\left(\omega \right),\omega \in \mathbb{R}.\hfill \end{array}$$

Then,

$$\widehat{\sigma }\left(x,\omega \right)=E\left(\omega \right)\widehat{\epsilon }(x,\omega ),\omega \in \mathbb{R},x\in \mathbb{R},$$

(12)

implies

$$\begin{array}{ccc}\hfill E\left(\omega \right)& =& \frac{1+b\frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}}{1+a\frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}}=1+\frac{\left(b-a\right)\frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}}{1+a\frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}}\hfill \\ & =& E_1\left(\omega \right)+i{E}_2\left(\omega \right),\omega \in \mathbb{R}.\hfill \end{array}$$

(13)

Let ${\left(i\omega +\mu \right)}^{1-\alpha }=\rho exp\left(i\theta \right),\rho >0,\theta \in \left[0.2\pi \right)$ so that $i\omega ={\rho }^{\frac{1}{1-\alpha }}exp\left(\frac{i\theta }{1-\alpha }\right)-\mu$. This leads to

$${\rho }^{\frac{1}{1-\alpha }}cos\left(\frac{\theta }{1-\alpha }\right)=\mu ,\omega ={\rho }^{\frac{1}{1-\alpha }}sin\left(\frac{\theta }{1-\alpha }\right).$$

(14)

Then, we have

$$\frac{i\omega }{1+a\frac{i\omega }{{\left(i\omega +\mu \right)}^{1-\alpha }}}=\frac{i{\rho }^{\frac{1}{1-\alpha }}sin\left(\frac{\theta }{1-\alpha }\right)}{\rho cos\theta +i\rho sin\theta +ai{\rho }^{\frac{1}{1-\alpha }}sin\left(\frac{\theta }{1-\alpha }\right)}=$$

$$\frac{{\rho }^{\frac{1}{1-\alpha }}sin\left(\frac{\theta }{1-\alpha }\right)\left[\rho sin\theta +{\rho }^{\frac{1}{1-\alpha }}sin\left(\frac{\theta }{1-\alpha }\right)\right]+i\rho cos\theta {\rho }^{\frac{1}{1-\alpha }}sin\left(\frac{\theta }{1-\alpha }\right)}{{(\rho cos\theta )}^2+{(\rho sin\theta +\omega )}^2},$$

(15)

with $\rho$ and $\omega$ given by (14). Now, it is easy to prove the next proposition.

Proposition 1.

Condition $b>a$ is a sufficient one for the components $E_1$ and $E_2$ of the complex dynamic modulus E defined by (12) to satisfy

$$E_1\left(\omega \right)=E_1\left(-\omega \right),E_2\left(\omega \right)=-E_2\left(-\omega \right),\omega \in \mathbb{R},E_2\left(\omega \right)\ge 0,\omega \in {\mathbb{R}}_{+},$$

$${\int }_0^{\infty }\frac{1}{\omega }\frac{E_2\left(\omega \right)}{{(1+{\omega }^2)}^{\frac{m}{2}}}d\omega ={\int }_0^{\infty }\frac{(b-a)\rho cos\theta }{({(\rho cos\theta )}^2+{(\rho sin\theta +\omega )}^2){(1+{\omega }^2)}^{\frac{m}{2}}}d\omega <\infty ,$$

for some $m>0$.

Proof. 

Relations (13)–(15) imply that

$$\begin{array}{ccc}\hfill E_1\left(\omega \right)& =& 1+\left(b-a\right)\frac{\omega \left[\rho sin\theta +\omega \right]}{{(\rho cos\theta )}^2+{(\rho sin\theta +\omega )}^2},\hfill \\ \hfill E_2\left(\omega \right)& =& \left(b-a\right)\frac{\omega \rho cos\theta }{{(\rho cos\theta )}^2+{(\rho sin\theta +\omega )}^2},\omega \in \mathbb{R}.\hfill \end{array}$$

Since

$${(\rho cos\theta )}^2+{(\rho sin\theta +\omega )}^2={(\rho cos\theta )}^2+{(\rho sin\theta +{\rho }^{\frac{1}{1-\alpha }}sin\frac{\theta }{1-\alpha })}^2,$$

and $sin\theta$ and $sin\frac{\theta }{1-\alpha }$ have the same sign + for $\theta \in (0,\pi /2)$ and sign − for $\theta \in (-\pi /2,0)$, we can easily conclude that

$$E_1\left(\omega \right)=E_1\left(-\omega \right),E_2\left(\omega \right)=-E_2\left(\omega \right),\omega \in \mathbb{R},$$

Now, we simply conclude that

$$E_2\left(\omega \right)\ge 0,\omega >0ifb>a.$$

The last part is clear. □

2.2. Restrictions for H

In this subsection, $K_1=K,H_1=H$. We present the restrictions for (10) with the assumptions

$$a\ge 0,b>0,0<\alpha \le \beta ,\beta \in \left[0,1\right].$$

In fact, we have the same formulation as in the case of Proposition 2,but now in a quite different context:

Proposition 2.

Condition $b>a$ is a sufficient one for the components $E_1$ and $E_2$ of the complex dynamic modulus E, defined by (12), to satisfy

$$E_1\left(\omega \right)=E_1\left(-\omega \right),E_2\left(\omega \right)=-E_2\left(-\omega \right),\omega \in \mathbb{R},E_2\left(\omega \right)\ge 0,\omega \in {\mathbb{R}}_{+},$$

$${\int }_0^{\infty }\frac{1}{\omega }\frac{E_2\left(\omega \right)}{{(1+{\omega }^2)}^{\frac{m}{2}}}d\omega ={\int }_0^{\infty }\frac{(b-a)\rho cos\theta }{({(\rho cos\theta )}^2+{(\rho sin\theta +\omega )}^2){(1+{\omega }^2)}^{\frac{m}{2}}}d\omega <\infty ,$$

for some $m>0$.

Proof. 

Applying the Fourier transform to (10) and using (7), we obtain

$$\mathcal{F}\left({}^C{D}_{\left(K_1\right)}^{t}f\left(t\right)\right)\left(\omega \right)=\frac{{\left(i\omega \right)}^{1+\alpha -\beta }}{{\left(i\omega \right)}^{\alpha }+1},\omega \in \mathbb{R},$$

so that, instead of (13), we have

$$E\left(\omega \right)=\frac{1+b\frac{{\left(i\omega \right)}^{1+\alpha -\beta }}{{\left(i\omega \right)}^{\alpha }+1}}{1+a\frac{{\left(i\omega \right)}^{1+\alpha -\beta }}{{\left(i\omega \right)}^{\alpha }+1}},\omega \in \mathbb{R}.$$

This implies that

$$E\left(\omega \right)=1+\left(b-a\right)\frac{{\left(i\omega \right)}^{1+\alpha -\beta }}{{\left(i\omega \right)}^{\alpha }+1+a{\left(i\omega \right)}^{1+\alpha -\beta }},\omega \in \mathbb{R}.$$

(16)

It will be clear from the proof that is to follow that the converse assumption $a>b$ does not imply the claims of the proposition. We will consider the case when $\omega >0$ since the case of $\omega <0$ is quite similar.

Let $\omega >0$ and $z=$${\left(i\omega \right)}^{1+\alpha -\beta }=Rexp\left(i\theta \right)$. Then, we have

$$i\omega =R^{\frac{1}{1+\alpha -\beta }}\left[cos\frac{\theta }{1+\alpha -\beta }+isin\frac{\theta }{1+\alpha -\beta }\right],$$

which implies that

$$cos\frac{\theta }{1+\alpha -\beta }=0 \mathrm{i.e.},\frac{\theta }{1+\alpha -\beta }=\frac{\pi }{2}\mathrm{and}sin\frac{\theta }{1+\alpha -\beta }>0,\mathrm{since}\omega >0.$$

Thus, $1+\alpha -\beta <1$ implies that

$$z={\left(i\omega \right)}^{1+\alpha -\beta }=Rexp\left(i(1+\alpha -\beta )\pi /2\right)$$

lies in the first quadrant. Rewrite (16) in the form

$$E\left(\omega \right)=1+\left(b-a\right)\frac{z}{z^{\frac{\alpha }{1+\alpha -\beta }}+1+az}.$$

Since

$$z^{\frac{\alpha }{1+\alpha -\beta }}=R^{\frac{\alpha }{1+\alpha -\beta }}(cos\frac{\alpha \pi }{2}+isin\frac{\alpha \pi }{2}),$$

we again have that $z^{\frac{\alpha }{1+\alpha -\beta }}$ belongs to the first quadrant.

Then, using simple geometric arguments (sum of three terms, each having an argument in $(0,\frac{\pi }{2})$), we conclude that

$$arg\left(\frac{z}{z^{\frac{\alpha }{1+\alpha -\beta }}+1+az}\right)={\theta }_E\in (0,\frac{\pi }{2}).$$

In the case of $\omega <0$, we obtain ${\theta }_E\in (-\pi /2,0)$. Rewriting E as

$$E\left(\omega \right)=1+(b-a)\rho e^{i{\theta }_E},$$

(17)

we obtain $E_1\left(\omega \right)\ge 0,E_2\left(\omega \right)\ge 0,\omega >0$.

Therefore, looking at the real and imaginary part of E in (17), we easily conclude that

$$E_1\left(\omega \right)=E_1\left(-\omega \right),E_2\left(\omega \right)=-E_2\left(-\omega \right),\omega \in \mathbb{R}.$$

For the last part, we note that $E_2\left(0\right)=0$, and it is differentiable in a neighborhood of zero so that $\frac{E_2\left(\omega \right)}{\omega }\ne 0$. □

2.3. Dissipation Inequality

Our main theorem of this section is to follow. Since the condition for cases H and Z are the same ($b>a$), the next theorem holds in both cases.

For simplicity, we assume that $\sigma$ and $\epsilon$ belong to $C^1[0,\infty )$.

Theorem 1.

With the quoted assumptions on σ and ε, and the constitutive Equation (13), condition $b>a$ implies the dissipation inequality (10).

Proof. 

The basic idea is the Bochner–Schwartz theorem for non-negative measures; see [33], p. 331.

We will rewrite (10) in another form.

$$\begin{array}{ccc}\hfill D\left(x\right)& =& {\int }_0^T\sigma (x,t)\frac{\partial \epsilon \left(x,t\right)}{\partial t}\left(t,x\right)dt={\int }_0^T\left({\mathcal{F}}^{-1}E\right)(x,t)\ast \epsilon (x,t)\frac{\partial \epsilon }{\partial t}\left(t,x\right)dt\hfill \\ & =& {\int }_0^T{\mathcal{F}}^{-1}\left(\frac{1}{i\omega }E\left(\omega \right)\right)(x,t)\ast \frac{\partial \epsilon }{\partial t}(x,t)\frac{\partial \epsilon }{\partial t}\left(t,x\right)dt,x\in \mathbb{R}.\hfill \end{array}$$

Let

$$F_1\left(\omega \right)=E_2\left(\omega \right)/\left(i\omega \right),and{F}_2\left(\omega \right)=E_1\left(\omega \right)/\left(i\omega \right),\omega \in \mathbb{R}.$$

Now, $F_1$ is even and $F_1\ge 0,\omega >0$, whereas $F_2$ is odd. Since $F_2$ is odd, as in our previous paper [19], we can show that

$${\int }_0^T{\mathcal{F}}^{-1}\left(F_2\left(\omega \right)\right)(x,t)\ast \frac{\partial \epsilon }{\partial t}(x,t)\frac{\partial \epsilon }{\partial t}\left(t,x\right)dt=0.$$

Thus, we have to prove that the non-negativity of $F_1\left(\omega \right)\ge 0,\omega \in \left(-\infty ,\infty \right),$ implies

$${\int }_0^T{\mathcal{F}}^{-1}\left(F_1\left(\omega \right)\right)(x,t)\ast \frac{\partial \epsilon }{\partial t}(x,t)\frac{\partial \epsilon }{\partial t}\left(t,x\right)dt\ge 0.$$

This will complete the proof of (11).

Let $\theta \in \mathcal{S}\left(\mathbb{R}\right)$. Then, by Theorem IX.10 in [33], the positivity of $F_1$ implies that ${\mathcal{F}}^{-1}\left(F_1\right)$ is positive definite; that is, for every $\theta \in \mathcal{S}\left(\mathbb{R}\right)$ and $supp\theta \in [0,T]$,

$${\int }_0^T{\int }_0^T({\mathcal{F}}^{-1}\left(F_1\right)(t-\tau )\theta \left(\tau \right)\theta \left(t\right)d\tau dt$$

$$=2{\int }_0^T{\int }_0^t{\mathcal{F}}^{-1}\left(F_1\left(\xi \right)\right)\left(t-\tau \right)\theta \left(\tau \right)\theta \left(t\right)d\tau dt\ge 0,$$

where we also use the assumption that $F_1$ is even.

Finally, since any function in $C[0,T]$ is a pointwise limit of a real-valued sequence ${\theta }_k,k\in \mathbb{N}$ in $\mathcal{S}\left(\mathbb{R}\right)$, let ${\theta }_k\to {\epsilon }^{\left(1\right)}$ pointwisely on open set $(0,T)$ as $k⟼\infty$. For such ${\theta }_k$,

$$\begin{array}{ccc}& & \underset{k\to \infty }{lim}{\int }_0^T{\int }_0^T{\mathcal{F}}^{-1}\left(F_1\left(\xi \right)\right)\left(t-\tau \right){\theta }_k\left(\tau \right){\theta }_k\left(t\right)d\tau dt\hfill \\ & =& {\int }_0^T{\int }_0^T{\mathcal{F}}^{-1}\left(F_1\left(\xi \right)\right)\left(t-\tau \right){\epsilon }^{\left(1\right)}\left(\tau \right){\epsilon }^{\left(1\right)}\left(t\right)d\tau dt\ge 0.\hfill \end{array}$$

The last expression is simply $D\ge 0$, i.e., (11) holds. □

Remark 1.

Conditions stated in Theorem 1 are obtained by a different approach in relation to the one used by Bagley and Torvik; see [34,35]. Their approach is based on the assumption that a sinusoidal stress, imposed on a viscoelastic body, after a transition period, implies that the strain has the same form but with a phase shift. The energy dissipation condition has to be satisfied during a deformation process starting from a virginal state and it is not required for this deformation to constitute a cycle. The approach used here is also used in [18,31,36].

3. Stress Relaxation and Wave Equation

We will show how condition $b>a$, which follows from the dissipation inequality, implies the solvability and the unicity of a solution for a real model, which will be described below. As we already noted, in the sequel, we consider only case Z.

The next lemma is needed.

Lemma 1.

(a) Let $s\in \mathbb{C},\Re s>0$. Then,

$$\Re s\sqrt{\frac{1+a\frac{s}{{(s+\mu )}^{1-\alpha }}}{1+b\frac{s}{{(s+\mu )}^{1-\alpha }}}}\ge 0.$$

(b) If $\Re s>c>0$, then there exists $c_0>0$ such that

$$\Re s\sqrt{\frac{1+a\frac{s}{{(s+\mu )}^{1-\alpha }}}{1+b\frac{s}{{(s+\mu )}^{1-\alpha }}}}\ge c_0$$

Proof. 

(a) The proof will be given in several simple steps. Let $\Re s>0,\Im s>0$. Then, using simple geometry, one can conclude that

$$z=\frac{s}{{(s+\mu )}^{1-\alpha }}=\rho e^{i\beta }\mathrm{satisfies}\arg \mathrm{z}=\beta \in (0,\pi /2).$$

The same reasoning shows that, for $\Re s>0$ and $\Im s<0$, there holds

$$\arg z=\beta \in (-\pi /2,0).$$

Thus, for $\Re s>0$, in both cases, we have $\Re z>0$.

Now, consider

$$Z=\frac{1+az}{1+bz}=r{e}^{i\gamma }\mathrm{with}\Re z>0.$$

Again, we conclude by elementary observation that $\gamma \in (0,\pi /2)$ implies

$$\arg Z=\gamma \in (-\pi /2,0),$$

since $\arg (1+bz)\ge arg(1+az)$. In addition, in a similar way, we conclude that $\gamma \in (-\pi /2,0)$ implies

$$\arg Z=\gamma \in (0,\pi /2),$$

because $arg(1+bz)\ge arg(1+az)$.

Now, we will determine the location of points of $(t,r)=A\left(z\right)$, where

$$A\left(z\right)=t+ir=\sqrt{\frac{1+az}{1+bz}}=\sqrt{\frac{a}{b}+\frac{b-a}{b(1+bz)}},\Re z>0.$$

We will decompose mapping A in several simple mappings, having in mind that $\Re z>0$.

1. $z\mapsto w_1=b(1+bz)$ transforms $\Re z>0$ into $\Re w_1>b$;

2. $w_1\mapsto w_2=\frac{b-a}{w_1},$ transforms $\Re w_1>b$ into the interior of the circle

$$(b-a)\frac{1}{b}[\left({(t-1/2)}^2\right)+r^2=1/4],$$

that is, to the set of points $w_2=t+ir$ belonging to the interior of the circle ${(t-1/2)}^2+r^2={\left(\sqrt{(b-a)/2b}\right)}^2$;

3. $w_2\mapsto w_3=\frac{a}{b}+w_2$ transforms the interior of the circle to the translated one,

$$C=\{(t,r):\frac{a}{b}+[{(t-1/2)}^2+r^2={\left(\sqrt{(b-a)/2b}\right)}^2\}.$$

(18)

With this analysis, we conclude that $\Re z>0$ is transformed by the mapping A into the interior of C given by (18).

We conclude:

If $\Re s>0$ and $\Im s>0$, then we know that z has the real part $>0$ and the imaginary part $>0$. Thus, $A\left(z\right)$ has a positive real part and negative imaginary part. With this, we conclude that $arg\left(sA\right(z\left)\right)$ must belong to $(-\pi /2,\pi /2)$. The same arguments show that the assumption $\Re s>0,\Im s<0$ imply that $arg\left(sA\right(z\left)\right)$ must belong to $(-\pi /2,\pi /2)$. This proves assertion (a).

(b) The first quadrant in the z-plane goes by the mapping A into the lower half of C in (18), whereas the fourth quadrant in the z-plane goes to the upper half of C. Thus, by the first part, we have that A maps points $\Re z>0$ into points $(t,r)$ of the complex plane so that $(t,r)$ belongs to the set of points

$$\{(t,r):\sqrt{\mathrm{points}\mathrm{in}\mathrm{the}\mathrm{interior}\mathrm{of}}C\}.$$

Now, the first part of the proof implies that, after the multiplication of s and $A\left(z\right)$, one must have $\Re \left(sA\left(z\right)\right)\ge \sqrt{a/b}$. This proves assertion (b). □

Remark 2.

Previous proof shows that, for $a>b$, the circle $(b-a)\frac{1}{b}[\left({(t-1/2)}^2\right)+r^2=1/4]$ lies in the left side of $\Re z=0$ and intersects both the left and right half of the complex plane. This contradicts a). Thus, the conclusion $b>a$ that we have obtained from the dissipation inequality appears as the essential one for the existence and the uniqueness of a solution of a wave equation in Theorem 2.

3.1. Stress Relaxation for Z in Case $\epsilon (x,t)=H\left(t\right){\mathbf{1}}_x$

We treat the creep problem. Applying the Laplace transform to (10), we obtain

$$\tilde{\sigma }\left(x,s\right)=\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\tilde{\epsilon }\left(x,s\right),x\in \mathbb{R},\Re s>0.$$

For the stress relaxation test, we take $\epsilon \left(x,t\right)=H\left(t\right){\mathbf{1}}_x$, where ${\mathbf{1}}_x$ is the characteristic function of $\mathbb{R}$ (${\mathbf{1}}_x=1,x\in \mathbb{R}$). Since $\widetilde{H\left(s\right)}=\frac{1}{s},\Re s>0$, we have

$$\tilde{\sigma }\left(x,s\right)=\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\tilde{\epsilon }\left(x,s\right)=\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\frac{1}{s},x\in \mathbb{R},\Re s>0.$$

(19)

Proposition 3.

Let $\epsilon \left(x,t\right)=H\left(t\right){\mathbf{1}}_x$ be the strain in (10). Then, the stress σ has the form

$$\sigma \left(x,t\right)=\delta \left(x\right)(\frac{b}{a}H\left(t\right)+g\left(t\right)),x\in \mathbb{R},t\in [0,\infty ),$$

where, for any $x_0>0$, which means that the integral does not depend on $x_0>0$,

$$g\left(t\right)=\frac{1}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }exp\left(ts\right)\frac{a-b}{bs}\frac{1}{1+a\frac{s}{{(s+\mu )}^{1-\alpha }}}ds,t\in {\mathbb{R}}_{+},$$

is a continuous function.

Proof. 

By Lemma 1, part (b), the integrand does not have zeros in the right half of $\mathbb{C}$. Moreover,

$$\frac{a-b}{bs}\frac{1}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\mathrm{behaves}C{s}^{-1-\alpha },\left|s\right|\to \infty ,$$

(20)

so that, using the Cauchy formula, the complex integral below does not depend on the choice of $x_0>0$. Thus,

$$\sigma \left(x,t\right)=\frac{1}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }exp\left(ts\right)\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\frac{1}{s}ds,x\in \mathbb{R},t\in {\mathbb{R}}_{+},$$

where $x_0>0$. Thus,

$$\sigma \left(x,t\right)=\frac{\delta \left(x\right)}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }exp\left(ts\right)\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\frac{1}{s}ds,$$

$$=\frac{\delta \left(x\right)}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }exp\left(ts\right)(\frac{b}{as}+\frac{a-b}{as}\frac{1}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}})ds$$

$$=\delta \left(x\right)(\frac{b}{a}H\left(t\right)+g\left(t\right)),$$

where

$$g\left(t\right)=\frac{1}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }exp\left(ts\right)\left(\frac{a-b}{bs}\frac{1}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\right)ds,t>0,$$

is a continuous function because of (20). □

The result of the numerical inversion of (19) is shown in Figure 1.

In Figure 1, $x\in \mathbb{R}$. The stress relaxation curve for the case of a “classical” fractional Zener model is given in [16], p. 64.

3.2. Wave Equation for Z

We present a wave equation for case Z. Note that the waves in the Zener model of the viscoelastic body with the Riemann–Liouville and Caputo derivative have been studied in many papers. We refer to the review articles [19,37]. The one-dimensional equation of motion, constitutive equation, and geometrical conditions in the dimensionless form are

$$\begin{array}{cc}\hfill & \frac{\partial }{\partial x}\sigma (x,t)+f(x,t)=\frac{{\partial }^2}{\partial t^2}u(x,t),\hfill \end{array}$$

(21)

$$\begin{array}{cc}\hfill a{}^C{D}_{\left(K\right)}^t\sigma \left(t,x\right)+\sigma \left(t,x\right)& =b{}^C{D}_{\left(K\right)}^t\epsilon \left(t,x\right)+\epsilon \left(t,x\right),x\in \mathbb{R},t\in {\mathbb{R}}_{+}\hfill \end{array}$$

(22)

$$\begin{array}{cc}\hfill & \epsilon (x,t)=\frac{\partial }{\partial x}u(x,t),x\in \mathbb{R},t\in {\mathbb{R}}_{+}\hfill \end{array}$$

(23)

where $f\left(t,x\right)$ denotes the body force. We associate to (21)–(23) the following initial conditions:

$$u\left(0,x\right)=u_0\left(x\right),\frac{\partial u\left(0,x\right)}{\partial t}=v_0\left(x\right),x\in \mathbb{R}.$$

(24)

where $u_0$ and $v_0$ are functions with properties that will be discussed in the main theorem of this section. The use of the Laplace transform, (1)${}_3$ and (9), give, for $\Re s>0$,

$$\begin{array}{ccc}\hfill \mathcal{L}\left({}^C{D}_{\left(K\right)}^{t}f\right)\left(s\right)& =& \widetilde{\left({}^C{D}_{\left(K\right)}^{t}f\right)}\left(s\right)={\tilde{K}}_2\left(s\right)\widetilde{f^{\left(1\right)}}\left(s\right)\hfill \\ & =& \frac{1}{{\left(s+\mu \right)}^{1-\alpha }}\left[s\tilde{f}\left(s\right)-f\left(0\right)\right]\hfill \\ & =& \frac{s}{{\left(s+\mu \right)}^{1-\alpha }}\tilde{f}\left(s\right)-\frac{f\left(0\right)}{{\left(s+\mu \right)}^{1-\alpha }}.\hfill \end{array}$$

Note that, for $\mu =0$,we recover the result presented in [38], p. 98. Applying the Laplace transform to (22), we obtain

$$\tilde{\sigma }\left(x,s\right)=\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\tilde{\epsilon }\left(x,s\right)=\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\widetilde{\frac{\partial u}{\partial x}}\left(x,s\right),x\in \mathbb{R},\Re s>0,$$

where we used $\sigma \left(x,0\right)=\epsilon \left(x,0\right)=0$. Let

$$L\left(t\right)={\mathcal{L}}^{-1}\left(\frac{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\right)\left(t\right),t>0,$$

so that

$$\sigma \left(t,x\right)=L\left(t\right){\ast }_t\frac{\partial u\left(t,x\right)}{\partial x}={\int }_0^{t}L\left(t-\tau \right)\frac{\partial u\left(\tau ,x\right)}{\partial x}d\tau ,t\in \left(0,\infty \right),x\in \mathbb{R}.$$

Then, (21)–(24) become

$$\frac{{\partial }^2}{\partial t^2}u(x,t)-L\left(t\right){\ast }_t\frac{{\partial }^{2}u\left(t,x\right)}{\partial x^2}=f\left(x,t\right)+u_0\left(x\right){\delta }^{\left(1\right)}\left(t\right)+v_0\left(x\right)\delta \left(t\right).$$

(25)

In order not to go into cumbersome detail, we will say that a continuous function $m\left(x\right),x\in \mathbb{R},$ rapidly decreases enough as $\left|x\right|\to \infty$ if it is convolvable with a given function b; that is, $x\mapsto {\int }_{\mathbb{R}}h(x-t)b\left(t\right)dt,x\in \mathbb{R},$ is a continuous function. Thus, we assume:

(i)

$f(·,t)$ and the first three derivatives $\frac{d}{dx}f(·,t)$, $\frac{d^2}{d{x}^2}f(·,t)$, and $\frac{d^3}{d{x}^3}f(·,t),t\ge 0,$ rapidly decrease enough with respect to $\left|x\right|\to \infty$ so that their convolutions in x with a bounded continuous function are continuous on $\mathbb{R}\times (0,\infty )$.

(ii)

$u_0$ and $v_0$ and their derivatives up to the third order rapidly decrease enough with respect to $\left|x\right|\to \infty$ so that their convolutions in x with a bounded continuous function are continuous on $\mathbb{R}\times (0,\infty )$.

Theorem 2.

Assume that f, $u_0$, and $v_0$ satisfy assumptions $\left(i\right)$ and $\left(ii\right)$ given above. Equations (21)–(24) have a unique solution given by

$$u(x,t)=P\left(x,t\right){\ast }_{t,x}B\left(x,t\right),x\in \mathbb{R},t>0,$$

where

$$P\left(x,t\right)={\mathcal{L}}^{-1}\left(\frac{1}{2}s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}exp\left(-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)\right)\left(x,t\right),$$

is a solution of the equation

$$\frac{{\partial }^2}{\partial t^2}u(x,t)-L\left(t\right){\ast }_t\frac{{\partial }^{2}u\left(t,x\right)}{\partial x^2}=\delta \left(x\right)\delta \left(t\right),$$

and

$$B\left(x,t\right)={\mathcal{L}}^{-1}\left[\frac{\tilde{f}(x,s)+u_0\left(x\right)s+v_0\left(x\right)}{s^2}\right]\left(x,t\right)$$

(26)

$$=f(x,t)t_{+}+u_{0}H\left(t\right)+v_0{t}_{+},x\in \mathbb{R},t>0.$$

$P(x,t)$ is a distribution so that it is a smooth bounded function out of a set of points $(x,t)\in (\mathbb{R}\setminus (-\epsilon ,\epsilon \left)\right)\times [0,\infty )$ for any $\epsilon >0$. The singularity is the point $x=0$ so that, in a neighbourhood of this point, $x\in (-a,a)$, one has that P is the third derivative with respect to x of a bounded continuous function over to $(-a,a)\times [0,\infty )$.

The solution $u(x,t),(x,t)\in \mathbb{R}\times [0,\infty )$ is a continuous function.

Proof. 

The Laplace transform ${\mathcal{L}}_{t\mapsto s}$ applied to (25), with $x\in \mathbb{R}$ and $\mathfrak{R}s>0$, gives

$$\frac{{\partial }^2}{\partial x^2}\tilde{u}(x,s)-s^2\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\tilde{u}(x,s)=-\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\tilde{B}(x,s).$$

(27)

Applying the Fourier transform ${\mathcal{F}}_{x\mapsto \xi }$ to (27), we obtain

$$-{\xi }^2\widehat{\tilde{u}}(\xi ,s)-s^2\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\widehat{\tilde{u}}(\xi ,s)=-\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}\widehat{\tilde{B}}(\xi ,s).$$

This implies that

$$\widehat{\tilde{u}}(\xi ,s)=\frac{s^2\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}{{\xi }^2+s^2\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}·\frac{\widehat{B}(\xi ,s)}{s^2},\xi \in \mathbb{R},\Re s>0.$$

(28)

The next step is to apply the inverse Fourier transform. We note that

$${\mathcal{F}}^{-1}\left(\frac{a}{{\xi }^2+a}\right)\left(x\right)=\frac{1}{2}\sqrt{a}{e}^{-\left|x\right|\sqrt{a}},x\in \mathbb{R},$$

if $\mathfrak{R}a>0$. In order to use this result, we need to check that $\mathfrak{R}s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}$ is positive if $\mathfrak{R}s>0$. This follows from Lemma 1, part (a). Thus, for $x\ne 0$, the function $P(x,t),t\ge 0$ is smooth. Thus, P is smooth over $(\mathbb{R}\setminus \{0\left\}\right)\times [0,\infty )$.

Next, (28) becomes

$$\widehat{\tilde{u}}(x,s)=s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}{\mathcal{F}}^{-1}\left(\frac{s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}}{{\xi }^2+s^2\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\frac{\widehat{\widehat{B}}(\xi ,s)}{s^2}\right)\left(x,s\right),$$

(29)

for $\xi \in \mathbb{R},\mathfrak{R}s>0$, and, thus,

$$\begin{array}{ccc}\hfill \tilde{u}\left(x,s\right)& =& s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}exp\left(-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)\hfill \\ & & {\ast }_x\frac{\tilde{f}(x,s)+u_0\left(x\right)s+v_0\left(x\right)}{s^2}.\hfill \end{array}$$

(30)

Put

$$\begin{array}{ccc}\hfill P\left(x,t\right)& =& {\mathcal{L}}^{-1}\left(\frac{1}{2}\left(s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)\right.\hfill \\ & & \times \left.exp\left(-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)\left(x,t\right)\right),\begin{array}{c}x\in \mathbb{R},t>0.\end{array}\hfill \end{array}$$

Lemma 1 implies that we can take any $x_0>0$ in the integrals below since we will perform the use of te Cauchy formula in the analysis of these integrals. Thus, for $x\in \mathbb{R}$, $t\in {\mathbb{R}}_{+}$,

$$P(x,t)=\frac{1}{4\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}exp\left(ts-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)ds,$$

(31)

we conclude that the point $x=0$ is singular. Thus, the integral (31) must be understood in the sense of dual pairing (it does not exist in the classical sense). We have to enter into the space of the distribution and rewrite P into the form

$$P\left(x,t\right)=$$

$$-\frac{1}{2}\frac{d^3}{d{x}^3}{\mathcal{L}}^{-1}\left({\left(s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)}^{-2}\times exp(-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}})\right)\left(x,t\right),$$

so that

$$\begin{array}{ccc}\hfill K(x,t)& =& \frac{1}{2}{\left(s\sqrt{\frac{1+a\frac{s}{{s+\mu )}^{1-\alpha }}}{1+b\frac{s}{{(s+\mu )}^{1-\alpha }}}}\right)}^{-2}exp\left(-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right),\hfill \\ \hfill x& \in & \mathbb{R},t\ge 0\hfill \end{array}$$

is integrable in a neighbourhood of $x=0,t\ge 0$. We have that

$$P(0,t)={\mathcal{L}}^{-1}\left(\frac{1}{2}{\left(s\sqrt{\frac{1+a\frac{s}{{s+\mu )}^{1-\alpha }}}{1+b\frac{s}{{(s+\mu )}^{1-\alpha }}}}\right)}^{-2}\right),t\ge 0,$$

is a bounded continuous function.

Let

$$P_0(x,t)=\frac{1}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }K(x,s)ds,x\in \mathbb{R},t>0.$$

Regarding the behavior of the integrand, we conclude that it is a bounded continuous function on $\mathbb{R}\times [0,\infty )$. Since (26) holds for B, we obtain

$$u(x,t)=-P_0(x,t){\ast }_{x,t}(\frac{{\partial }^3}{\partial x^3}f(x,t){\ast }_t{t}_{+}+u_0^{\left(3\right)}\left(x\right){\ast }_{t}H\left(t\right)+v_0^{\left(3\right)}{\ast }_t{t}_{+}),x\in \mathbb{R},t>0.$$

(32)

This, with the assumptions on $f,u_0$ and $v_0$, implies that $u(x,t),x\in \mathbb{R}\times [0,\infty ),$ is a continuous function. This completes the proof of the theorem. □

Corollary 1.

If $u_0=0$ and f and $v_0$ and their derivatives up to order 5 satisfy assumptions $\left(i\right)$ and $\left(ii\right)$, then the solution u is a classical one.

Proof. 

The assumptions imply that the second derivative of $u(x,t)$ is continuous in x.We note that the member with $u_0$ disappears and that the second derivative of $t_{+}$ is $\delta \left(t\right)$. Thus, in the convolution of B with respect to t with $P_0$, one obtains for the solution that the second derivative $\frac{{\partial }^2}{\partial t^2}u$ is continuous for $t\in [0,\infty )$. □

Corollary 2.

If $f,u_0$ and $v_0$ and their derivatives up to order 6 satisfy assumptions $\left(i\right)$ and $\left(ii\right)$, then the solution u is a classical one.

Proof. 

Note that $u_0\ne 0$ gives an additional member on the right hand side of $P_0(x,t){\ast }_{x,t}\left(u_0^{\left(3\right)}\left(x\right){\delta }^{\prime }\left(t\right)\right)$. Thus, in order to improve the properties of $P_0$ with respect to t, one should put in the definition of P, $P(x,t)=\frac{{\partial }^4}{\partial x^4}....$. Then, u given by (32) will have a new form: a convolution of a new $P_0$ and the fourth derivatives of $f,u_0$ and $v_0$ With this, one has that $\frac{{\partial }^2}{\partial t^2}P(x,t)$ is also continuous with respect to $x\in \mathbb{R}$ and $t\ge 0$. This completes the proof. □

Example

As a specific example, we consider Equations (21)–(24) with

$$f\left(x,t\right)=0,v_0\left(x\right)=0,u_0\left(x\right)=\delta \left(x\right).$$

The solution is given by (30). It has the form

$$u\left(x,t\right)=-\frac{1}{2\pi i}{\int }_{x_0-i\infty }^{x_0+i\infty }\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}exp\left(-\left|x\right|s\sqrt{\frac{1+a\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}{1+b\frac{s}{{\left(s+\mu \right)}^{1-\alpha }}}}\right)ds$$

(33)

It has a strong singularity at $x=0$ for any $t>0$, so, in our numerical experiment, we consider it for $x\ne 0$. In Figure 2, we show (33) for several values of parameters. The results presented in Figure 2 show that, for selected values of $a=0.2,b=1,\alpha =0.3$, by increasing the parameter $\mu$, the speed of the propagation of the maximum decreases and the amplitude of the maximum decreases.

Figure 3 shows the solution when $\mu$ increases further. It is seen that the speed of the propagation of the maximum decreases further, whereas the amplitude of the maximum increases. This interesting property needs to be studied further.

4. Conclusions

We considered a general form of a fractional derivation and fractional integral, suggested by several leading experts in the field. Two Sonin pairs, called cases $\mathbf{H}$ of Hanyga and $\mathbf{Z}$ of Zacher, were analyzed through the dissipation inequality for the Zener-type constitutive equation. Our approach for the proof of the dissipation inequality was based on the Bochner–Schwartz theorem as in our recent papers [18,36]. The limitation of our approach is that it can only be applied to linear constitutive equations.

Our framework includes the space of generalized functions. This enables us to use the strong results of the Schwartz theory. However, our results are closely connected with the formulations of classical analysis.

In the case $\mathbf{Z}$, we present results related to a stress relaxation and a wave propagation for a Zener-type viscoelastic body for which the proposed thermodynamical restriction for coefficients in the constitutive equation guarantee the existence and uniqueness.

The analysis of a solution for a wave propagation is the main part of our analysis. It appears that some estimates based on the condition $b>a$, followed from the dissipation inequality, are necessary for the use of the Fourier and Laplace transforms and their inverses. Additional assumptions on initial data $u_0,v_0$ and the perturbation f imply the existence of the classical solution given in the integral form.

It is shown in examples that a parameter $\mu$ appearing in the definition of the generalized derivative in a Sonin pair proposed by Zacher allows for a better fit of experimental results. The numerical results show the properties of a solution out of a neighbourhood of $x=0$.

Finally, we comment the possible extension of the present work. The introduction of stochastic terms in the constitutive equation for the case Z when $\mu \ne 0$ (cf. (8)) seems possible. However, the first step in such an analysis is the restriction on the coefficients that are obtained in this paper.