1 Introduction

There has been a great deal of interest in the fractional calculus for the last decades. The reasons of this interest are the findings which came to hand when some of the researchers modeled some dynamic systems making use of fractional operators [1–7]. One of the most interesting results among them is the one obtained by Bagley and Tovirk [8, 9] who studied the viscoelasticity structures and the behavior of materials using fractional derivatives. The equation used by these two scientists is called the Bagley–Torvik equation and has the form [3]

$$ \lambda _{2} u''(t)+ \lambda _{1} D^{\nu }u(t)+\lambda _{0} u(t)=f(t),\quad t>0, $$
(2)

where \(\lambda _{i}\), \(i=0,1,2\) are constants, \(\lambda _{2}\neq 0\), f is a real-valued function and \(D^{\nu }\) is the fractional derivative \(\nu =\frac{1}{2}\) or \(\nu =\frac{3}{2}\). Because of the applications of the Bagley–Torvik equation, many researches tackled the problem of finding the analytic and numerical solution of this equation [3, 10–21].

On the other side, discrete fractional calculus has also attracted the attention of many researchers. This type of calculus dealing with the sums and differences with non-integer quantities has also many applications in variety of fields [22–25]. Motivated by the above, we intended to find the analytic and the numerical solution of a certain discrete version of the Bagley–Torvik equation in this article. The considered equation contains the h-difference which as \(h\rightarrow 0\) gives the classical derivative. To the best of our knowledge no one [26, 27] has discussed this before.

The article is arranged as follows.

The second and third sections discuss the theory needed to handle the equation under consideration. The fourth section proposes the solutions of the discrete Bagley–Torvik equation. The fifth section present numerical solutions of some special cases of the mentioned equation. The sixth section is devoted to the conclusion.

2 Preliminaries

In this section, some basic definitions and results which will be used further are presented.

Definition 2.1

Let \(u(t),\ t\in [0,\infty )\), be a real- or complex-valued function and \(h>0\) be a fixed shift value. Then the forward difference operator on \(h\mathbb{Z}\) is defined as

$$ \Delta _{h} u(t) =\frac{ u(t+h) - u(t)}{h}, $$
(3)

and the backward difference operator on \(h\mathbb{Z}\) is defined as

$$ \nabla _{h} u(t) =\frac{ u(t) - u(t-h)}{h}. $$
(4)

For \(h=1\), this gives \(\Delta u(t) = u(t+1) - u(t)\) and \(\nabla u(t) = u(t) - u(t-1)\), respectively.

The forward jumping operator on the time scale \(h\mathbb{Z}\) is \(\sigma _{h}(t)=t+h\) and the backward jumping operator is \(\rho _{h}(t)=t-h\). For \(a,b\in R\) and \(h>0\), we use the notation \(\mathbb{N}_{a,h}=\{a, a+h, a+2h,\ldots,\}\) and \(_{b,h}\mathbb{N}=\{b,b-h,b-2h,\ldots\}\).

Definition 2.2

For \(h>0\) and \(\mu \in R\), the increasing h-polynomial factorial function is defined as

$$ t_{h}^{\bar{\mu }}=h^{\mu } \frac{\Gamma (\frac{t}{h}+\mu )}{\Gamma (\frac{t}{h})}, $$
(5)

where \(t_{h}^{[{0}]}=1\), Γ is the Euler gamma function and \(\frac{t}{h}+\mu ,\frac{t}{h},\notin \{0,-1,-2,-3,\ldots\}\), as the division at a pole yields zero.

If μ is a positive integer, then

$$ t_{h}^{\bar{\mu }} = t(t+h) (t+2h)\cdots \bigl(t+(\mu -1)h \bigr). $$
(6)

Remark 2.3

Applying the nabla operator on (5), then

$$ \nabla _{h} t_{h}^{\bar{\mu }} = \mu \ t_{h}^{\overline{\mu -1}}. $$
(7)

Proposition 2.4

([28] The relation between nabla \(h-RL\) fractional difference and h-Caputo fractional difference)

  1. (i)

    \({}_{a}^{C}\nabla _{h}^{\nu }u(t)={}_{a}\nabla _{h}^{\nu }u(t)- \frac{(t-a)_{h}^{\overline{-\nu }}u(a)}{\Gamma (1-\nu )}\).

  2. (ii)

    \({}_{h}^{C}\nabla _{b}^{\nu }u(t)={}_{h}\nabla _{b}^{\nu }u(t)- \frac{(b-t)_{h}^{\overline{-\nu }}u(b)}{\Gamma (1-\nu )}\).

