#
# Numerical dsolve is used with Euler's method.  The numeric solutions  
# are graphed and compared with the analytic solution.
#
# Maple Worksheet written by Denise Kirschner
# 
> with(DEtools): with(plots):  
--------------------------------------------------------------------------------
# 
> p1:=DEplot(diff(y(x),x)=2*y(x)-1,[x,y],0..0.5,{[0,1]},stepsize=.1,method=`Euler`,thickness=1);

p1 := PLOT(

CURVES(

    [[0, 1.], [.1, 1.1], [.2, 1.22], [.3, 1.364], [.4, 1.5368], [.5, 1.74416]],

    STYLE(LINE), COLOUR(RGB, 1.0, 1.0, 1.0)),

POINTS(op([[0, 1.]], COLOUR(RGB, 0, 0, 1.0))),

VIEW(-.02500000000 .. .5250000000, .9627920000 .. 1.781368000), THICKNESS(1),

AXESLABELS(x, y))
--------------------------------------------------------------------------------
# 
> p2:=DEplot(diff(y(x),x)=2*y(x)-1,[x,y],0..0.5,{[0,1]},stepsize=.05,method=`Euler`,thickness=2):
--------------------------------------------------------------------------------
> p3:=DEplot(diff(y(x),x)=2*y(x)-1,[x,y],0..0.5,{[0,1]},stepsize=.01,method=`Euler`,thickness=3):
# 
> diffeq1:=diff(y(x),x)=2*y(x)-1; inits:=y(0)=1;

                                     d
                        diffeq1 := ---- y(x) = 2 y(x) - 1
                                    dx

                                inits := y(0) = 1
> sol1:=dsolve({diffeq1,inits},y(x));

                        sol1 := y(x) = 1/2 + 1/2 exp(2 x)
> p4:=plot(rhs(sol1),x=0..0.4):
> 
> display([p1,p2,p3,p4]);
> p5:=DEplot(diff(y(x),x)=y(x)^2+x^2,[x,y],0..0.5,{[0,1]},stepsize=.1,method=`Euler`,thickness=1):
> p6:=DEplot(diff(y(x),x)=y(x)^2+x^2,[x,y],0..0.5,{[0,1]},stepsize=.05,method=`Euler`,thickness=2):
> p7:=DEplot(diff(y(x),x)=y(x)^2+x^2,[x,y],0..0.5,{[0,1]},stepsize=.01,method=`Euler`,thickness=3):
> diffeq2:=diff(y(x),x)=y(x)^2+x^2; inits:=y(0)=1;

                                     d             2    2
                        diffeq2 := ---- y(x) = y(x)  + x
                                    dx

                                inits := y(0) = 1
> sol2:=dsolve({diffeq2,inits},y(x));

sol2 := y(x) =

                                    2                    2
                        x GAMMA(3/4)  BesselY(-3/4, 1/2 x )
  ------------------------------------------------------------------------------
  /            2                   2                        \
  |  GAMMA(3/4)  BesselY(1/4, 1/2 x )                     2 |            2
  |- -------------------------------- + BesselJ(1/4, 1/2 x )| (GAMMA(3/4)  - Pi)
  |                    2                                    |
  \          GAMMA(3/4)  - Pi                               /

                                           2
                      x BesselJ(-3/4, 1/2 x )
   - ---------------------------------------------------------
                 2                   2
       GAMMA(3/4)  BesselY(1/4, 1/2 x )                     2
     - -------------------------------- + BesselJ(1/4, 1/2 x )
                         2
               GAMMA(3/4)  - Pi
> p8:=plot(rhs(sol2),x=0..0.5):
> display([p5,p6,p7,p8]);
--------------------------------------------------------------------------------
>