How to solve 4th order ODE by MAXIMA?

Ronald
2013-11-21
2013-11-28
  • Ronald

    Ronald - 2013-11-21

    Hi,
    I would like to solve a 4th order ODE : y""-5y"+4y=0, and I tried to use "ode2" and "contrib_ode" to solve it but it seems not work:
    ode2 method
    de : 'diff(y,x,4)-5('diff(y,x,2))+4y;
    gsoln : ode2(de, y, x);
    msg1
    false
    contrib_ode method
    load(contrib_ode)$
    contrib_ode(de,y,x);
    false
    Can someone help me to solve this question ? thanks!!

     
  • Albrecht Mueller

    According to the description the function ode2 solves ordinary differential equations of first or second order only. You can transform your equation into a system of linear ordinary differential equations and use the desolve function. General idea:

    y4(x)=diff(y3(x),x), ..., y1(x)=diff(y(x),x), y4(x) - 5*y2(x) + 4*y(x) = 0
    
     
  • Aleksas

    Aleksas - 2013-11-23

    I have created "odes" package - 20 new functions for solving
    ordinary differential equations. One function is odeL.

    odeL
    Function: odeL(eqn, dvar, ivar)
    The function odeL solves an linear ODEs with constant coefficients.

    Examples:
    (%i1) load(odes)$ load(contrib_ode)$
    1. y''' - 2y'' + y' = 0.
    (%i3) eq:'diff(y,x,3)-2'diff(y,x,2)+'diff(y,x) = 0$
    (%i4) odeL(eq,y,x);
    (%o4) y=x
    %e^xC3+%e^xC2+C1

    1. y""-5y"+4y=0
      (%i5) eq:'diff(y,x,4)-5'diff(y,x,2)+4y=0$
      (%i6) odeL(eq,y,x);
      (%o6) y=%e^(2x)C4+%e^xC3+%e^(-x)C2+%e^(-2x)C1
      (%i7) fs(eq,y,x);
      (%o7) [%e^(-2x),%e^(-x),%e^x,%e^(2x)]
      (%i8) solve(k^4-5*k^2+4=0,k);
      (%o8) [k=-2,k=2,k=-1,k=1]

    2. y'''' + 8y'' + 16y = x ^7+sin(x)^3
      (%i9) eq:'diff(y,x,4)+8'diff(y,x,2)+16y=x^7+sin(x)^3$
      (%i10) sol:odeL(eq,y,x);
      (%o10) y=xsin(2x)C4+sin(2x)C3+xcos(2x)C2+cos(2x)C1-(24sin(3x)-200sin(x)-150x^7+3150x^5-23625x^3+47250*x)/2400
      (%i11) ode_check(eq,sol);
      (%o11) 0

    3. y'''- 3y'' + y = xsin(x)exp(x).
      (%i12) eq:'diff(y,x,3)-3'diff(y,x,2)+y=xsin(x)exp(x);
      (%o12) 'diff(y,x,3)-3
      ('diff(y,x,2))+y=x%e^xsin(x)
      (%i13) solvet(k^3-3k^2+1=0,k);
      (%o13) [k=2cos(%pi/9)+1,k=2cos((5%pi)/9)+1,k=2cos((7*%pi)/9)+1]
      (%i14) sol:odeL(eq,y,x);
      (%o14) y=%e^(2
      cos((7%pi)/9)x+x)C3+%e^(2cos((5%pi)/9)x+x)C2+%e^(2cos(%pi/9)x+x)C1-(%e^x(x(17sin(x)-68cos(x))+90sin(x)+48cos(x)))/289
      (%i15) ode_check(eq,sol);
      (%o15) 0

    best

    Aleksas Domarkas

     
  • Ronald

    Ronald - 2013-11-25

    Thanks for all your response !

    Dear Aleksas,
    Thank you for your information. Could you give me the ode package you created just for course teaching (school) ?
    By the way, if the initial condition for 4th oder ODE : y'"(0)=5, y"(0)=3, y'(0)=0, y(0)=1, how to solve it by MAXIMA (I think the function "ic2" does not work on it).
    Thanks for your help!!

     
    • Joseph Cusumano

      Joseph Cusumano - 2013-11-25

      It's nice to have a package to use (I do it all of the time), but you
      might find it useful to learn how to solve these equations "by hand"
      using Maxima. It's not that hard, assuming you know how to do it with
      pencil and paper.

      Here's an example, for your problem, that can be easily modified for
      other ODEs.

      /* coefficients for 4th order ODE */
      c:[1,0,-5,0,4];
      /* assemble ODE */
      ODE:sum(c[k]*'diff(y,x,k-1),k,1,5);
      
      /* Find the characteristic polynomial */
      subst(y=exp(s*x),ODE);
      charPoly:coeff(factor(ev(%,eval,diff)),exp(s*x));
      
      /* Find roots to charPoly, and assemble general solution */
      sSols:solve(charPoly,s);
      yGen:sum(concat('A,k)*exp(rhs(sSols[k]*x)),k,1,4);
      
      /* Initial values of derivatives at x=0 */
      ICvals:[1,0,3,5];
      /* IC equations */
      ICs:makelist(subst(x=0,diff(yGen,x,k-1))=ICvals[k],k,1,4);
      
      /* Solve for coefficients in general solution */
      ASols:solve(ICs,[A1,A2,A3,A4])[1];
      
      /* Obtain desired solution to ODE */
      ySol:subst(ASols,yGen);
      

      Best regards,

      Joe


      Joseph Cusumano
      Professor of Engineering Science & Mechanics
      Penn State University
      212 Earth & Engineering Sciences Building
      University Park, PA 16802
      814.865.3179 (o)
      814.865.9974 (f)

       
  • Aleksas

    Aleksas - 2013-11-25

    (%i1) load(odes);
    (%o1) "C:/Users/Aleksas/maxima/odes.mac"

    y'''' - 5y'' + 4y = 0, y(0)=1, y'(0)=0, y''(0)=3, y'''(0)=5

    For this I have "odeL_ic":

    (%i2) eq:'diff(y,x,4)-5'diff(y,x,2)+4y=0$
    (%i3) sol:odeL_ic(eq,y,x,[0,1,0,3,5]);
    (%o3) y=(3%e^(2x))/4-(2%e^x)/3+%e^(-x)-%e^(-2x)/12

    Documentation see in
    http://uosis.mif.vu.lt/~aleksas/DL/DLsuMaxima/odes-doc.pdf

    best

    Aleksas Domarkas

     
  • Ronald

    Ronald - 2013-11-25

    Dear Aleksas,
    I tried to find the function "odes", but it showed the following message:
    load(odes);
    file_search1: odes not found in file_search_maxima,file_search_lisp
    -- an error. To debug this try: debugmode(true);
    Did I do anything wrong ?

     
  • Aleksas

    Aleksas - 2013-11-25

    On windows you need put odes.mac similar to

    "C:/Users/Aleksas/maxima/odes.mac"

    in

    "C:/Users/user_name/maxima/"

    or in maxima instalation directory share/...

    best

    Aleksas D

     
  • Loukas Zachilas

    Loukas Zachilas - 2013-11-27

    Dear Aleksas,
    I am using (and teaching) differential equations too. How could I get (like Ronald asked) the package: ODES.MAC?
    Would you, please, send it to me?

    Best regards,
    Loukas

     
  • Loukas Zachilas

    Loukas Zachilas - 2013-11-28

    Dear Aleksas,
    thanks for the link of the Maxima Files, and a lot of congratulations for the very good educational material that you have created on Differential Equations.

    best regards,
    Loukas

     
  • Loukas Zachilas

    Loukas Zachilas - 2013-11-28

    Dear Aleksas,
    I tried to load ODES.MAC, but I have got the following message:

    (%i1) load("/Users/loukas/Downloads/odes.mac")$
    Maxima encountered a Lisp error:
    decoding error on stream
    #<SB-SYS:FD-STREAM for="" "file="" Users="" loukas="" Downloads="" odes.mac"="" {14879499}="">
    (:EXTERNAL-FORMAT :ASCII):
    the octet sequence (251) cannot be decoded.
    Automatically continuing.
    To enable the Lisp debugger set debugger-hook to nil.

    Any suggestions?

    Loukas

     
  • Aleksas

    Aleksas - 2013-11-28

    with wxMaxima 13.04.2

    Available directories you can see with
    (%i1) file_search_maxima;
    (%o1) ["C:/Users/Aleksas/maxima/$$$.{mac,mc}","C:\PROGRA~2\MAXIMA~1.2/share/maxima/5.31.2/share/$$$.{mac,mc}","C:\PROGRA~2\MAXIMA~1.2/share/maxima/5.31.2/share/{affine,algebra,algebra/charsets,algebra/solver,amatrix,bernstein,calculus,cobyla,cobyla/ex,cobyla/lisp,colnew,colnew/ex1,colnew/ex2,colnew/ex3,colnew/ex4,colnew/lisp,combinatorics,contrib,contrib/Grobner,contrib/Zeilberger,contrib/altsimp,contrib/bitwise,contrib/boolsimp,contrib/diffequations,contrib/diffequations/tests,contrib/format,contrib/fresnel,contrib/gentran,contrib/gentran/man,contrib/gentran/test,contrib/gf,contrib/integration,contrib/levin,contrib/lurkmathml,contrib/maximaMathML,contrib/mcclim,contrib/namespaces,contrib/noninteractive,contrib/prim,contrib/rand,contrib/rkf45,contrib/sarag,contrib/state,contrib/unit,contrib/vector3d,descriptive,diff_form,diffequations,distrib,draw,dynamics,ezunits,finance,fourier_elim,fractals,graphs,hypergeometric,integequations,integer_sequence,integration,lapack,lapack/blas,lapack/lapack,lbfgs,linearalgebra,lsquares,macro,matrix,minpack,minpack/lisp,misc,mnewton,multiadditive,numeric,numericalio,orthopoly,pdiff,physics,simplex,simplex/Tests,simplification,solve_rat_ineq,solve_rec,stats,stringproc,sym,tensor,to_poly_solve,trigonometry,utils,vector,z_transform}/$$$.{mac,mc}"]
    (%i2) load(odes);
    (%o2) "C:/Users/Aleksas/maxima/odes.mac"

    Add directory you can if from wxMaxima menu select:
    Maxima>Add to Path..> Galutinis

    (%i3) file_search_maxima : cons(sconcat("C:/Users/Aleksas/Desktop/Galutinis/###.{lisp,mac,mc}"), file_search_maxima)$

    (%i4) load(odes);
    (%o4) "C:/Users/Aleksas/Desktop/Galutinis/odes.mac"

    (%i5) file_search_maxima;
    (%o5) ["C:/Users/Aleksas/Desktop/Galutinis/###.{lisp,mac,mc}",
    "C:/Users/Aleksas/maxima/$$$.{mac,mc}",
    "C:\PROGRA~2\MAXIMA~1.2/share/maxima/5.31.2/share/$$$.{mac,mc}","C:\PROGRA~2\MAXIMA~1.2/share/maxima/5.31.2/share/{affine,algebra,algebra/charsets,algebra/solver,amatrix,bernstein,calculus,cobyla,cobyla/ex,cobyla/lisp,colnew,colnew/ex1,colnew/ex2,colnew/ex3,colnew/ex4,colnew/lisp,combinatorics,contrib,contrib/Grobner,contrib/Zeilberger,contrib/altsimp,contrib/bitwise,contrib/boolsimp,contrib/diffequations,contrib/diffequations/tests,contrib/format,contrib/fresnel,contrib/gentran,contrib/gentran/man,contrib/gentran/test,contrib/gf,contrib/integration,contrib/levin,contrib/lurkmathml,contrib/maximaMathML,contrib/mcclim,contrib/namespaces,contrib/noninteractive,contrib/prim,contrib/rand,contrib/rkf45,contrib/sarag,contrib/state,contrib/unit,contrib/vector3d,descriptive,diff_form,diffequations,distrib,draw,dynamics,ezunits,finance,fourier_elim,fractals,graphs,hypergeometric,integequations,integer_sequence,integration,lapack,lapack/blas,lapack/lapack,lbfgs,linearalgebra,lsquares,macro,matrix,minpack,minpack/lisp,misc,mnewton,multiadditive,numeric,numericalio,orthopoly,pdiff,physics,simplex,simplex/Tests,simplification,solve_rat_ineq,solve_rec,stats,stringproc,sym,tensor,to_poly_solve,trigonometry,utils,vector,z_transform}/$$$.{mac,mc}"]

    best
    Aleksas

     
  • Loukas Zachilas

    Loukas Zachilas - 2013-11-28

    Dear Aleksas,
    I don't believe it is a matter of wrong directory.
    I have copied the file odes.mac in the right directory:

    (%i1) file_search_maxima;
    (%o1) ["/Users/loukas/.maxima/$$$.{mac,mc}","/Applications/Maxima.app/Contents/Resources/maxima//share/maxima/5.30.0/share/$$$.{mac,mc}","/Applications/Maxima.app/Contents/Resources/maxima//share/maxima/5.30.0/share/{affine,algebra,algebra/charsets,algebra/solver,amatrix,bernstein,calculus,cobyla,cobyla/ex,cobyla/lisp,colnew,colnew/lisp,combinatorics,contrib,contrib/Grobner,contrib/Zeilberger,contrib/altsimp,contrib/bitwise,contrib/boolsimp,contrib/diffequations,contrib/diffequations/tests,contrib/format,contrib/fresnel,contrib/gentran,contrib/gentran/man,contrib/gentran/test,contrib/gf,contrib/integration,contrib/levin,contrib/lurkmathml,contrib/maximaMathML,contrib/mcclim,contrib/namespaces,contrib/noninteractive,contrib/prim,contrib/rand,contrib/rkf45,contrib/sarag,contrib/state,contrib/unit,contrib/vector3d,descriptive,diff_form,diffequations,distrib,draw,dynamics,ezunits,finance,fourier_elim,fractals,graphs,hypergeometric,integequations,integer_sequence,integration,lapack,lapack/blas,lapack/lapack,lbfgs,linearalgebra,lsquares,macro,matrix,minpack,minpack/lisp,misc,mnewton,multiadditive,numeric,numericalio,orthopoly,pdiff,physics,simplex,simplex/Tests,simplification,solve_rat_ineq,solve_rec,stats,stringproc,sym,tensor,to_poly_solve,trigonometry,utils,vector,z_transform}/$$$.{mac,mc}"]
    (%i2) load(odes);
    Maxima encountered a Lisp error:
    decoding error on stream
    #<SB-SYS:FD-STREAM for="" "file="" Applications="" Maxima.app="" Contents="" Resources="" maxima="" share="" maxima="" 5.30.0="" share="" diffequations="" odes.mac"="" {149634B9}="">
    (:EXTERNAL-FORMAT :ASCII):
    the octet sequence (251) cannot be decoded.
    Automatically continuing.
    To enable the Lisp debugger set debugger-hook to nil.

    I use a Mac Pro (computer) and I am sure that I have copied the file in the right folder (as you may see when I try to load it).
    Is it a problem of ASCII format?
    I download it from the link that you have sent me (Save the download link as...)
    And I choose the keep the extension (.mac).
    Is there any case that a character in the file cannot be decoded?
    Any suggestions from you or from Andrejv??

    Best regards,
    Loukas

     
  • Loukas Zachilas

    Loukas Zachilas - 2013-11-28

    Dear Aleksas,
    Ok. Trouble is solved!
    I have just rewritten the commands (line-by-line) from the scratch!
    I created the file odes.mac, saved it in the right place and everything is working!
    Thanks again.

    best regards,
    Loukas

     

Log in to post a comment.

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:

JavaScript is required for this form.





No, thanks