function [ gx, gy, Cx, Cy ] = closed_link_constraints ...
( theta1, theta2, theta3, theta4, omega, l1, l2, l3, l4, x1, y1, x3, y3, alpha )
% •ÂƒŠƒ“ƒN‹@�\ �§–ñˆÀ’è‰»–@
    %theta1 = theta(1); theta2 = theta(2); theta3 = theta(3); theta4 = theta(4);
    C1 = cos(theta1); S1 = sin(theta1);
    C12 = cos(theta1+theta2); S12 = sin(theta1+theta2);
    C3 = cos(theta3); S3 = sin(theta3);
    C34 = cos(theta3+theta4); S34 = sin(theta3+theta4);
    
    X = l1*C1+l2*C12+x1 - (l3*C3+l4*C34+x3);
    Y = l1*S1+l2*S12+y1 - (l3*S3+l4*S34+y3);
    gx = [ -l1*S1-l2*S12; -l2*S12; ...
            l3*S3+l4*S34;  l4*S34 ];
    gy = [  l1*C1+l2*C12;  l2*C12; ...
           -l3*C3-l4*C34; -l4*C34 ];
    Hx = [ -l1*C1-l2*C12, -l2*C12, 0, 0; ...
                 -l2*C12, -l2*C12, 0, 0; ...
           0, 0, l3*C3+l4*C34, l4*C34; ...
           0, 0,       l4*C34, l4*C34 ];
    Hy = [ -l1*S1-l2*S12, -l2*S12, 0, 0; ...
                 -l2*S12, -l2*S12, 0, 0; ...
           0, 0, l3*S3+l4*S34, l4*S34; ...
           0, 0,       l4*S34, l4*S34 ];
    Cx = omega'*Hx*omega + 2*alpha*gx'*omega + alpha^2*X;
    Cy = omega'*Hy*omega + 2*alpha*gy'*omega + alpha^2*Y;
end
