#### Abstract

The linearly damped oscillator equation is considered with the damping term generalized to a Caputo fractional derivative. The order of the derivative being considered is . At the lower end the equation represents an undamped oscillator and at the upper end the ordinary linearly damped oscillator equation is recovered. A solution is found analytically, and a comparison with the ordinary linearly damped oscillator is made. It is found that there are nine distinct cases as opposed to the usual three for the ordinary equation (damped, over-damped, and critically damped). For three of these cases it is shown that the frequency of oscillation actually increases with increasing damping order before eventually falling to the limiting value given by the ordinary damped oscillator equation. For the other six cases the behavior is as expected, the frequency of oscillation decreases with increasing order of the derivative (damping term).

#### 1. Introduction

In this paper the linearly damped oscillator equation is considered with the damping term replaced by a fractional derivative [1] whose order, , will be restricted to, ,

Burov and Barkai [2, 3] examined such an equation in connection with critical behavior. They were able to determine a solution in terms of generalized Mittag-Leffler functions. Nonlinear fractional oscillators have been studied numerically by Zaslavsky et al. [4]. He was primarily interested in chaotic behavior. It is hoped that a careful study of the analytic solution to the linear fractionally damped equation will help shed light on properties of the nonlinear equation and be of use for direct applications of fractionally damped oscillations (see, e.g., [5, 6]).

In this paper the Caputo formulation of the fractional derivative will be used. The Caputo derivative is preferred over the Riemann-Liouville derivative for physical reasons. Consider the Laplace transform of the two formulations of the fractional derivative for

The constant term arising from the Laplace transform of the Caputo derivative is merely the initial value of the function. For the Riemann-Liouville derivative this is not the case. The constant term arising from the Laplace transform currently has no simple physical interpretation. Hence the Caputo fractional derivative seems to be more useful for modeling physical systems.

If the order of the fractional damping term is allowed to become (outside the range of values considered in this paper) the equation is usually referred to as the Bagley-Torvik equation (see, e.g., [1, 7, 8]). The solution of this equation exhibits damped oscillatory behavior similar to what we expect to find for the equation studied in this paper. The Bagley-Torvik equation was originally derived to study the motion of a rigid plate in a Newtonian fluid [7].

The analytic solution to the fractionally damped equation is found by means of Laplace transform. For the sake of clarity, and for pointing out some unique difficulties with the factional equation, a comparison with the Laplace transform method as applied to the nonfractional case is made. It is found that there are nine distinct cases for the fractionally damped equation as opposed to the usual three cases for the Nonfractional equation. In six of the nine cases the results are as expected; increasing the order of the factional derivative increases the effects of the damping (i.e., the frequency of the damping slows as the order of the derivative increases). However, in three cases, the frequency of the damping actually increases as the order of the fractional derivative increases until a peak value is reached after which the frequency falls to its Nonfractional limit. The physical reason for this increase in the oscillation frequency is not yet clear.

#### 2. The Nonfractional Case

Before a solution to the linear fractionally damped oscillator equation is constructed it will be useful to review the Laplace transform method of solution for the linearly damped oscillator equation

The constants and are taken to be real and positive. is the damping force per unit mass. is the restoring force per unit mass. In both cases (fractional and Nonfractional) the following initial conditions will be used:

Transforming (2.1) together with the initial conditions (2.2) gives the following:

or

Equation (2.4) can be inverted using tables, however, to shed light on a problem that will happen later, (2.4) will be inverted via the complex inversion integral. The exponents on the variable in both terms are whole numbers, hence there will not be a branch cut in the contour integral and the Bromwich contour can be used

Recall that the Bromwich contour begins at goes vertically up to (where is chosen so that all poles will lie to the left of the vertical contour line and thus all poles will be captured within the contour) and then travels in a half circle (to the left, counter clockwise) back to . For this problem there is no contribution from the contour integral. The only contribution comes from the residue. The residue is generated from the roots of the following quadratic equation

There are three different cases.

; 2 unequal real roots that are negative

; 2 repeated real roots that are negative

; 2 complex roots whose real parts are negative

See Figure 1 for a graphical representation of the location of the roots.