Definition 2.5

([29])

Assume that u is defined on \(\mathbb{N}_{a,h}\). Then the h-discrete Laplace transform of u is defined by

$$ \mathcal{N}_{a,h}\bigl[u(t)\bigr](s)= \int _{a}^{\infty }(1-hs)^{ \frac{t-a-h}{h}}u(t)\nabla _{h}t=h\sum_{t=a/h+1}^{\infty }(1-hs)^{t-a/h-1}u(ht). $$
(8)

When \(a=0\), this gives

$$ \mathcal{N}_{0,h}\bigl[u(t)\bigr](s)= \int _{0}^{\infty }(1-hs)^{ \frac{t-h}{h}}u(t)\nabla _{h}t=h\sum_{t=1}^{\infty }(1-hs)^{t-1}u(ht). $$
(9)

The following results are the h-discrete Laplace transform for the Caputo fractional difference and also for the integer difference operator.

Lemma 2.6

For the function \(u(t)\) defined on \(\mathbb{N}_{a,h}\) and \(n-1<\nu \leq n\), then

$$ \mathcal{N}_{a,h}\bigl[{}_{a}^{C} \nabla _{h}^{\nu }u(t)\bigr](s)=s^{\nu } \mathcal{N}_{a,h}\bigl[u(t)\bigr](s)-\sum_{k=0}^{n-1} s^{\nu -1-k} \nabla _{h}^{k}u(a). $$
(10)

For the positive integer n,

$$ \mathcal{N}_{a,h}\bigl[{}_{a}\nabla _{h}^{n} u(t)\bigr](s)=s^{n} \mathcal{N}_{a,h}\bigl[u(t)\bigr](s)- \sum _{k=0}^{n-1} s^{n-1-k}\nabla _{h}^{k}u(a). $$
(11)

Definition 2.7

(Nabla h-discrete Mittag-Leffler)

For \(\lambda \in R\), \(|\lambda |<1\) and \(\theta ,\beta ,\rho ,t\in C\) with \(\operatorname{Re}(\theta )>0\), the nabla h-discrete Mittag Leffler functions are defined by

$$ {}_{h}{E}_{\overline{\theta ,\beta }}^{\rho }( \lambda ,t)=\sum_{k=0}^{ \infty }\lambda ^{k} \frac{t_{h}^{\overline{k\theta +\beta -1}}(\rho )_{k}}{\Gamma {(\theta k+\beta )}k!}. $$
(12)

For \(h=\beta =\rho =1\), one can write

$$ {E}_{\overline{\theta }}(\lambda ,t)=\sum _{k=0}^{\infty } \lambda ^{k} \frac{t^{\overline{k\theta }}}{\Gamma {(\theta k+1)}}, $$
(13)

where \((\rho )_{k}=\rho (\rho +1)\cdots (\rho +k-1)\) and \((1)_{k}=k!\).

Lemma 2.8

(Finite inverse principle law)

Let \(t>0\), \(h>0\) and m be a positive integer. Then, for the equation \(\nabla _{h} v(t)=u(t)\), \(v(t)=\nabla _{h}^{-1}u(t)\) obeys the finite inverse principle law as

$$ v(t)-v(t-mh)=h\sum_{r=0}^{m-1}u(t-rh). $$
(14)

Proof

Take \(\nabla _{h} v(t)=u(t)\), now applying Eq. (4), then

$$\begin{aligned}& \frac{v(t)-v(t-h)}{h}=u(t)\quad \Rightarrow\quad v(t)-v(t-h)=hu(t) \\& v(t)=hu(t)+v(t-h). \end{aligned}$$
(15)

Replace \(v(t)\) by \(v(t-h)\) in (15) and resubstitute in (15), then

$$ v(t)=hu(t)+hu(t-h)+v(t-2h)\quad \Rightarrow\quad v(t)-v(t-2h)=h\sum _{r=0}^{1}u(t-rh). $$
(16)

Continuing like this gives (14). □

Lemma 2.9

Let \(t\in R\), \(a,h>0\). Then

