Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

[MATLAB] Modeling response of a particle in an anharmonic potential

  1. Apr 8, 2015 #1
    I'm trying to model the response of a particle in an anharmonic potential in MATLAB. I know that writing position, x as a function of potential energy U is,

    upload_2015-4-8_16-54-10.png

    where T is the period, m is mass, and E is the energy of the system. How would I get a function for T(E) since it's not independent of amplitude obviously. As a starting point, I'm using the anharmonic potential of

    upload_2015-4-8_17-0-7.png
    Where a and b are constants. Does anybody know how I could program this in matlab?

    Thanks
     
  2. jcsd
  3. Apr 9, 2015 #2

    kreil

    User Avatar
    Gold Member

  4. Apr 9, 2015 #3
    Yes my goal is to plot a potential vs. displacement curve. Is this even possible considering the number of functions in the first equation above?
    Here's some of my code: For simplicity, I assume period T(E) = 1


    %Define Functions

    %For a,b = 1

    %U = ax^2 + bx^3

    %Assuming T(E) = 1



    U = @(x) x^2+x^3

    fun = @(E,U) T/(sqrt(U-E))^-1

    int = @(E) integral(fun,0,U)

    ezplot('(-x+(sqrt(2)*2*pi)^-1)*int')
     
  5. Apr 9, 2015 #4

    kreil

    User Avatar
    Gold Member

    I think something like this is a step closer, but I'm not convinced it's correct yet.

    Code (Text):

    a=1;
    b=1;
    T=1;
    U = @(x) a*x.^2 + b*x.^3;
    A = (2*pi*sqrt(2))^-1;
    I = [];

    for n =0:0.1:10
      fun = @(E) T./(sqrt(U(n)-E));
      I = [I, A*integral(fun,0,U(n))];
    end

    plot(0:0.1:10, I)
     
    The difference here is that I loop through possible values for U and evaluate the integral once for each one.
     
  6. Apr 9, 2015 #5
    Thank you I'll try it. I'm wondering if it would be more beneficial/ shed more light if I plotted x as a function of time. In that case I solved lagrange
    L = K - U = .5mv^2 - x^2 - x^3 since U = x^2 + x^3
    the associated equation of motion is then ma = -3x^2 - 2x
    I would then need to solve this differential equation to plot x(t). However I put this equation into wolfram alpha, and got a massive expression out and it makes me thing that it might not be possible to model this response.

    Is this the case?
     
  7. Apr 10, 2015 #6

    kreil

    User Avatar
    Gold Member

    You can solve the second order differential equation numerically in MATLAB. Since the ODE functions only work with first order equations, you'll have to rewrite this second order equation

    [tex]
    m \frac{d^2 x}{dt^2} = -3x^2-2x
    [/tex]

    As an equivalent system of first order equations,

    [tex]
    m \dot{y_2} = -3y_1^2-2y_1
    [/tex]
    [tex]
    y_2 = \dot{x}
    [/tex]
    [tex]
    y_1 = x
    [/tex]

    You'll also need to define the initial conditions of the problem.

    Take a look at the examples here:
    http://www.mathworks.com/help/matlab/ref/ode45.html

    Save the following function file, which defines the equations above (assume m=1):
    Code (Text):

    function dy = anharmonic(t,y)
      dy = zeros(2,1);
      dy(1) = y(2);
      dy(2) = -3*y(1).^2 - 2*y(1);
    end
     
    Then set the time span and initial conditions before calling ode45:
    Code (Text):

    tspan = [0 15];
    y0 = [0; 0.25];
    [t,y] = ode45(@anharmonic, tspan, y0);
    plot(t,y(:,1),t,y(:,2))
    legend('Position','Velocity')
     
    ?temp_hash=9239f5d362fbcfce69e50037841c7f57.png

    You might want to check my work, because the plot looks more harmonic than anharmonic to me.
     

    Attached Files:

    Last edited: Apr 10, 2015
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook




Similar Discussions: [MATLAB] Modeling response of a particle in an anharmonic potential
  1. Modeling and Matlab (Replies: 3)

  2. 3D modelling in MATLAB (Replies: 2)

Loading...