function m = mfs_mce(elem)

# usage: m = mfs_mce(elem)
#
# Input  elem     Structure containing data of one element
#
# Output k        Element stiffness matrix
#
# The function computes the stiffness matrix of a 2-dimensional beam 
# element.
#
# -------------------------------------------------------------------

# Geometry

  d = elem.coor(2, :) - elem.coor(1, :);
  L = norm(d);
  c = d(1) / L;
  s = d(2) / L;

# Mass matrix in element coordinate system

  if (isfield(elem.mat, "rho"))
     m0 = elem.mat.rho * elem.geom.A * L / 420;
  else
     m = sparse(6, 6);
     return;
  end

  m11 = 140;
  m14 =  70;

  m22 = 156;
  m23 =  22 * L;
  m25 =  54;
  m26 = -13 * L;
  m33 =  4 * L^2;
  m35 = -m26;
  m36 = - 3 * L^2;
  m55 = 156;
  m56 = -m23;
  m66 = m33;
  
  mE = m0 * [ m11,   0,   0, m14,   0,   0; 
                0, m22, m23,   0, m25, m26;
                0, m23, m33,   0, m35, m36;
              m14,   0,   0, m11,   0,   0;
                0, m25, m35,   0, m55, m56;
                0, m26, m36,   0, m56, m66 ];

# Transformation matrix to structure coordinate system

  T33 = [c, s, 0; -s, c, 0; 0, 0, 1];
  T   = zeros(6, 6, "double");
  T(1 : 3, 1 : 3) = T33;
  T(4 : 6, 4 : 6) = T33;

# Mass matrix in structure coordinate system

  m = T' * mE * T;

end
