%
disp(' ');
disp(' long_force_sine.m   ver 1.2   April 11, 2013 ');
disp(' ');
disp(' by Tom Irvine ');
disp(' ');
disp(' Longitudinal vibration of a fixed-free rod, ');
disp(' subjected to an axial sinusoidal force. ');
disp(' ');
%
disp(' Select units:  1=English  2=metric ');
iu=input(' ');
%
[E,I,rho,c,L,A]=geometry_materials_units(iu);
%
%
damp=input(' Enter the modal damping ');
%
if(iu==1)
    P=input(' Enter the force (lbf) ');
else
    P=input(' Enter the force (N) ');    
end
f=input(' Enter the forcing frequency (Hz) ');
%
Q=2*P/(rho*A*L);
%
omega=2*pi*f;
%
j=sqrt(-1);
%
u=0;
strain=0;
disp(' ');
disp(' Natural Frequency (Hz) ');
for(k=1:4)
%    
    i=2*k-1;
%    
    omegan=i*pi*c/(2*L);
    omegad=omegan*sqrt(1-damp^2);
%
    out1=sprintf(' %8.4g ',omegan/(2*pi));
    disp(out1);
%
    den=omegan^2-omega^2 + j*2*damp*omega*omegan;
    arg1=omegan*L/c;
%    
    x=L;
    arg2=omegan*x/c;
    num=sin(arg1)*sin(arg2);
    u=u+num/den;
%
    x=0;
    arg2=omegan*x/c;
    num=(omegan/c)*sin(arg1)*cos(arg2);
    strain=strain+num/den;
%
end
disp(' ');
u=u*Q;
mu=max(abs(u));
%
strain=strain*Q;
ms=max(abs(strain));
R=E*A*ms;
%
disp(' ');
disp(' Maximum Response Values ');
disp(' ');
%
out3=sprintf('  micro strain at x=0   %8.4g \n',ms*1.0e+06);
%
if(iu==1)
    out1=sprintf('  displacement at x=L   %8.4g in ',mu);
    out2=sprintf('  velocity     at x=L   %8.4g in/sec \n',mu*omega);    
    out4=sprintf('  reaction force at x=0   %8.4g lbf \n',R);
else
    out1=sprintf('  displacement at x=L   %8.4g mm ',mu*1000);
    out2=sprintf('  velocity     at x=L   %8.4g m/sec \n',mu*omega);    
    out4=sprintf('  reaction force at x=0   %8.4g N \n',R);
end
%
disp(out1);
disp(out2);
disp(out3);
disp(out4);