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

### Howard Wilson (view profile)

14 Oct 2002 (Updated )

Companion Software (amamhlib)

shakestr
```function shakestr
% Example:  shakestr
% ~~~~~~~~~~~~~~~~~~
% This program illustrates wave motion in
% a string which is initially at rest when
% the right end is shaken with a nearly
% resonant harmonic oscillation.

% User m functions required:
%    shkstrng, ploteasy

fprintf('\nFORCED MOTION OF A VIBRATING ')
fprintf('STRING\n'), close
wf=1.98*pi; tmax=(2*pi)*.8;
[y,t,x]=shkstrng(wf,80,0,tmax,75,0,1,51);
surf(x,t,y); ylabel('time')
view([30,30]), xlabel('x axis')
zlabel('deflection'), % colormap([1 1 1]);
title(['Motion of a String with the Right', ...
' End Shaken Harmonically'])
fprintf('\nPress [Enter] for the\n')
fprintf('deflection when t=0.5\n')
figure(gcf)
pause
print -deps strngsrf
[yp5,tp5,xp5]=shkstrng(wf,80,.5,.5,1,0,1,51);
ploteasy(xp5,yp5,'horizontal direction', ...
'transverse deflection', ...
'String Deflection When t=0.5')
fprintf('Press [Enter] for the deflection\n')
fprintf('history at x=0.25\n')
print -deps dflatep5
pause

[yxc,txc,xc]= ...
shkstrng(wf,80,0,tmax,151,.25,.25,1);
ploteasy(txc,yxc,'dimensionless time', ...
'transverse deflection', ...
'Motion at x=0.25 in the String');
print -deps motnqrtp
fprintf('Press [Enter] to animate')
fprintf('\nthe string motion\n')
pause
motion(x,y)
disp(' '), disp('All Done')

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

function [y,t,x]= ...
shkstrng(w,nsum,t1,t2,nt,x1,x2,nx)
%
% [y,t,x]=shkstrng(w,nsum,t1,t2,nt,x1,x2,nx)
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Simulation of the motion of a string having
% one end fixed and the other end shaken
% harmonically.
%
% w     - forcing frequency
% t1,t2 - minimum and maximum times
% nt    - number of time values
% x1,x2 - minimum and maximum x values
%         lying between zero and one
% nx    - number of x values
%
% t,x   - vectors of time and position values
% y     - matrix of transverse deflection
%         values having nt rows and nx
%         columns
%
% User m functions called:  none
%----------------------------------------------

t=linspace(t1,t2,nt)'; x=linspace(x1,x2,nx);
np=pi*(1:nsum); y=sin(w*t)/sin(w)*sin(w*x);
a=2*w*ones(nt,1)*(cos(np)./(np.^2-w^2));
y=y+a.*sin(t*np)*sin(np'*x);

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

function motion(x,y,inct,trac)
%
% motion(x,y,inct,trac)
% ~~~~~~~~~~~~~~~~~~~~~
% This function animates the motion history
% of the string.
%
% x    - horizontal position coordinates
%        corresponding to various columns
%        of matrix y
% y    - matrix with row j specifying the
%        string position at the j'th time
%        value
% inct - the number of row increments used
%        to select positions for plotting.
%        Using inct=2 would plot every other
%        row of y. inct=1 is the default value.
% trac - if this parameter is present,
%        successive plot images are left on
%        the screen. Otherwise, each
%        configuration is shown and removed
%        before the next image is shown. The
%        default choice is to remove
%        successive images.
%
% User m functions called:  none
%----------------------------------------------

if nargin ==2, inct=1; trac=0; end
if nargin ==3, trac=0; end
if inct > 1
[nt,nx]=size(y); y=y(1:inct:nx,:);
end

xmin=min(x); xmax=max(x);
ymin=min(y(:)); ymax=max(y(:)); clf;
axis([xmin,xmax,2*ymin,2*ymax]);
[nt,nx]=size(y); axis off; hold on

for j=1:nt-1
plot(x,y(j,:)); drawnow; figure(gcf);
if trac ==0, cla; end
pause(.1)
end
plot(x,y(nt,:)); figure(gcf); hold off;

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

function ploteasy(x,y,xlabl,ylabl,titl)
%
% ploteasy(x,y,xlabl,ylabl,titl)
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Easy plot function with a simple
% argument list
%
% x,y   - data to be plotted
% xlabl - horizontal axis label for the graph
% ylabl - vertical axis label for the graph
% titl  - title for the graph
%
% User m functions called:  none
%----------------------------------------------

plot(x,y);
if nargin==2, figure(gcf), return, end
if nargin>2, xlabel(xlabl), end
if nargin>3, ylabel(ylabl), end
if nargin>4, title(titl), end
figure(gcf)```