# Advanced Mathematics and Mechanics Applications Using MATLAB, 3rd Edition

### Howard Wilson (view profile)

14 Oct 2002 (Updated )

Companion Software (amamhlib)

strwave(rundemo)
```function strwave(rundemo)
% strwave(rundemo)
% This function illustrates propagation of
% waves in a tightly stretched string having
% given initial deflection. A function call of
% strngrun inputs data interactively. Otherwise,
% strngrun(1) executes a sample data case.

pltsav=1; % flag to save or not save graphs

disp(' ')
disp('WAVE PROPAGATION IN A STRING'), disp(' ')
if nargin==0 % Input data interactively
L='['; R=']';
v=input(['Input wave speed (a) and ',...
'string length (len) >? '],'s');
v=eval([L,v,R]); a=v(1); len=v(2);
disp(' ')
xd=input('Input the xdata values >? ','s');
xd=eval([L,xd,R]); disp(' ')
yd=input('Input the ydata values >? ','s');
yd=eval([L,yd,R]); disp(' ')
disp('Input tmax and the number of time steps')
d=input('(Try len/a and 40) >? ','s');
d=eval([L,d,R]); tmax=d(1); nt=d(2);
disp(' ')
disp('Specify position x=x0 where the time')
x0=input('history is to be evaluated (try len/4) >? ');
disp(' ')
disp('Specify time t=t0 when the deflection')
t0=input('curve is to be plotted >? ');
disp(' ')
titl=input('Input a graph title >? ','s');

else % Example for triangular initial deflection
a=1; len=1; tmax=len/a; nt=40;
xd=[0,.33,.5,.67,1]*len; yd=[0,0,-1,0,0];

% Different example for a truncated sine curve
% xd=linspace(0,len,351); yd=sin(3*pi/len*xd);
% k=find(yd<=0); xd=xd(k); yd=yd(k);

x0=0.25*len; t0=0.33*len/a;
titl='TRANSLATING WAVE OVER HALF A PERIOD';
end

nx=80; x=0:len/nx:len; t=0:tmax/nt:tmax;

h=max(abs(yd)); xplot=linspace(0,len,201);
tplot=linspace(0,max(t),251)';

[Y,X,T]=strngwav(xd,yd,x,t); plot3(X',T',Y','k');
xlabel('x axis'), ylabel('time'), zlabel('y(x,t)')
title(titl)
if pltsav, print(gcf,'-deps','strngplot3'); end
drawnow, shg

disp(['when t = ',num2str(t0)]), pause

[yt0,xx,tt]=strngwav(xd,yd,xplot,t0,len,a);
close; plot(xx(:),yt0(:),'k')
xlabel('x axis'), ylabel('y(x,t0)')
title(['DEFLECTION WHEN T = ',num2str(t0)])
axis([min(xx),max(xx),-h,h])
if pltsav, print(gcf,'-deps','strngyxt0'); end
drawnow, shg

disp(' ')
disp(['at x = ',num2str(x0)]), pause

yx0=strngwav(xd,yd,x0,tplot,len,a);
plot(tplot,yx0,'k')
xlabel('time'), ylabel('y(x0,t)')
title(...
['DEFLECTION HISTORY AT X = ',num2str(x0)])
axis([0,max(t),-h,h])
if pltsav, print(gcf,'-deps','strngyx0t'); end
drawnow, shg

disp(' ')
disp('over two periods of motion'), pause
x=linspace(0,len,81); t=linspace(0,4*len/a,121);
[Y,X,T]=strngwav(xd,yd,x,t);
titl='MOTION OVER TWO PERIODS';
animate(X(1,:),Y',titl,.1), pause(2)

if pltsav, print(gcf,'-deps','strnganim'); end

disp(' '), disp('All Done')

%===============================================

function [Y,X,T]=strngwav(xd,yd,x,t,len,a)
% [Y,X,T]=strngwav(xd,yd,x,t,len,a)
% xd,yd - data vectors defining the initial
%         deflection as a piecewise linear
%         function. xd values should be increasing
%         and lie between 0 and len
% x,t   - position and time vectors for which the
%         solution is evaluated
% len,a - string length and wave speed

if nargin<6, a=1; end; if nargin <5, len=1; end
xd=xd(:); yd=yd(:);  p=2*len;

% If end values are not zero, add these points
if xd(end)~=len, xd=[xd;len]; yd=[yd;0]; end
if xd(1)~=0, xd=[0;xd]; yd=[0;yd]; end
nd=length(xd);

% Extend the data definition for len < x < 2*len
xd=[xd;p-xd(nd-1:-1:1)]; yd=[yd;-yd(nd-1:-1:1)];
[X,T]=meshgrid(x,t); xp=X+a*T; xm=X-a*T;
shape=size(xp); xp=xp(:); xm=xm(:);

% Compute the general solution for a piecewise
% linear initial deflection
Y=(sign(xp).*interp1(xd,yd,rem(abs(xp),p))...
+sign(xm).*interp1(xd,yd,rem(abs(xm),p)))/2;
Y=reshape(Y,shape);

%===============================================

function animate(x,y,titl,tim,trace)
% animsmpl(x,y,titl,tim,trace)
% This function performs animation of a 2D curve
% x,y - arrays with columns containing curve positions
%       for successive times. x can also be a single
%       vector if x values do not change. The animation
%       is done by plotting (x(:,j),y(:,j)) for
%       j=1:size(y,2).
% titl- title for the graph
% tim - the time in seconds between successive plots

if nargin<5, trace=0; else, trace=1; end;
if nargin<4, tim=.05; end
if nargin<3, trac=''; end; [np,nt]=size(y);
if min(size(x))==1, j=ones(1,nt); x=x(:);
else, j=1:nt; end; ax=newplot;
if trace, XOR='none'; else, XOR='xor'; end
r=[min(x(:)),max(x(:)),min(y(:)),max(y(:))];
%axis('equal') % Needed for an undistorted plot
axis(r), % axis('off')
curve = line('color','k','linestyle','-',...
'erase',XOR, 'xdata',[],'ydata',[]);
xlabel('x axis'), ylabel('y axis'), title(titl)
for k = 1:nt
set(curve,'xdata',x(:,j(k)),'ydata',y(:,k))
if tim>0, pause(tim), end, drawnow, shg
end
```