# Motion described by $m \frac{\mathrm{d}^2 x}{\mathrm{d}t^2}=-k\frac{\mathrm{d}^{\frac12 }x}{\mathrm{d}t^{\frac12}}$

What kind of motion would a (preferably dimensionless for simplicity) body do if the force acted on it was proportional to the semi-derivative of displacement, i.e.

$$m \frac{\mathrm{d}^2 x}{\mathrm{d}t^2}=-k\frac{\mathrm{d}^{\frac12}x}{\mathrm{d}t^{\frac12}} \, \, ?$$

It would be helpful if someone with a copy of Mathematica plotted this for various values of the constants.

• Might be relevant: physics.stackexchange.com/q/4005 – Kyle Kanos Oct 8 '14 at 18:56
• More on fractional derivatives. – Qmechanic Oct 8 '14 at 19:05
• Hi user148432. If you haven't already done so, please take a minute to read the definition of when to use the homework tag, and the Phys.SE policy for homework-like problems. – Qmechanic Oct 8 '14 at 19:13
• @user148432 you might find math.stackexchange.com a good place for detailed help. – user6972 Oct 8 '14 at 20:26
• It is good practice to state what you have done so far instead of letting answerers start from scratch. At its current state, your (interesting) question does unfortunately sound quite like a "hey, anyone please work this out for me"... – Tobias Kienzler Oct 15 '14 at 11:37

If $$D^n$$ denotes the $$n$$th derivative and $$D^{-n}$$ the $$n$$th integral, then we have that,

$$D^n f(t) = D^m[D^{-(m-n)}f(t)]$$

providing $$m \geq \lceil{n}\rceil$$. For our half derivative, we choose $$n=1/2$$, and $$m=2$$, in which case we have,

$$D^{1/2}f(t) = D^2[D^{-(3/2)}f(t)]$$

There is a general formula for the $$n$$th integral of a function, one of my favorite results of Cauchy:

$$f^{-(n)}(t) = \frac{1}{\Gamma(n)}\int_{0}^t (t-u)^{n-1}f(u) \, du$$

which is essentially a convolution $$f(t) \ast t^{n-1}$$. Applying it, we find,

$$D^{1/2}f(t) = \frac{d^2}{dt^2} \left[ \frac{2}{\sqrt{\pi}}\int_0^t (t-u)^{1/2}f(u) \, du\right]$$

Given the differential equation,

$$\frac{d^2 x(t)}{dt^2} = -\frac{k}{m} \frac{d^{1/2} x(t)}{dt^{1/2}}$$

we can substitute in our definition of $$D^{1/2}x(t)$$, and conclude,

$$x(t) = -\frac{2k}{m\sqrt{\pi}}\int_{0}^t (t-u)^{1/2}x(u)\, du + c_1t +c_2$$

for $$c_1,c_2 \in \mathbb{R}$$ which is an integral equation. If we can assume $$x(t)$$ is supported on $$[0,\infty)$$ only, then the integral is a convolution $$x(t) \ast \sqrt{t}$$ and taking the Laplace transform, we find,

$$X(s) = \left( 1+ \frac{k}{ms^{3/2}}\right)^{-1} \left( \frac{c_1}{s^2} + \frac{c_2}{s} \right) = \frac{m(c_1 + c_2 s)}{k\sqrt{s}+ms^2}$$

The solution $$x(t)$$ is then the inverse Laplace transform of $$X(s)$$. Formaly, this is given by,

$$x(t) = \frac{1}{2\pi i} \int_{\Gamma} e^{st} \frac{m(c_1 + c_2 s)}{k\sqrt{s}+ms^2} \, ds$$

where the contour $$\Gamma$$ is in the complex plane; it is a vertical line of infinite length with all poles of $$F(s)$$ to its left. In practice, we close the contour with an additional contour, ensure the second integral tends to zero (e.g. by the estimation lemma), and use the residue theorem.

The integrand, which we denote $$F(s)$$, has three poles located at $$s^3 = k^2/m^2$$, or equivalently,

$$s_1 = \omega^{4/3}_0, \quad s_2 = \frac{1}{2}(1+i\sqrt{3})\omega^{4/3}_0, \quad s_3 = \frac{1}{2}(i\sqrt{3}-1)\omega^{4/3}_0$$

as well as at $$s_0= 0$$, where we define $$\omega^2_0 := k/m$$. The vertical contour should begin after $$s_1$$ so all poles are to the left. However, doing so analytically is somewhat tedious. I chose to use a numerical method for the evaluation of inverse Laplace transforms due to H.E. Salzer which uses aquadrature formula. With Mathematica, I managed to reconstruct $$x(t)$$ partially:

in the simplified case when $$c_1 = c_2 = k/m = 1$$. It seems, by visual inspection, the solution resembles that of damped harmonic motion, such as when one introduces a damping term $$\gamma \dot{x}$$ in the equations of motion of a standard harmonic oscillator.

