% Sensitivity of load-dependent compliance w.r.t nodes: plates
function dFdxU = assmdFdXU(U,nodes,plates,loadmag,th,massden)
nplates = size(plates,1);
dFdxU = zeros(size(nodes,1)*2,1);
% Loop over plates
for i = 1:nplates
    n1 = plates(i,1); n2 = plates(i,2); n3 = plates(i,3);
    plate_X = [nodes(n1,1) nodes(n2,1) nodes(n3,1)];
    plate_Y = [nodes(n1,2) nodes(n2,2) nodes(n3,2)];
    x1 = plate_X(1);  x2 = plate_X(2); x3 = plate_X(3);
    y1 = plate_Y(1);  y2 = plate_Y(2); y3 = plate_Y(3);
    dAdx1 = 0.5*(y2-y3); dAdy1 = 0.5*(x3-x2);
    dAdx2 = 0.5*(y3-y1); dAdy2 = 0.5*(x1-x3);
    dAdx3 = 0.5*(y1-y2); dAdy3 = 0.5*(x2-x1);
    loaddof = [6*n1-3 6*n2-3 6*n3-3];
    dFU = -(loadmag+th*massden)/3*sum(U(loaddof,1));
    dFdxU(2*n1-1) = dFdxU(2*n1-1) + dFU*dAdx1;
    dFdxU(2*n1) = dFdxU(2*n1) + dFU*dAdy1;
    dFdxU(2*n2-1) = dFdxU(2*n2-1) + dFU*dAdx2;
    dFdxU(2*n2) = dFdxU(2*n2) + dFU*dAdy2; 
    dFdxU(2*n3-1) = dFdxU(2*n3-1) + dFU*dAdx3;
    dFdxU(2*n3) = dFdxU(2*n3) + dFU*dAdy3;
end