$$ \nabla _{h}^{-1}(1-hs)^{\frac{t}{h}-1}=- \frac{(1-hs)^{\frac{t}{h}}}{s}. $$
(17)

Proof

Take \(u(t)=(1-hs)^{\frac{t}{h}-1}\) in (4), which gives

$$ \nabla _{h} (1-hs)^{\frac{t}{h}-1} = \frac{(1-hs)^{\frac{t}{h}-1}-(1-hs)^{\frac{t-h}{h}-1}}{h}= \frac{(1-hs)^{\frac{t}{h}-1}(1-(1-hs)^{-1})}{h}. $$
(18)

Now, the proof of (17) follows by taking \(\nabla _{h}^{-1}\) on both sides. □

Corollary 2.10

Let \(t\in (-\infty ,\infty )\) and \(h,s>0\), then

$$ -\frac{(1-hs)^{\frac{t}{h}}}{s}+\frac{(1-hs)^{\frac{t-mh}{h}}}{s}=h \sum _{r=0}^{m-1}(1-hs)^{\frac{t-rh}{h}-1}. $$
(19)

Proof

The proof follows by using the finite inverse principle law in (17). □

Example 2.11

For the particular values of \(h=2\), \(s=5\), \(t=6\) and \(m=100\), Eq. (19) is verified by MATLAB and it turns out that its numerical value is 145.8.

Lemma 2.12

Let \(h>0\) and u, w be real-valued bounded functions. Then

$$ \indent \nabla _{h}^{-1}\bigl(u(t)w(t) \bigr)=u(t)\nabla _{h}^{-1}w(t)-\nabla _{h}^{-1} \bigl(\nabla _{h}^{-1}w(t-h) \nabla _{h} u(t)\bigr). $$
(20)

Proof

Applying the nabla operator on the function \(u(t)v(t)\) gives

$$\begin{aligned}& \begin{aligned} \nabla _{h}\bigl[u(t)v(t)\bigr]&=\frac{u(t)v(t)-u(t-h)v(t-h)}{h} \\ & =\frac{u(t)v(t)-u(t)v(t-h)+u(t)v(t-h)-u(t-h)v(t-h)}{h}, \end{aligned} \\& \nabla _{h}\bigl[u(t)v(t)\bigr] = u(t)\nabla _{h} v(t) +v(t-h)\nabla _{h} u(t), \\& u(t)\nabla _{h} v(t) = \nabla _{h}\bigl[u(t)v(t) \bigr] - v(t-h)\nabla _{h} u(t). \end{aligned}$$

Now considering \(w(t)=\nabla _{h} v(t)\) and taking \(\nabla _{h}^{-1}\) on both sides gives (20). □

3 Generalized nabla discrete h-Laplace transform and its convolution

Following the time scale calculus, one gave the following definition for the nabla discrete Laplace transform on \(\mathbb{N}_{a,h}\) modifying Lemma 2.6 using the closed form(inverse difference operator).

Definition 3.1

Assume that \(u(t)\) is defined on \(\mathbb{N}_{a,h}\). Then the generalized nabla discrete Laplace transform of u is defined by

$$\begin{aligned} \mathcal{N}_{a,h}\bigl\{ u(t)\bigr\} (s) =& \int _{a}^{\infty }{}_{h} \tilde{e}_{\ominus s}^{\rho }(t,a)u(t)\nabla _{h} t \\ =& \int _{a}^{\infty }\frac{{}_{h}\tilde{e}_{\ominus s}(t,a)}{1-hs}u(t)\nabla _{h} t = \int _{a}^{\infty } ({1-hs} )^{\frac{t-a-h}{h}}u(t) \nabla _{h} t. \end{aligned}$$
(21)

Using the closed and summation form solution, the above equation can be written as

$$ \mathcal{N}_{a,h}\bigl\{ u(t)\bigr\} (s)={}_{a}\nabla _{h}^{-1} u(t) ({1-hs} )^{\frac{t-a-h}{h}}|_{a}^{\infty } =h\sum _{i=a/h+1}^{ \infty }u(ih) ({1-hs} )^{i-\frac{a}{h}-1}. $$
(22)

Remark 3.2

(i) In the case \(a=0\),

