How to solve highly oscillating differential equation

[karlzr]

The equation looks like: ##x''(t)+b x'(t)+cx(t)+ d x^3(t)=0##. This is the motion of a particle in a potential ##cx^2/2+d x^4/4## with friction force ## b x'##. In my case, the friction term is very small and the particle will oscillate billions of times before the magnitude decreases significantly. So how do we solve this kind of equations to a very late final moment when the magnitude almost dies out?
Mathematica or Matlab doesn't work because of the stiffness problem, i.e. the particle oscillates so fast that we need to set the step size extremely small in order to have reliable numerical result. Another difficulty is the potential energy is not symmetric around the minimum. Thus I guess we cannot use the approximate form for the solution :##y(x)=A(t) cos (\omega(t)t)##.
Thanks a lot for any kind of suggestion!

 

[mfb]

Another difficulty is the potential energy is not symmetric around the minimum.

Your potential looks very symmetric.
If the x⁴ term is relevant, then the simple cosine won't work any more.

The small b makes it possible to separate the undamped motion and the damping (will decrease the amplitude over time).
Something like A(t+T) as function of A(t) should be possible to find, where A(t+T) means the amplitude one oscillation after A(t).

Do you have to get the phase right after so many oscillations, or just the amplitude and shape of oscillation?

 

[jasonRF]

As a first guess, I would try to use multiple scale analysis:
http://en.wikipedia.org/wiki/Multiple-scale_analysis

That may yield an approximate solution for you.

jason

 

[karlzr]

Your potential looks very symmetric.
If the x⁴ term is relevant, then the simple cosine won't work any more.

The small b makes it possible to separate the undamped motion and the damping (will decrease the amplitude over time).
Something like A(t+T) as function of A(t) should be possible to find, where A(t+T) means the amplitude one oscillation after A(t).

Do you have to get the phase right after so many oscillations, or just the amplitude and shape of oscillation?

Phase is not important in my case. So how do we find ##A(t+T)## given ##A(t)##?
Actually ##c<0##, ##d>0##

 

[Strum]

What values for d and c do you want to investigate? I just quickly wrote some code which may work and I like to check it before I show you( if you want to see it that is).

 

[JJacquelin]

The second order ODE is a particular case of Duffing differential equation :
http://mathworld.wolfram.com/DuffingDifferentialEquation.html
It can be reduced to a first order ODE of the Abel's kind (see attachement)
Except for particular values of the parameters, it is not analytically solvable :
http://eqworld.ipmnet.ru/en/solutions/ode/ode0124.pdf
Nevertheless, the numerical methods of solving might be simpler for the first order Abel's ODE than for the second order Duffing's ODE.

Sorry I don't succed to joint the attachment.
I do not have enough time to rewrite it in Latex. Roughly, the method is based on the change of function : dx/dt = -b y(x)

 

[karlzr]

What values for d and c do you want to investigate? I just quickly wrote some code which may work and I like to check it before I show you( if you want to see it that is).

We can assume ##c=-4000##, ##d=0.03##, ##b=10^{-13}##. The initial condition can be ##x(0)=10##, ##x'(0)=0##. Thank you so much!

 

[karlzr]

Roughly, the method is based on the change of function : dx/dt = -b y(x)

Do you mean we solve the first order equation with x as the independent variable? Then it can only being solved in only half oscillation, right? since x does not change monotonically.

 

[mfb]

Phase is not important in my case. So how do we find ##A(t+T)## given ##A(t)##?
Actually ##c<0##, ##d>0##

Something that will give a good approximation: evaluate one oscillation each numerically with b=0 and a few different starting amplitudes. For each starting amplitude also calculate the total energy (which is just the potential at start).

b acts as friction - determine the total energy that gets lost in one oscillation if x(t) follows the previous calculation but with a non-zero b (should be the integral over ##bx'^2##).

Determine energy loss as function of energy (alternatively: amplitude change as function of amplitude). This is now independent of the individual oscillations, you can set up a differential equation and solve this.

 

[Strum]

We can assume ##c=-4000##, ##d=0.03##, ##b=10^{-13}##. The initial condition can be ##x(0)=10##, ##x'(0)=0##. Thank you so much!

I think you should choose some other parameters :).

 

[mfb]

Maybe he likes working with large x values? The x⁴ term makes the state bound.
It does make the choice of x(0)=10 a bit unusual, right. Could be a typo.

 

[karlzr]

I think you should choose some other parameters :).

What is wrong with the numbers? ##x=0## is unstable in this case and the oscillation will be around a minimum ##x \approx 170##. We can start from any ##x(0)##between 0 and 250. Here is the plot of the potential energy given the parameters:

 

[Strum]

Well in principle you can start from anything really. But choosing your important parameters with16 orders of magnitude difference can be complicated for the numerics.

 

[mfb]

Well the small b was the whole point of the thread I think. It allows (and makes it necessary) to treat it as small perturbation.

 

[Strum]

Haha I completely missed that!

 

[karlzr]

Well the small b was the whole point of the thread I think. It allows (and makes it necessary) to treat it as small perturbation.

Oh sorry, I didn't make the most important point clear.
The physical time is ##\tau=10^{-24}t ## rather than ##t##, so the equation should be ##x''(\tau)+bx'(\tau)+cx(\tau)+d x^3(\tau)=0##. There is the fast oscillating problem only in the physical time ##\tau##. With the parameters given before, the period of oscillation is about ##10^{-24} s## , which is far less than the scale ##10^{-10} s## I am interested in. I need to find the magnitude of oscillation at ## \tau_f=10^{-10}## given initial conditions at ##\tau_0=0##. So it is not likely to be solved numerically with finite number of steps.
I am so sorry for this inconvenience!

 

[karlzr]

Determine energy loss as function of energy (alternatively: amplitude change as function of amplitude). This is now independent of the individual oscillations, you can set up a differential equation and solve this.

I don't understand the last part. So I can get the energy loss a couple of times numerically given different input energy. How do I generalize those data into a differential equation?

 

[mfb]

I would expect some reasonable relation between energy and energy loss per oscillation - probably close to a quadratic function for small amplitudes. Certainly something you can approximate well with multiple data points.

EW symmetry breaking at the big bang?