function out = gamli(theta)
    gamma = pi / 2
    l1 = 2
    l2 = sqrt(2)
    l3 = sqrt(2)
    p1 = sqrt(5)
    p2 = sqrt(5)
    p3 = sqrt(5)
    x1 = 4
    x2 = 0
    y1 = 0
    y2 = 4

    a2 = l3.*cos(theta) - x1
    b2 = l3.*sin(theta)
    a3 = l2.*cos(theta + gamma) - x2
    b3 = l2.*sin(theta + gamma) - y2

    p1.^2 = x.^2 + y.^2
    p2.^2 = (x + a2).^2 + (y + b2).^2
    p3.^2 = (x + a3).^2 + (y + b3).^2

    p2.^2 = x.^2 + y.^2 + a2.*x.*2 + b2.*y.*2 + a2.^2 + b2.^2
    p3.^2 = x.^2 + y.^2 + a3.*x.*2 + b3.*y.*2 + a3.^2 + b3.^2

    d = 2 * (a2.*b3 - b2.*a3)
    n1 = b3.*(p2.^2 - p1.^2 - a2.^2 - b2.^2) - b2.*(p3.^2 - p1.^2 - a3.^2 - b3.^2)
    n2 = -1*a3.*(p2.^2 - p1.^2 - a2.^2 - b2.^2) + a2.*(p3.^2 - p1.^2 - a3.^2 - b3.^2)

    out = n1.^2 + n2.^2 - p1.^2 * d.^2