$$ \mathcal{N}_{0,h}\bigl\{ u(t)\bigr\} (s)=\nabla _{h}^{-1} u(t) ({1-hs} )^{\frac{t-h}{h}}|_{0}^{\infty } =h\sum_{i=1}^{ \infty }u(ih) ({1-hs} )^{i-1}. $$
(23)

(ii) In the special case \(h=1\),

$$ \mathcal{N}_{a}\bigl\{ u(t)\bigr\} (s)={}_{a}\nabla ^{-1} u(t) ({1-s} )^{{t-a-1}} |_{a}^{\infty } =\sum_{i=a+1}^{\infty }u(ih) ({1-hs} )^{i-{a}-1}. $$
(24)

Theorem 3.3

For \(t\in \mathbb{N}_{a,h}\), \(h,\mu >0\) and \(s\neq 0\),

$$ \nabla _{h}^{-1} \bigl[t_{h}^{\overline{\mu }}(1-hs)^{\frac{t}{h}-1} \bigr]=-\sum_{i=1}^{\mu +1} \frac{\mu ^{(i-1)}t_{h}^{\overline{\mu +1-i}}(1-hs)^{\frac{t}{h}}}{s^{i}}. $$
(25)

Proof

Taking \(u(t)=t\) and \(w(t)=(1-hs)^{\frac{t}{h}-1}\) in (20), using (7) and (17), then

$$ \nabla _{h}^{-1}t(1-hs)^{\frac{t}{h}-1}=- \frac{t(1-hs)^{\frac{t}{h}}}{s}-\frac{(1-hs)^{\frac{t}{h}}}{s^{2}}. $$
(26)

Again taking \(u(t)=t_{h}^{\overline{2}}\) and \(w(t)=(1-hs)^{\frac{t}{h}-1}\) in (20) gives

$$ \nabla _{h}^{-1}t_{h}^{\overline{2}}(1-hs)^{\frac{t}{h}-1}=- \frac{t_{h}^{\overline{2}}(1-hs)^{\frac{t}{h}}}{s}+\frac{2}{s} \nabla _{h}^{-1} \bigl[t{(1-hs)^{\frac{t}{h}}} \bigr]. $$
(27)

Now applying (26) and simplifying give

$$ \nabla _{h}^{-1}t_{h}^{\overline{2}}(1-hs)^{\frac{t}{h}-1}=- \frac{t_{h}^{\overline{2}}(1-hs)^{\frac{t}{h}}}{s}- \frac{2t_{h}^{\overline{1}}(1-hs)^{\frac{t}{h}}}{s^{2}}- \frac{2(1-hs)^{\frac{t}{h}}}{s^{3}}, $$

which can be rewritten as

$$ \nabla _{h}^{-1} \bigl[t_{h}^{\overline{2}}(1-hs)^{\frac{t}{h}-1} \bigr]=-\sum_{i=1}^{3} \frac{2^{(i-1)}t_{h}^{\overline{3-i}}(1-hs)^{\frac{t}{h}}}{s^{i}}. $$

By proceeding the above process up to μ times one finds (25). □

Lemma 3.4

Let \(\mu ,h>0\) and \(s\neq 0\), then

$$ \mathcal{N}_{h} \bigl[t_{h}^{\overline{\mu -1}} \bigr]= \frac{(\mu -1)!}{s^{\mu }}. $$
(28)

Proof

The proof follows by applying the limits 0 to ∞ in (25) and using (24). □

Example 3.5

By equating (24) and (28),

$$ \mathcal{N}_{h} \bigl[t_{h}^{\overline{\mu -1}} \bigr]=h \sum_{i=0}^{\infty }(ih)_{h}^{\overline{\mu -1}}(1-hs)^{i-1}= \frac{(\mu -1)!}{s^{\mu }}, $$

which is verified by MATLAB for the particular values of \(h=2\), \(s=1/3\) and \(\mu =3\) being numerically equal as regards both closed and summation form solution as 54.

Remark 3.6

For the fraction ν, one can write (28) as

$$ \mathcal{N}_{h} \bigl[t_{h}^{\overline{\nu -1}} \bigr]= \frac{\Gamma (\nu )}{s^{\nu }}. $$
(29)

Definition 3.7

([29])

Let \(s\in R\), \(0<\nu <1\) and \(u,v:\mathbb{N}_{a,h}\rightarrow R\) be a function. The nabla h-discrete convolution of u with v is defined by

