import math


def simplify_U(theta, phi, lam):
    """Return the gate u1, u2, or u3 implementing U with the fewest pulses.

    U(theta, phi, lam) is the input gate.

    The returned gate implements U exactly, not up to a global phase.

    Return (gate_string, params, "OpenQASM string") where gate_string is one of
    "u1", "u2", "u3", "id" and params is a 3-tuple of parameter values. The
    OpenQASM string is the name of the gate with parameters substituted.
    """
    epsilon = 1e-13
    name = "u3"
    params = (theta, phi, lam)
    qasm = "u3(%.15f,%.15f,%.15f)" % params
    # Y rotation is 0 mod 2*pi, so the gate is a u1
    if abs(params[0] % (2.0 * math.pi)) < epsilon:
        name = "u1"
        params = (0.0, 0.0, params[1] + params[2] + params[0])
        qasm = "u1(%.15f)" % params[2]
    # Y rotation is pi/2 or -pi/2 mod 2*pi, so the gate is a u2
    if name == "u3":
        # theta = pi/2 + 2*k*pi
        if abs((params[0] - math.pi / 2) % (2.0 * math.pi)) < epsilon:
            name = "u2"
            params = (math.pi / 2, params[1],
                      params[2] + (params[0] - math.pi / 2))
            qasm = "u2(%.15f,%.15f)" % (params[1], params[2])
        # theta = -pi/2 + 2*k*pi
        if abs((params[0] + math.pi / 2) % (2.0 * math.pi)) < epsilon:
            name = "u2"
            params = (math.pi / 2, params[1] + math.pi,
                      params[2] - math.pi + (params[0] + math.pi / 2))
            qasm = "u2(%.15f,%.15f)" % (params[1], params[2])
    # u1 and lambda is 0 mod 4*pi so gate is nop
    if name == "u1" and abs(params[2] % (4.0 * math.pi)) < epsilon:
        name = "id"
        params = (0.0, 0.0, 0.0)
        qasm = "id"
    return name, params, qasm