# D'Alembert's Method
# (Section 3.4)
# Example 1  (See Example 1, Section 3.4 for details.)
#
# To be able to use d'Alembert's method, we need to define the odd
# extension of a given function f(x), where f(x) denotes the initial
# shape of the string. # Here is the odd extension on the interval 
# [-1,1].

f:= x->piecewise(0<x and x<1/3, 3*x/10, 1/3<=x and x<=1,(3*(1-x))/20) ;

fodd:=x->piecewise(0<=x and x<=1, f(x), -1<=x and x<0,-f(-x));
plot(fodd(x),x=-2..2);

# Now we wish to continue the portionof the graph seen between x = -1;
# and x = 1; periodically to the whole real line. We will use the
# formula shown in exercise 20 on p.22. For Maple, floor is the same as
# the greatest integer function mentioned in the text. We chose the
# function name extend to remind us that we are extending part of a
# graph.

extend:=(x,p)->x-2*p*floor( (x+p)/(2*p));
plot( extend(x,1),x=-4..4);

# This gives a periodic extension of the basic graph y = x; on the
# interval [-1,1]. To get extensions of other functions, compose them
# with extend. For example, we will periodically extend our odd
# extension fodd. In fact, we could make that our main function fstar.

fstar:=x->fodd(extend(x,1));
plot(fstar(x) , x=-3..3);

# We can now use D'Alembert's method to get the shape of the string at
# various value of t; by simply translating then averaging the graphs of
# the odd periodic extension fstar(x);.  For example, since c = 1/Pi;,
# we have u(x, t) = (fstar(x-t/Pi)+fstar(x+t/Pi))*(1/2);. We define
# this and plot on the interval [0,1] for t = (1/3)*Pi;.

u:=(x,t)-> 1/2*(fstar(x-1/Pi*t)+fstar(x+1/Pi*t));
plot( u(x,Pi/3),x=0..1);

# We can compare this shape to what we would get by using the series
# solution that we derived in the previous section.  Recall that a good
# approximation of this function was denoted S(10,x,t) in the worksheet
# for Section 3.3 entitled string.mws.  We will recall from this
# notebook the necessary data for the analytical solution.  So,

b:=n->2*int(f(x)*sin(n*Pi*x),x=0..1);
assume('n',integer);
b(n);
S:=(N,x,t)->sum( b(n)*sin(n*Pi*x)*cos(n*t),n=1..N);
plot( [S(10,x,Pi/3),u(x,Pi/3)],x=0..1);

# It is important to note that the plot that we obtained from
d'Alembert's solution is the exact shape of the string. # Let us check
the shape of the string when t = 2*Pi;. This shows that the string has
returned ot it's original shape.

plot( u(x,2*Pi),x=0..1);

# At this point, it is easy to animate several graphs as time increases.

with(plots):
animate( u(x,t),x=0..1,t=0..2*Pi,frames=12);

# What you see in the above pictures are the exact shapes of the string
# at the stated times.  There is no approximation here.  D'Alembert's
# method gives you the precise shape of the string.  We can arrange the
# plots in an array as follows withthe shortcut % assuming that these
# commands are executed immediately after the animate command..

display(%,axes=none);

#