$$ (u*v) (t)= \int _{a}^{t}u(s)v\bigl(t-\rho (s)+a\bigr) \nabla _{h} s=h\sum_{k=a/h+1}^{t/h}u(kh)v \bigl(t-\rho (kh)+a\bigr). $$
(30)

Theorem 3.8

([29] The h-convolution theorem)

For any \(\nu \in R\slash \{\ldots,-2,-1,0\}\), \(s\in R\) and u, v defined on \(\mathbb{N}_{a,h}\), we have

$$ \mathcal{N}_{a,h}\bigl[(u*v) (t)\bigr](s)=\mathcal{N}_{a,h} \bigl[u(t)\bigr](s)\times \mathcal{N}_{a,h}\bigl[v(t)\bigr](s). $$
(31)

Lemma 3.9

For \(\lambda \in R\), \(|\lambda |<1\) and \(\theta ,\beta ,\nu ,t\in C\) with \(\operatorname{Re}(\theta )>0\),

$$ \mathcal{N}_{h} \bigl[{}_{h}{E}_{\overline{\theta ,\beta }}^{\nu }( \lambda ,t) \bigr]= \frac{s^{\alpha \nu -\beta }}{(s^{\alpha }-\lambda )^{\nu }}. $$
(32)

Proof

Applying the Laplace transform in (12) and using (29) give

$$\begin{aligned}& \begin{aligned} \mathcal{N}_{h} \bigl[{}_{h}{E}_{\overline{\theta ,\beta }}^{\nu }( \lambda ,t) \bigr]& =\sum_{k=0}^{\infty } \frac{\lambda ^{k}(\nu )_{k}}{\Gamma {(\theta k+\beta )}k!} \mathcal{N}_{a,h} \bigl[t_{h}^{\overline{k\theta +\beta -1}} \bigr] =\sum_{k=0}^{\infty } \frac{\lambda ^{k}}{s^{\theta k+\beta }}\frac{(\nu )_{k}}{k!} \\ & = s^{-\beta }\sum_{k=0}^{\infty } \biggl( \frac{\lambda }{s^{\theta }} \biggr)^{k}\frac{(\nu )_{k}}{k!} = s^{-\beta } \biggl(\frac{s^{\theta }}{s^{\theta }-\lambda } \biggr)^{\nu }, \end{aligned} \\& \mathcal{N}_{h} \bigl[{}_{h}{E}_{\overline{\theta ,\beta }}^{\nu }( \lambda ,t) \bigr] = \frac{s^{\alpha \nu -\beta }}{(s^{\alpha }-\lambda )^{\nu }}. \end{aligned}$$

 □

4 Solution of discrete Bagley–Torvik equation

In this section, we find the analytic solution of the discrete fractional Bagley–Torvik equation given as

$$ \nabla _{h}^{2} u(t)+A{}^{C} \nabla _{h}^{\nu }u(t)+Bu(t)=f(t),\quad t>0, $$
(33)

where \(0<\nu <1\) or \(1<\nu <2\), subject to \(u(0)=a\) and \(\nabla _{h} u(0)=b\), with a and b being real numbers.

Here, one can solve the above fractional equation in the two cases of the particular values of \(\nu =\frac{1}{2}\) and \(\frac{3}{2}\) with numerical analysis.

4.1 Case 1: \(\nu =\frac{1}{2}\)

Here, the researchers solve the discrete fractional Bagley–Torvik equation for \(\nu =\frac{1}{2}\) by employing the Laplace transform.

Theorem 4.1

The discrete fractional Bagley–Torvik equation

$$ \nabla _{h}^{2} u(t)+A{}^{C} \nabla _{h}^{1/2} u(t)+Bu(t)=f(t),\quad t>0 $$
(34)

where \(0<\nu <1\), subject to \(u(0)=a\) and \(\nabla _{h} u(0)=b\), has the solution

$$\begin{aligned} u(t) =&a\sum_{n=0}^{\infty }(-1)^{n} A^{n} {}_{h} E_{ \overline{2,\frac{3n}{2}}}^{n+1}(-B,t)+b \sum_{n=0}^{\infty }(-1)^{n} A^{n} {}_{h} E_{\overline{2,\frac{3n}{2}+2}}^{n+1}(-B,t) \\ &{}+ \sum_{n=0}^{\infty }(-1)^{n} A^{n} \sum_{\tau =1}^{t}{}_{h} E_{\overline{2,\frac{3n}{2}+2}}^{n+1}(-B,\tau )f(t-\rho {\tau }),\quad \biggl\vert \frac{B}{s^{2}} \biggr\vert < 1. \end{aligned}$$
(35)

