Tools for verifying equations with 6+ variables

Question

TLDR: looking for a tool to check equations with 6+ parameters

Hello, this may already have an answer, but here goes:
I am working on solving geometric equations for Ray-Ellipse intersections. a ray here would be a line that ends at a given point. the Equations I am using are as follows: $$Ellipse:\frac{(x-h)²}{a²} + \frac{(y - k)²}{b²} = 1$$ $$Ray:\begin{cases}x = o_1 + m v_1 \\\\y = o_2 + mv_2\end{cases}, m < 0$$ $$O = (o_1, o_2), V = \begin{bmatrix} v_1 \\\\ v_2 \end{bmatrix}$$ a = eccentricity in X, b = eccentricity in Y, h and k are distance from origin in X and Y.
To simplify the algebra: $$y = cx+d , c = \frac{v_2}{v_1}, d = o_2 - c*o_1$$ Solving for X and Y, I have found an equation I won't bother writing here(mostly due to the current version likely being wrong), I get an equation with 6 unknown parameters excluding x and y. And no, I cannot replace the parameters due to the equations being used in code with code-generated parameters.

I have already tried using Wolfram Mathworld, but the generated equations are useless for my application. Similarly, when using Geogebra to compare generated coordinates with the intersect function, The program acts up and is very hard to use.
Is there any alternative, or better way of doing this?
The equations need to have a degree of error that is as low as possible, due to the error from the calculations accumulating through tens and possibly hundreds of iterations.
Sorry for the text wall, and a good day to whoever reads this.

Answer 1

How do you have six unknowns, not including $x,y$? Do you not know what ellipse and what ray you are intersecting? If so, then what DO you know?

I am going to assume that you are misunderstanding what "unknown" means here. If a parameter is something that will be provided to you when you actually need to calculate a number, then it is a known parameter, even if you do not know right now what that parameter value will be.

If your goal is: given an ellipse and a ray, to find the points of intersection between the two, then $(h, k), a, b, (o_1, o_2), (v_1, v_2)$ are all known parameters. They are the inputs to your function, defining which ellipse and which ray you are intersecting. $t$ (which you label $m$, but I am changing it to something more traditional) is the unknown you need to solve for, and the actual intersection point(s) $(x,y)$ are obtained from the value(s) found for $t$.

You have the three equations $$\dfrac{(x - h)^2}{a^2} +\dfrac{(y - k)^2}{b^2} = 1\\x = o_1 + tv_1\\y = o_2 + tv_2$$ We can eliminate $x,y$ by substituting the values in the second and third equation into the first:

$$\dfrac{(v_1t + o_1 - h)^2}{a^2} + \dfrac{(v_2t + o_2 - k)^2}{b^2} = 1$$

We can simplify this first by dividing the denominators into each value: let $$V_1 = \frac{v_1}a, O_1 = \frac{h-o_1}a\\V_2 = \frac{v_2}b, O_2 = \frac{k-o_2}b$$ and the equation becomes $$(V_1t - O_1)^2 + (V_2t - O_2)^2 = 1\\ (V_1^2 + V_2^2)t^2 - 2(V_1O_1 + V_2O_2)t + (O_1^2 + O_2^2 - 1) = 0$$ Next define $$A = V_1^2 + V_2^2\\B = V_1O_1 + V_2O_2\\C = O_1^2 + O_2^2 - 1$$ and the equation becomes $$At^2 - 2Bt + C = 0$$ which we can throw into the quadratic formula to get $$t = \dfrac{2B \pm \sqrt{4B^2 - 4AC}}{2A} = \dfrac{B \pm \sqrt{B^2 - AC}}A$$

Now you do some checks:

• If $B^2 - AC < 0$, then the ellipse and the ray do not intersect.

Otherwise calculate the two values of $t$: $$t_1 = \frac BA + \frac{\sqrt{B^2 - AC}}A\\t_2 = \frac BA - \frac{\sqrt{B^2 - AC}}A$$

If $B^2 = AC$, there will be a single value: $t_1 = t_2$. This indicates that the line of the ray is tangent to the ellipse (though possibly not on the side of the ray). Otherwise there will be two values. Since the ray extends in the direction $t < 0$, if either of $t_1$ or $t_2$ is $\ge 0$, it should be discarded.

• if both $t_1, t_2 \ge 0$, there is no point of intersection.
• if both $t_1, t_2 < 0$, there will be two points of intersection if $t_1 \ne t_2$, or one point of tangent intersection if $t_1 = t_2$.
• if only one of $t_1, t_2$ is $< 0$, then $(o_1, o_2)$ is inside the ellipse and there is a single point of intersection.

For each valid value of $t = t_1, t_2$, you calculate the point of intersection as $(x,y) = (v_1t + o_1, v_2t+o_2)$

Here is some pseudocode:

class Vector {
   double X      \\x-coordinate
   double Y      \\y-coordinate
}

class Ray {
   Vector Origin 
   Vector AntiDirection

   function Vector PointAt(double t) {
      return new Vector {X = Origin.X + t * AntiDirection.X, Y = Origin.Y + t * AntiDirection.Y 
   }
}

class Ellipse {
   Vector Center       \\(h, k)
   double Xradius      \\a
   double Yradius      \\b
}

function Vector[] IntersectionPoints(Ray R, Ellipse E) {
   double xcoef = R.AntiDirection.X / E.Xradius
   double ycoef = R.AntiDirection.Y / E.Yradius
   double xoffset = (E.Center.X - R.Origin.X) / E.Xradius
   double yoffset = (E.Center.Y - R.Origin.Y) / E.Yradius
   double A = xcoef * xcoef + ycoef * ycoef
   double B = xcoef * yoffset + ycoef * xoffset
   double C = xoffset * xoffset + yoffset * yoffset - 1

   Vector[] OutList = new Vector[]
   
   double discriminant = B * B - A * C
   if (discriminant < 0) return OutList  \\empty list
   
   double t_1 = (B + Sqr(discriminant)) / A
   double t_2 = (B - Sqr(discriminant)) / A
   if (t_1 < 0) OutList.Add( R.PointAt(t_1) )
   if (t_2 < 0 and t_2 != t_1) OutList.Add( R.PointAt(t_2) )
   return OutList
}