What kind of software can be used to solves systems of equations?

Question

For example, I have to solve the following equations:

$$\left\{\begin{align*} &x^2 + y^2 + z^2 = 1\\ &Ax + By + Cz = 0 \end{align*}\right.$$

for $y$ and $z$, where $x$, $A$, $B$ and $C$ are known values. I need this for production. I guess I could solve this by hand, but there are many other of those and, if did all manually, my production would stop. Is there a software that solves this kind of situation?

Answer 1

Hot off the Mathematica presses:

> Solve[{x^2 + y^2 + z^2 - 1 == 0, A x + B y + C z == 0}, {y, z}] // Simplify

$$\begin{align} \Bigg\{\Bigg\{y & \to -\frac{A B x+\sqrt{-C^2 \left(A^2 x^2+B^2 \left(-1+x^2\right)+C^2 \left(-1+x^2\right)\right)}}{B^2+C^2}, \\ z & \to \frac{-A C^2 x+B \sqrt{-C^2 \left(A^2 x^2+B^2 \left(-1+x^2\right)+C^2 \left(-1+x^2\right)\right)}}{C \left(B^2+C^2\right)}\Bigg\}, \\ \Bigg\{y & \to \frac{-A B x+\sqrt{-C^2 \left(A^2 x^2+B^2 \left(-1+x^2\right)+C^2 \left(-1+x^2\right)\right)}}{B^2+C^2}, \\ z & \to -\frac{A C^2 x+B \sqrt{-C^2 \left(A^2 x^2+B^2 \left(-1+x^2\right)+C^2 \left(-1+x^2\right)\right)}}{C \left(B^2+C^2\right)}\Bigg\}\Bigg\} \\ \end{align}$$

I have intentionally not simplified these equations, so as to give you a feel what Mathematica solutions are like.

Answer 2

Based on Mario Carneiro's comment, you can do just this:

Let $V = (A,B,C)$ be a vector. Your equation is asking for the set of points $(x,y,z)$ on the unit sphere that is orthogonal to $V$.

1. If $A = B = 0$ we have that the circle is just the circle in the $x, y$ plane. So you can parametrise as normal: $\theta \mapsto (\cos\theta, \sin\theta, 0)$ would give you points on this circle by angle.

2. If at least one of $A,B$ are not zero you can rotate the unit circle from the $x,y$ plane via an orthogonal transformation. We do this by constructing an orthogonal matrix that sends the $z$ axis to the subspace spanned by the vector $V$. Since orthogonal transformations preserve length and angles, the image of the unit circle from the $x,y$ plane under this orthogonal transformation would be the circle you want. There are lots of ways to do this. We choose one here that is algebraically simpler.

   Set the vectors $W_3 = V$, $W_2 = (B, -A, 0)$, $W_1 = (-AC,-BC,B^2+A^2)$. A direct computation will verify to you that the three are orthogonal. Set $U_3 = W_3 / \|W_3\|$ where $\|\cdot\|$ is the vector 2-norm, so $U_3 = W_3 / \sqrt{A^2 + B^2 + C^2}$. Set $U_2 = W_2 / \|W_2\|$ and set $U_1 = W_1 / \|W_1\|$. Now the matrix $O$ whose columns are given by the vectors $U_1, U_2, U_3$ is orthogonal. And it maps the $z$ axis to the subspace spanned by $U_3$, which is also the subspace spanned by $V$.

   So the following mapping would give you a parametrisation of the required circle: $$ \theta \mapsto \cos\theta U_1 + \sin\theta U_2 $$ for completeness we expand the definition of $U_1$ and $U_2$ to get $$ \theta \mapsto \cos\theta \frac{(-AC, -BC, A^2 + B^2)}{\sqrt{ (A^2 + B^2+ C^2)(A^2 + B^2)}} + \sin\theta \frac{(B,-A,0)}{\sqrt{A^2 + B^2}} $$ Note that when at least one of $A,B$ is non-zero, the denominators are non-zero and so the mapping is well defined.

One caveat: because of the method of construction, the family of maps is discontinuous in $(A,B,C)$ as $A,B$ both approach 0. This is why we have to fill in the edge case specifically in step 1.

---

In view of Dokkat's comment below, let me quickly describe the construction of $O$ and the choice of $W_1,W_2,W_3$. The first thing to remember from linear algebra is that

If $O$ is the matrix corresponding to an orthogonal transformation, then its column vectors form an orthonormal basis of the vector space, and vice versa.

Furthermore, if the column vectors of $O$ are $U_1, U_2, U_3$ in order, then the vector $(0,0,1)$ would be mapped to $U_3$ by the transformation $O$ by matrix multiplication.

So to get our desired $O$ (which sends the $z$ axis to the subspace spanned by $V$), we need that $U_3$ is the unit vector corresponding to $V$. Then we just need to find $U_1, U_2$ such that $U_1, U_2, U_3$ together is an orthonormal basis.

Given that the direction of $U_3$ is fixed, we first fix the directions of $U_1$ and $U_2$. After we fix the directions, we can just take the unit vectors in those directions to fix the length.

To fix the directions we observe the following two facts:

• Suppose $A \neq 0$, the vector $(A,B,C)$ is perpendicular to the vector $(B,-A,0)$. This is an analogue of "rotating by 90 degrees in the plane" when we think about vectors in a plane. Note that we keep only the $x,y$ components and drop the $z$ component. You can do this in also higher dimensional spaces. We set $W_2$ to be this vector.

• In 3D space, we have the cross product. Which has the nice property that given vectors $\vec{v}, \vec{w}$, their cross product $\vec{v}\times \vec{w}$ is perpendicular to either $\vec{v}$ or $\vec{w}$. Using this property we set $W_1 = W_2 \times W_3$. (Or maybe the negative of that, I forgot which sign I took.)

Then the $W_1, W_2, W_3$ will be mutually orthogonal nonvanishing vectors.