Note that in cases one and three the poles will be of order one and in case two the pole will be of order two. Note also that if there were no damping the poles would be on the imaginary axis at .

Computing the residue for case one gives

or

As and are both negative this solution will decay exponentially. This is usually referred to as the over-damped case.

Case three is computed the same way as case one. Now the poles are complex so the exponential function can be expressed using sine and cosine with an over-all exponential damping factor

where and Notice that the presence of damping causes the effective angular frequency, , to be smaller than the un-damped angular frequency; that is, the oscillations go slower, as one might expect if there were damping to impede the motion. This is usually referred to as the under-damped case. By comparison, cases one and two could be viewed as having a zero frequency or an infinite period.

For case two the pole is of order 2, and the residue is given by

Recall that , this allows the denominator to be factored. The limit then becomes

or

This is usually called the critically damped case. Graphs of sample solutions to these three cases can be found in any introductory book on differential equations.

#### 3. The Fractional Case

Now consider (1.1). In this case, has units of time raised to the power . Hence the over all units of the second term remain the same as in (2.1). There are two cases to consider, and . We will consider the former in this paper. The Laplace transform of (1.1) is

or

If (3.2) is inverted using a contour integral a branch cut is needed on the negative real axis due to the fractional exponents on the complex variable . Hence, a Hankel contour will be used. This contour starts at and goes vertically up to (where is again chosen so that all poles will lie to the left of the vertical contour line) then travels in a quarter circle arc (to the left) to just above the negative real axis (i.e., ). The contour then has a cut that goes into the origin (following the negative real axis), around the origin in a clockwise sense (to just below the negative real axis) and then back out to . The contour then has another quarter circle arc to .

Now the question is, where are the poles? This is a somewhat more involved question than in the standard linearly damped model. To find the poles the following equation needs to be solved:

Which, for an arbitrary , is not a trivial problem. To determine if there are solutions, and if so how many, let then (3.3) breaks into 2 equations, a real and an imaginary part

Could there be a solution on the positive real axis? No, in this case and the first equation of (3.4) would be the sum of three positive nonzero terms, which would never be zero. Could there be a solution on the negative real axis? No, in this case and the second term of the second equation of (3.4) would never be zero. Using similar arguments we can show that there are no solutions on the positive or negative imaginary axes, recall . It can also be shown that no solutions are in the right half plane (both terms of the second equation would always be positive). If there are solutions, they should be in pairs, complex conjugates, with and . To attempt to find a solution first solve the second equation of (3.4) for *r* and substitute this into the first equation (only look for positive values first)

The reader may be worried about the negative sign and the fractional exponent in (3.5), however, for the restricted angular range being considered, , is always negative. So, the argument of the root will always be positive.

Given values for it would appear to be impossible to solve (3.5) for . Equation (3.5) can be simplified to a more aesthetically pleasing form

Now it needs to be seen if there is a value that will satisfy (3.6). For this equation to be true needs to be positive. This will only happen on the restricted domain . Now the question becomes, on this restricted domain can we pick a value that will make the left-hand side of (3.6) as large or as small as we wish? Thus ensuring that no matter what the values of , , and we are given we can always find a value that will satisfy (3.6). Consider the two limits

Since the left-hand side of (3.6) is continuous in and we have the two limits above, (3.7), it is guaranteed that there will be at least one solution to (3.6) and hence there will be at least two poles for the residue calculation. If we can show that the left hand side of (3.6) decreases monotonically in over the restricted domain then we know that there will be only one solution to (3.6), and thus only two poles in the residue calculation. To show that the left hand side of (3.6) decreases monotonically in we need to show that the derivative of the left hand side of (3.6) with respect to is always negative, that is,

Computing the derivative, doing some algebra, and throwing away over all factors that are always positive we have

On the restricted domain , , and . This reduces (3.9) to

Equation (3.10) can now be factored into a perfect square and prove the assertion made in(3.8)

Hence the left hand side of (3.6) will decrease monotonically on the restricted domain with the upper bound being and the lower bound being 0. To summarize, it has just been shown that there is always one solution, with a positive angle, to (3.6) and this solution must be such that . Consequently there will be two poles for the residue calculation and they will be complex conjugates of each other. Notice that for the fractionally damped equation repeated roots are not possible. Repeated roots can only happen when the order of the derivative becomes one. See Figure 2 for a graphical representation of the location of the roots.