• Not trying to be rude or anything, I'm grateful you answered so I upvoted, by I am not trying to find a way to solve the differential equation analytically, so that's why I don't accept your answer. I posted this on PSE so I can get some intuition on the kind of motion the body will do under the influence of a force proportional to that fractional derivative. I'm looking forward to a solution other than hypergeometric but I fear it will just be an alternative representation of hypergeometric functions, but we never know what intuition a different solution might provide us. – user148432 Oct 8 '14 at 20:35
• @user148432: Yes, I understood your question, but I started fiddling with the equation, and found a cute result, so I thought I'd share it as it is partially relevant, and may help you with the problem in general. – JamalS Oct 8 '14 at 20:37
• @JamalS. Didn't you forget the initial condition while applying the Laplace transform ? The Laplace transform of $\ddot x$ is $s^2X(s)-s x(0)-\dot x(0)$. – Tom-Tom Oct 9 '14 at 8:32
• @Tom-Tom: No, I did not. I applied the Laplace transform to the equation $x(t) = x(t) \ast t^{1/2}$ which has no derivatives; it is a convolution equation. – JamalS Oct 9 '14 at 8:37
• @fibonatic: No, it's not the square root of a derivative, it is literally the application of a differentiation half times, it is like half an integral and half a derivative. – JamalS Oct 9 '14 at 14:19

I am no mathematician and am a little afraid that my answer is too simple to be true, but here goes:

I use Fourier transforms to define the fractional derivative. $x(\omega)$ is defined such that

$$x(t) = \int_{-\infty}^\infty \, \frac{\text{d}\omega}{2\pi} \text{e}^{i \omega t} \, x(\omega) \, .$$

Then any integer derivatives is

$$\frac{\text{d}^n}{\text{d} t^n} x(t) = \int_{-\infty}^\infty \, \frac{\text{d}\omega}{2\pi} \text{e}^{i \omega t} \, (i \omega )^n \, x(\omega) \, ,$$

so that

$$\left(\frac{\text{d}^n x}{\text{d} t^n}\right)(\omega) = (i \omega )^n \, x(\omega) \, .$$

Then this can be generalised to any real number $n$. In particular,

$$\frac{\text{d}^{1/2}}{\text{d} t^{1/2}} x(t) = \int_{-\infty}^\infty \, \frac{\text{d}\omega}{2\pi} \text{e}^{i \omega t} \, \sqrt{i \omega } \, x(\omega) \, .$$

Since your equation is linear it can be solved separately for each Fourier modes. In Fourier space it becomes

$$\left(- m \omega^2 + k \sqrt{i \omega}\right) x(\omega) = 0 \, .$$

This tells us that either $x(\omega) = 0$ or $\omega^2 = \frac{k}{m} \sqrt{i \omega}$. The first solution is trivial. It is equivalent to $x(t)= 0$. The second however has four solutions:

\begin{align} & \omega_1 = \left(\frac{k}{m}\right)^{2/3} \text{e}^{i\,\pi/6} \, ,\\ & \omega_2 = \left(\frac{k}{m}\right)^{2/3} \text{e}^{i\,\pi/2} \, , \\ & \omega_3 = \left(\frac{k}{m}\right)^{2/3} \text{e}^{i \, 5\pi/6} \, , \\ & \omega_4 = 0 \, . \end{align}

Then the most general solution to your problem is

$$x(t) = C_1 \text{e}^{i\omega_1 t} + C_2 \text{e}^{i\omega_2 t} + C_3 \text{e}^{i\omega_3 t} + C_4\, .$$

$C_i$ are integration constants. Explicitly one finds,

$$x(t) = C_1 \, \text{e}^{i \,t \, (k/m)^{2/3} \sqrt{3}/2} \, \text{e}^{-(k/m)^{2/3} t/2} + C_2 \, \text{e}^{- (k/m)^{2/3} t} + C_3 \, \text{e}^{-i \, t \, (k/m)^{2/3} \sqrt{3}/2} \, \text{e}^{-(k/m)^{2/3} t/2} + C_4 \, .$$

We find three solutions that decay exponentially and one constant. Two of the decaying solutions oscillate as well.

In order to make a plot, I set C_4 = 0 because this amounts to a simple vertical shift of $x(t)$, I choose $C_1 = C_2^* = (c_1 + i c_2)/2$ and $C_3 = C_3^*$ so that $x(t)$ is real, I rescale the time according to $\tau = t \, (k/m)^{2/3}/2$ and rescale $x(t)$ according to $y = x \, C_3$. Then my solution becomes,

$$y(\tau) = \text{e}^{-\tau} \left[ c_1 \, \cos(\sqrt{3} \, \tau) + c_2 \sin(\sqrt{3} \, \tau) \right] + \text{e}^{-2 \tau} \, .$$