Proof

By applying the Laplace transform on (34) and using Lemma 2.6, one finds

$$ s^{2}\mathcal{N}_{h}\bigl[u(t) \bigr](s)+\sum_{k=0}^{1} s^{1-k}\nabla _{h}^{k} u(0) +As^{1/2}\mathcal{N}_{h}\bigl[u(t)\bigr](s)+B \mathcal{N}_{h}\bigl[u(t)\bigr](s)= \bar{f}(s), $$
(36)

simplifying and applying the initial conditions lead to

$$\begin{aligned} \begin{gathered} \bigl(s^{2}+As^{1/2}+B\bigr) \mathcal{N}_{h}\bigl[u(t)\bigr](s) =as+b+\overline{f(s)} \\ \mathcal{N}_{h}\bigl[u(t)\bigr](s)=\frac{1}{(s^{2}+As^{1/2}+B)} \bigl(as+b+ \overline{f(s)}\bigr). \end{gathered} \end{aligned}$$
(37)

Now,

$$ \begin{gathered} \frac{1}{(s^{2}+As^{1/2}+B)}=\frac{1}{s^{2}+B} \frac{1}{1+\frac{As^{1/2}}{s^{2}+B}}=\frac{1}{s^{2}+B}\sum_{n=0}^{\infty }(-1)^{n}\frac{A^{n} s^{n/2}}{(s^{2}+B)^{n}},\\ \frac{1}{(s^{2}+As^{1/2}+B)}= \sum_{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{n/2}}{(s^{2}+B)^{n+1}}, \quad \biggl\vert \frac{As^{1/2}}{s^{2}+B} \biggr\vert < 1. \end{gathered} $$
(38)

Now using (38) in (37),

$$\begin{aligned} \mathcal{N}_{h}\bigl[u(t)\bigr](s) =&a \sum _{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{n/2+1}}{(s^{2}+B)^{n+1}}+b\sum_{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{n/2}}{(s^{2}+B)^{n+1}} \\ &{}+\sum_{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{n/2}}{(s^{2}+B)^{n+1}}\overline{f(s)}. \end{aligned}$$
(39)

Now applying the inverse Laplace transform and using (32) give (35). □

4.2 Case 2: \(\nu =\frac{3}{2}\)

Theorem 4.2

The discrete fractional Bagley–Torvik equation

$$ \nabla _{h}^{2} u(t)+A{}^{C} \nabla _{h}^{3/2} u(t)+Bu(t)=f(t),\quad t>0, $$
(40)

where \(1<\nu <2\), subject to \(u(0)=a\) and \(\nabla _{h} u(0)=b\), has the solution

$$\begin{aligned} u(t) =&a\sum_{n=0}^{\infty }(-1)^{n} A^{n} {}_{h} E_{ \overline{2,\frac{n}{2}+1}}^{n+1}(-B,t)+a \sum_{n=0}^{\infty }(-1)^{n} A^{n+1} {}_{h} E_{\overline{2,\frac{n}{2}+\frac{3}{2}}}^{n+1}(-B,t) \\ &{}+b \sum_{n=0}^{\infty }(-1)^{n} A^{n} {}_{h} E_{ \overline{2,\frac{n}{2}+2}}^{n+1}(-B,t)+ \sum_{n=0}^{\infty }(-1)^{n} A^{n} \sum_{\tau =1}^{t}{}_{h} E_{\overline{2,\frac{n}{2}+2}}^{n+1}(-B, \tau )f(t-\rho {\tau }). \end{aligned}$$
(41)

Proof

By applying the Laplace transform on (40) and using Lemma 2.6,

$$\begin{aligned}& s^{2}\mathcal{N}_{h}\bigl[u(t) \bigr](s)-\sum_{k=0}^{1} s^{1-k}\nabla _{h}^{k} u(0)+A \Biggl[s^{1/2}\mathcal{N}_{h}\bigl[u(t)\bigr](s)-\sum _{k=0}^{1} s^{1/2-k} \nabla _{h}^{k} u(0)\Biggr] \\& \quad {}+B\mathcal{N}_{h} \bigl[u(t)\bigr](s)=\bar{f}(s), \end{aligned}$$
(42)

