...

function [Wx,Vx,Mx,Thetax,Deltax]= ...
  seval(xloc,Nterms,Nloads,Power,A, ...
        W,V,M,Theta,Delta);
%
% [Wx,Vx,Mx,Thetax,Deltax]= ...
%   seval(xloc,Nterms,Nloads,Power,A, ...
%         W,V,M,Theta,Delta);
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% This function evaluates the discontinuity
% functions at a specified value of x.
%
% xloc   - location along beam
% Nterms - vector containing the number 
%          of terms required to model 
%          each type of loading
% Nloads - number of loadings
% Power  - matrix of exponents for each
%          load term
% A      - matrix of a values for the
%          bracket part of function for
%          each load term
% W      - matrix of load coefficients 
%          for each bracket term
% V      - matrix of shear coefficients
%          for each bracket term
% M      - matrix of moment coefficients
%          for each bracket term
% Theta  - matrix of rotation coefficients
%          for each bracket term
% Delta  - matrix of deflection coefficients
%          for each bracket term
% 
% Wx     - load at specified x
% Vx     - shear at specified x
% Mx     - moment at specified x
% Thetax - rotation at specified x
% Deltax - deflection at specified x
%
% User m functions called:  none.
%----------------------------------------------

Wx=0; Vx=0; Mx=0; Thetax=0; Deltax=0;
for j=1:Nloads
  for k=1:Nterms(j)
    bracket=xloc-A(j,k);
    if bracket >= 0
      if Power(j,k) >= 0
        %...W-x
        tmp=1;
        if Power(j,k) > 0
          tmp=bracket^Power(j,k);
        end
        Wx=Wx+W(j,k)*tmp;
      end
      if Power(j,k)+1 >= 0
        %...V-x
        tmp=1;
        if Power(j,k)+1 > 0
          tmp=bracket^(Power(j,k)+1);
        end
        Vx=Vx+V(j,k)*tmp;
      end
      %...M-x
      tmp=1;
      if Power(j,k)+2 > 0
        tmp=bracket^(Power(j,k)+2);
      end
      Mx=Mx+M(j,k)*tmp;
      %...Theta-x
      Thetax=Thetax+ ...
             Theta(j,k)*bracket^(Power(j,k)+3);
      %...Delta-x
      Deltax=Deltax+ ...
             Delta(j,k)*bracket^(Power(j,k)+4);
    end
  end
end