Here is a plot of $y(\tau)$ for different values of $c_{1,2}$:

• Sorry for the multiple algebra mistakes. I checked it and re-checked. I hope that it's Ok now. Please tell me if you find something strange. – Steven Mathey Oct 8 '14 at 23:32
• How does this generalization of derivatives to fractional derivatives compare to the definition JamalS uses? – CuriousKev Oct 16 '14 at 1:44
• I think this was an interesting approach. But from comments in Tobias's answer, this seems to be a different generalization of derivatives, since the fractional derivatives would commute now. – CuriousKev Oct 16 '14 at 5:17

You can simply take the semi-derivative of your equation again, which yields

\begin{align} m\frac{d^2}{dt^2}\underbrace{\frac{d^{\tfrac12}x}{dt^{\tfrac12}}}_{=-\frac mk\frac{d^2x}{dt^2}} &= -k\frac{dx}{dt} \\\Rightarrow m^2\frac{d^4x}{dt^4} &= k^2\frac{dx}{dt} \tag{*} \end{align}

and then solve that ODE. But, similarly to squaring an algebraic equation to eliminate roots, you (probably) have to discard halve of the solutions to satisfy the original equation.

As hinted at by user121330's comment, I assumed commutativity of $D^2$ (where $D:=d/dt$) and $D^{\tfrac12}$ which is not granted in general. I suspect that basically correlates with my above statement that some solutions of $(*)$ must actually be discarded. I'll try and think about this some more, but for now please be aware that there may be flaws or even uselessness in this answer...

• Any thoughts on initial conditions? This seems to require more conditions than usual, but then again the half-derivative is not time-local. – Emilio Pisanty Oct 15 '14 at 12:02
• @EmilioPisanty I honestly don't have much experience with fractional derivatives, but my guess is that similarly to squaring an equation (or taking derivatives of normal ODEs) you get more solutions than are valid for the original equation, i.e. two initial conditions should suffice. But as you state, the non-locality may yield curious effects. However, JamalS' answer (and others) clearly suggest two initial conditions just like for a usual 2nd degree ODE are required. – Tobias Kienzler Oct 15 '14 at 12:06
• I get a quite similar answer to Steven when I solve that, but I don't know what to discard as I don't know what restrictions I really have. – user148432 Oct 15 '14 at 12:07
• $m D^{1/2} D^2 x \neq m D^2 D^{1/2} x$ I desperately wanted an easy solution like this, but that commutation is pretty non-trivial. – user121330 Oct 15 '14 at 16:37
• I think this and Steven's answer are probably wrong, but I still found them very helpful for thinking about the problem. Maybe just add a note about why this approach doesn't work and keep it? It is still useful information. +1 for an interesting idea and realizing the mistake. – CuriousKev Oct 16 '14 at 5:13

