function out = beam_em_dd(lmb,l,xp)
%-------------------------------------------------------------------------
% Function name: beam_em
% Description: calculates SECOND DERIVATIVE OF THE two-dimensional beam/plate dry natural modes of vibration
% Note: natural modes are normalized in order to have inner product equal
%       to Kronecker symbol
%-------------------------------------------------------------------------

% Variables
% lmb   :   mode eigenvalue
% l     :   beam/plate half-length
% xp    :   x coord. (calculation points) vector
% out   =   value of the normalized eigenfunction.


out = zeros(length(lmb),length(xp));

for j=1:length(lmb)
    for k=1:length(xp)
        if j-1 == 0               % 1st (rigid) mode dd (dd=double derivative)
            out(j,k) = 0;
        elseif j-1 == 1           % 2nd (rigid) mode dd
            out(j,k) = 0;
        elseif rem(j-1,2) == 0    % symmetric modes
            out(j,k) = 1/sqrt(2*l)*(-lmb(j)^2*cos(lmb(j)*xp(k))/cos(lmb(j)*l)+lmb(j)^2*cosh(lmb(j)*xp(k))/cosh(lmb(j)*l));
        else                      % skew-symmetric modes
            out(j,k) = 1/sqrt(2*l)*(-lmb(j)^2*sin(lmb(j)*xp(k))/sin(lmb(j)*l)+lmb(j)^2*sinh(lmb(j)*xp(k))/sinh(lmb(j)*l));
        end
    end
end