and simplifying and applying the initial conditions,

$$\begin{aligned} \begin{gathered} \bigl(s^{2}+As^{3/2}+B\bigr) \mathcal{N}_{h}\bigl[u(t)\bigr](s) =a\bigl(s+As^{1/2} \bigr)+b+ \overline{f(s)}, \\ \mathcal{N}_{h}\bigl[u(t)\bigr](s)=\frac{1}{(s^{2}+As^{3/2}+B)}\bigl[a \bigl(s+As^{1/2}\bigr)+b+ \overline{f(s)}\bigr]. \end{gathered} \end{aligned}$$
(43)

Now,

$$ \begin{gathered} \frac{1}{(s^{2}+As^{3/2}+B)}=\frac{1}{s^{2}+B} \frac{1}{1+\frac{As^{3/2}}{s^{2}+B}}=\frac{1}{s^{2}+B}\sum_{n=0}^{\infty }(-1)^{n}\frac{A^{n} s^{3n/2}}{(s^{2}+B)^{n}},\\ \frac{1}{(s^{2}+As^{1/2}+B)}= \sum_{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{3n/2}}{(s^{2}+B)^{n+1}}, \quad \biggl\vert \frac{As^{3/2}}{s^{2}+B} \biggr\vert < 1. \end{gathered} $$
(44)

Now using (44) in (43) gives

$$\begin{aligned} \mathcal{N}_{h}\bigl[u(t)\bigr](s) =&a \sum _{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{3n/2+1}}{(s^{2}+B)^{n+1}}+a \sum_{n=0}^{ \infty }(-1)^{n} \frac{A^{n} s^{(3n/2+1)/2}}{(s^{2}+B)^{n+1}} \\ &{}+b\sum_{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{3n/2}}{(s^{2}+B)^{n+1}}+\sum_{n=0}^{\infty }(-1)^{n} \frac{A^{n} s^{3n/2}}{(s^{2}+B)^{n+1}}\overline{f(s)}. \end{aligned}$$
(45)

Now applying the inverse Laplace transform and using (32) give (41). □

5 Results and discussion

In this section, we give numerical solutions and graphical illustrations for the considered Bagley–Torvik equation for particular values of some parameters.

For the particular values of \(a=b=A=B=h=0.1\) and \(\nu =1/2\), \(f(t)=1\), the solution (35) is graphically shown in Fig. 1.

Figure 1
figure 1

Solution of Bagley–Torvik equation \(f(t)=1\), \(\nu =1/2\)

Full size image

Again for the particular values of \(a=b=0.0001\), \(A=B=0.001\), \(h=0.1\) and \(\nu =1/2\), \(f(t)=t\), the solution (35) is graphically shown in Fig. 2.

Figure 2
figure 2

Solution of Bagley–Torvik equation \(f(t)=t\), \(\nu =1/2\)

Full size image

When the non-integer order is \(\nu =3/2\) and for the values \(a=0.01\), \(b=0.02\), \(A=0.1\), \(B=0.2\), \(h=0.15\) and \(f(t)=t^{2}+1\) the solution (35) graphically is shown in Fig. 3.

Figure 3
figure 3

Solution of Bagley–Torvik equation \(f(t)=t^{2}+1\), \(\nu =3/2\)

Full size image

Finally, for the values \(a=0.125\), \(b=0.15\), \(A=B=0.25\), \(h=0.05\) and \(f(t)=0\) the solution (35) is graphically shown in Fig. 4,

Figure 4
figure 4

Solution of Bagley–Torvik equation \(f(t)=0\), \(\nu =3/2\)

Full size image

6 Conclusion

In this article, the authors discussed a certain version of the discrete Bagley–Torvik equation involving a nabla h-fractional Caputo difference. The researchers obtained the analytical solutions favorably associated to the discrete Laplace transform and the discrete Mittag-Leffler functions. The researchers presented the numerical solutions for specific values of initial values, parameters and the right hand side of the equation. The nabla difference considered can be replaced by the delta difference operator. In this case, one should not think that the analytical solutions can be easily obtained. On the other hand, researchers may also replace the h-fractional Caputo difference by newly defined fractional differences involving non-singular kernels.