One way to try to solve the equation is transforming it in an ODE. Apply the fractional derivative $D^{1/2}$ again to the equation to find $$D^{1/2}[D^2x(t)]=D^{5/2}x(t)-C_1t^{-3/2}-C_2t^{-5/2}-C_3t^{-7/2},$$ and $$D^{1/2}[D^{1/2}x(t)]=Dx(t)-C_4t^{-3/2}$$ Hence we got $$D^{5/2}x(t)=-\frac{k}{m}Dx(t)+C_1t^{-3/2}+C_2t^{-5/2}+C_3t^{-7/2}$$ But, we also have that $$D^{5/2}x(t)=D^2[D^{1/2}x(t)]=-\frac{m}{k}D^4x(t)$$ Hence we get the following ODE: $$x^{(4)}(t)-\omega^3 x'(t)=C_1t^{-3/2}+C_2t^{-5/2}+C_3t^{-7/2},$$ where $\omega^3=\dfrac{k^2}{m^2}$ Let $x'(t)=v(t)$ to get $$v'''(t)-\omega^3 v(t)=C_1t^{-3/2}+C_2t^{-5/2}+C_3t^{-7/2}$$ The general solution of the homogeneous equation for $v(t)$ is $$v_0(t)=c_1 e^{\omega t}+e^{-\omega t/2}\left(c_2\cos\dfrac{\sqrt3}{2}\omega t+c_3\sin\dfrac{\sqrt3}{2}\omega t\right)$$ The general solution of the homogeneous equation $x(t)$ is then $$x_0(t)=c_1 e^{\omega t}+e^{-\omega t/2}\left(c_2\cos\dfrac{\sqrt3}{2}\omega t+c_3\sin\dfrac{\sqrt3}{2}\omega t\right)+c_4,$$ where the constants are different in the two last equations. So, the general solution is of the form $$x(t)=x_0(t)+x_p(t),$$ where $x_p(t)$ is the particular solution.

• How did you get the first two equations? – CuriousKev Oct 16 '14 at 5:21
• By the rules of composition of fractional derivatives. The constants $C_i$ depend on the initial values of the function and fractional derivatives. The rules can be seen here in section 3.5.3: vutbr.cz/www_base/zav_prace_soubor_verejne.php?file_id=10060 – Mateus Sampaio Oct 16 '14 at 11:57

So, what you have here is, as others have mentioned, a fractional order differential equation. Since others have provided graphs, there seems to be little point to adding one. Also, you seem more interested in the qualitative aspect than actual analytical solutions, as you've mentioned.

In essence, what you have here is some second order system with some element that behaves like something between a spring and a damper. There is some dissipative behavior that will come into play here, but also some oscillatory behavior, as you see from the other graphs.

A special result of this being a fractional order system is that that specific derivative is also non-local, meaning that the current dynamics of the system are not only dependent on the current state, but also past state. This is called historicity, or a system with memory.

There are many ways to represent that fractional order derivative mathematically, if you want to get into that. Most of the people here seem to have used some integral definition (whether explicitly, or in the case of the Laplace transform, implicitly), but you can also represent it as a series.

If you wanted a full on analytical solution, you might run into some non-elementary functions, but you can derive series representations from them using some of the methods that Odibat as developed as a result of his work on generalized Taylor series expansions.

If you have any more questions, let me know.

An attempt for a more explicit solution. Using the definition of the half-derivative given by JamalS, one can transform the differential equation using the Laplace transform and get $$s^2X(s)-sx(0)-x'(0)=-\gamma^3{\sqrt s}\;X(s)$$ where $\gamma=\sqrt[3]{\frac km}$ is a positive constant. Solving for $X$ gives $$X(s)=\frac{sx(0)+\dot x(0)}{s^2+\gamma^3{\sqrt s}}.$$ Let us call $G(s)=1/(s^{3}+\gamma^3)$. The Laplace inverse of $s^{-1/2}G(s^{1/2})$ is given by (Erdelyi, Table of integral transforms, chapter IV, equation 4.1.33) $$z(t)=\frac{1}{\sqrt{\pi t}}\int_0^\infty g(u)\mathrm e^{-u^2/4t}\mathrm du.$$ Therefore we have $x(t)=\dot x(0)z(t)+x(0)\dot z(t)$.

Computation of $z(t)$ Let us write $\mathrm j=\mathrm e^{\mathrm i 2\pi/3}$. The inverse $g(t)$ of $G(s)$ is a simple Laplace inversion and yields $$g(t)=\frac1{3\gamma^2}\left[\mathrm e^{-\gamma t}+\mathrm e^{-\mathrm j\gamma t}+\mathrm e^{-\mathrm j^2\gamma t}\right]$$ Let us now write $$\Phi(a,t)=\frac{1}{\sqrt{\pi t}}\int_0^\infty\mathrm e^{-a u}\mathrm e^{-u^2/(4 t)}\mathrm du=\mathrm e^{a^2 t}\mathrm{erfc}(a\sqrt t)=\varphi(a^2t).$$

EDIT

A remarkable property of $\varphi$ is $\varphi'(x)=\varphi(x)-\frac1{\sqrt{\pi x}}$. Reporting in the definition of $z$ we find $$z(t)=\frac1{3\gamma^2}\left[\varphi(\gamma^2 t)+\varphi(\mathrm j\gamma^2t)+\varphi(\mathrm j^2\gamma^2 t)\right].$$ $$\dot z(t)=\frac13\left[\varphi(\gamma^2 t)+\mathrm j\varphi(\mathrm j\gamma^2 t)+\mathrm j^2\varphi(\mathrm j^2\gamma^2t)\right]-\frac2{3\sqrt{\pi\gamma^2t}}$$

Mathematica made the following plots of $z(t)$ and $\dot z(t)$ for $\gamma=1$, showing that $z$ is quite severely damped.

• It's somewhat curious that you obtain something radically different than what JamalS and Steven Mathey obtained. – Kyle Kanos Oct 15 '14 at 14:07
• @KyleKanos. The function ensemble in which the equations are solved differ, this is why, the final solutions are different. JamalS had troubles with the integration, the method I propose also uses a difficult integration but as the functions are analytical in the half complex plane instead of the imaginary line, it is possible to perform it. We actually used the same method, only with different initial conditions. Steven Mathey's answer is different, because he solved the problem without boundary conditions. Also keep in mind that $z(t)$ is not the solution, it is useful to build the solution. – Tom-Tom Oct 16 '14 at 7:30
• @KyleKanos. I have gone back to Mathematica to plot $\dot z$ together with $z$ and I did not get the same plot. I checked carefully and I am quite sure the new version is the correct one. There must have been a typo in the last version. The functions actually oscillate in my solution too. – Tom-Tom Oct 16 '14 at 8:05