6
\$\begingroup\$

I recently wrote a function to calculate the intersection point(s) of two parabolas. Working it out on paper was a nightmare, but I eventually did it (I think it's correct), and ended up with:

\$ x = \frac{\sqrt{c_2 - c_1 + \frac{b_1^2}{4a_1} - \frac{b_2^2}{4a_2} + E}}{\sqrt{a_2 - a_1}} - \frac{b1 - b2}{2(a1 - a2)} \$

with

\$ E = \frac{b_1^2a_2^3(a_2 - 2a_1^2) - a_1^3b_2^2(a_1 - 2a_2^2) - 2a_1^2a_2^2b_2(b_2 - b_1)}{(4a_1^2a_2^2)(a_1 - a_2)^2} \$

For any two parabolas \$ a_1x^2 + b_1x + c_1 = 0\$ and \$a_2x^2 + b_2x + c_2 = 0\$.

That can't be the fastest way to do it, I must be missing something. Does anyone know of any better ways I could implement parabola to parabola intersection point calculations?

The function I wrote is below:

def parabola_to_parabola_poi(a1, b1, c1, a2, b2, c2):
    """
    Calculate the intersection point(s) of two parabolas.
    """
    a1_sqrd, a2_sqrd = a1*a1, a2*a2
    b1_sqrd, b2_sqrd = b1*b1, b2*b2
    c1_sqrd, c2_sqrd = c1*c1, c2*c2

    E1 = b1_sqrd*a2_sqrd*a2*(a2 - 2*a1_sqrd)
    E2 = a1_sqrd*a1*b2_sqrd*(a1 - 2*a2_sqrd)
    E3 = 2*a1_sqrd*a2_sqrd*b2*(b2 - b1)
    E4 = (4*a1_sqrd*a2_sqrd)*(a1 - a2)*(a1 - a2)

    E = (E1 - E2 - E3)/E4

    F = (c2 - c1 + b1_sqrd/(4*a1) - b2_sqrd/(4*a2) + E)/(a2 - a1)
    G = (b1 - b2)/(2*(a1 - a2))

    if not F:
        px = math.sqrt(F) - G
        py = a1*px*px + b1*px + c1
        return [
            (px, py)
            ]
    elif F < 0:
        return []
    sqrt_F = math.sqrt(F)

    px1 = sqrt_F - G
    px2 = -sqrt_F - G

    py1 = a1*px1*px1 + b1*px1 + c1
    py2 = a1*px2*px2 + b1*px2 + c1

    return [
        (px1, py1),
        (px2, py2)
        ]
\$\endgroup\$
4
\$\begingroup\$

Well it seems as though I SIGNIFICANTLY overcomplicated the equation for some reason. I just did the equation it a second time and I got something much nicer:

For any two parabolas \$ ax^2 + bx + c = 0\$ and \$dx^2 + ex + f = 0\$, we get:

\$ ax^2 - dx^2 + bx - ex + c - f = 0 \\ x^2(a - d) + x(b - e) = f - c \$

Now it's just a matter of completing the square.

\$ x^2(a - d) + x(b - e) + \frac{(b - e)^2}{4(a - d)} = f - c + \frac{(b - e)^2}{4(a - d)} \\ (x\sqrt{a - d} + \frac{b - e}{2\sqrt{a - d}})^2 = f - c + \frac{(b - e)^2}{4(a - d)} \\ (a - d)(x + \frac{b - e}{2(a - d)})^2 = f - c + \frac{(b - e)^2}{4(a - d)} \\ x + \frac{b - e}{2(a - d)} = \sqrt{\frac{f - c + \frac{(b - e)^2}{a - d}}{a - d}} \\ x = \sqrt{\frac{f - c + \frac{(b - e)^2}{a - d}}{a - d}} - \frac{b - e}{2(a - d)} \\ \$

I don't know how I ended up with the monster I originally posted, but this seems much more reasonable.

\$\endgroup\$
  • 2
    \$\begingroup\$ Even simpler: apply the quadratic formula to \$(a-d)x^2 + (b-e)x + (c-f)\$. \$\endgroup\$ – Anonymous May 18 '14 at 3:04
  • \$\begingroup\$ You lost a factor of four in the denominator of \$(b-e)^2\$ between the third-to-last equation and the second-to-last equation. \$\endgroup\$ – David K May 19 '14 at 13:12
  • \$\begingroup\$ Also missed a \$\pm\$ sign. \$\endgroup\$ – David K May 19 '14 at 13:34
4
\$\begingroup\$

Many of the other answers have pointed out that evaluation the Quadratic Forumula is a simpler way solve this equation.

Here is my solution:

import math

def find_intersect(a1, b1, c1, a2, b2, c2):
    a = a1-a2
    b = b1-b2
    c = c1-c2

    inner_calc = b**2 - 4*a*c

    # Check if `inner_cal` is negative. If so, there are no real solutions.
    # Thus, return the empty set.
    if inner_calc < 0:
        return set()

    square = math.sqrt(inner_calc)
    double_a = 2*a
    answers = [(-b + square)/double_a, (-b - square)/double_a]

    # Using `set()` removes any possible duplicates.
    return set(answers)

If you were interested in returning all possible answers (real and imaginary), you could use the cmath library which automatically handles negative numbers in its sqrt() function. All you need to change is remove the if-statement and import cmath instead of math:

import cmath

def find_intersect(a1, b1, c1, a2, b2, c2):
    a = a1-a2
    b = b1-b2
    c = c1-c2

    square = cmath.sqrt(b**2 - 4*a*c)
    double_a = 2*a
    answers = [(-b + square)/double_a, (-b - square)/double_a]

    # Using `set()` removes any possible duplicates.
    return set(answers)
\$\endgroup\$
  • \$\begingroup\$ Yes, you have to test b**2 - (4*a*c) for a negative value before calling math.sqrt. In the code above, you go ahead and take the square root (using cmath.sqrt), then do a bunch of other operations, and then in the end you still need a conditional (if result.imag == 0.0) and if it's false, you throw away the results of all your previous calculations. It seems more efficient to me to just stop after computing b**2 - (4*a*c) when the result is negative. \$\endgroup\$ – David K May 19 '14 at 17:26
  • \$\begingroup\$ Good catch. I have updated my answer, thanks! \$\endgroup\$ – BeetDemGuise May 19 '14 at 17:51
  • \$\begingroup\$ I'd still prefer d = b*b - 4*a*c followed by if d < 0:, and in the case where d >= 0, then compute the square root. \$\endgroup\$ – David K May 19 '14 at 18:10
  • \$\begingroup\$ In most trivial cases, I would tend to disagree with you. I would rather do an extra trivial operation than store another local variable. However, in this case, cmath.sqrt is probably not a trivial operation. Thanks for all the input! \$\endgroup\$ – BeetDemGuise May 19 '14 at 18:33
  • \$\begingroup\$ Yes, I suppose I tend to think of square root as a heavy-weight operation. If I really wanted to minimize operations (at the expense of extra variables) I could have stored the results of b*b and 4*a*c in temporary variables and following this with something like if temp1 < temp2 ... but that is so extreme that I didn't even think of it until you made the tradeoff clear. I may also have been influenced by the fact that \$b^2 - 4ac\$ is such an important quantity in math that we have a name for it (discriminant), which makes it easier to think of naming it something in the code. \$\endgroup\$ – David K May 19 '14 at 20:12
3
\$\begingroup\$

Given two parabolas \$ f_1(x) = a_1 x^2 + b_1 x + c_1\$ and \$ f_2(x) = a_2 x^2 + b_2 x + c_2\$, you are asking to find \$x\$ such that \$f_1(x) = f_2(x)\$.

That is, you want \$f_1(x) - f_2(x) = 0\$. But if you write this as a difference of two polynomials in the usual way, you are trying to solve

\[f_1(x) - f_2(x) = (a_1 - a_2) x^2 + (b_1 - b_2) x + (c_1 - c_2) = 0.\]

Apply the usual quadratic formula with \$ a = a_1 - a_2\$, \$ b = b_1 - b_2\$, and \$ c = c_1 - c_2\$. The result is

\[ x = \frac{-(b_1 - b_2) \pm \sqrt{(b_1 - b_2)^2 - 4(a_1 - a_2)(c_1 - c_2)}} {2(a_1 - a_2)}.\]

You could try to manipulate this further, but I don't think you'll get anything easier to compute. Note that there are either two, one, or no possible solutions for \$x\$, depending on the value of the expression under the square root. Also note that we could write \$b_2 - b_1\$ instead of \$-(b_1 - b_2)\$, but I think writing it this way is easier to check against the standard formula and is also suggestive of the best code to compute the result.

\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.