Now the question of the poles that has been settled the solution to (1.1) can be generated. Denote the two poles as

where and are determined from the *r* and values that satisfy (3.6) in the usual way, and . Note that is negative, the two solutions are in the second and third quadrants and is the complex conjugate of . Note also that plays the role of from the Nonfractional case. The poles are of order one and the residue is given by

where has been replaced by . After some algebra this can be reduced to

where denote to .

For the contour integral the only contributions come from the paths along the negative real axis

The solution to (1.1) is then (3.15) subtracted from (3.14). This may look overly complicated but the solution does have the general form of

Notice that the decay function, (3.15), goes to zero if goes to zero or one; that is, if (1.1) goes to its Nonfractional limits the decay function goes away, as expected. The damping factor is similar to the damping factor for the Nonfractional case, . Notice that since the poles, for the residue calculation, have nonzero imaginary and nonzero real parts we will not have the same three distinct cases as we did for the Nonfractional case (critically damped, over-damped, and under-damped).

#### 4. The Oscillation Frequency

Consider the frequency of the oscillation component of the solution, . One question we might ask is: how does the frequency change as we change the order of the fractional damping? When is set to zero we have an un-damped oscillator with frequency

When is set to one we have the three cases given in Section 2: over-damped, critically-damped, and under-damped. So, the frequency may be zero or nonzero, that is,

For there will always be a nonzero frequency. Note that .

In the Nonfractional case increasing causes the frequency of oscillation to become smaller, monotonically, until the critical cases are reached and the oscillation period becomes infinite (these are the critical and over-damped cases). In the fractional case the frequency of oscillation, , now depends on the order of the derivative, , as well as and . To try to determine how depends on these three parameters consider *s* to be a function of , on , implicitly defined by

for fixed values of and (both being positive). Let us restrict our attention to the upper half plane for *s*. As such *s* will be one to one on . Due to the fractional exponent causing a branch cut on the negative real axis *s* will not be one to one at .

Now the question arises, does fall monotonically with respect to ? To get the answer to this consider the derivative of (4.3) with respect to and isolate (remember, and are being held fixed)

The imaginary part of this equation is

Specifically, consider this equation at

This gives three initial slopes for the rate of change of with respect to .

The frequency initially increases with increasing damping order.The frequency initially is not changing with increasing damping order. The frequency initially decreases with increasing damping order.This is not entirely what might have been expected. In the first case the oscillation frequency actually increases before falling. Hence there will be some values of for which the fractional damping will actually cause the oscillations to go faster than the un-damped oscillator (the damping will still cause the amplitude to decrease). Each of the above three cases can become any of the three Nonfractionally damped cases by letting ((2.7), (2.8), and (2.9)). Hence, there are nine cases for the linear fractionally damped oscillator.

There are some graphs of solutions to the imaginary part of (4.3) (the oscillation frequency) for various values of , , and . In all three graphs the oscillation frequency is on the vertical axis and the order of the derivative is on the horizontal axis. The three graphs for each case correspond to what would be under-damped, critically damped, and over-damped for a damped oscillator with whole order derivatives.

Figure 3 is a representative graph of case one, . The green graph is for , the black graph is for , and the red graph is for .

Figure 4 is a representative graph from case two, , a flat start. The red graph is for , the green graph is for , and the black graph is for .

Figure 5 is a representative graph for case three, a decreasing start. The green graph is for , the black graph is for , and the red graph is for .

#### 5. Conclusion

In this paper the linear fractionally damped oscillator equation was solved analytically. It was found that the solution is very similar to the Nonfractional case (decayed oscillations but with the inclusion of an additional decay function). It was found that there are nine distinct cases, as opposed to the usual three for the ordinary damped oscillator. An unexpected result was that for three of the cases the oscillation frequency actually increases with increasing order of derivative of the damping term (till a peak value is reached, then the frequency decreases as expected). The physical reason for this increase in oscillation frequency is not yet clear.