Computing Square Corners

2021-02-25

How to compute square corner vertices? This post looks into two different approaches where the first method goes on to apply matrix multiplication within a algorithmic context that results in a C++ application. The second approach develops the idea of the first method by decomposing rotation matrix into its eigenvalue counterpart, which is then utilized in two different algorithms and in a second C++ program.

The question at hand can be formalized as follows.

Assume a square \(S\) with a center point at origin of a coordinate system as is shown in Figure 1,

Square at its center point

and a requirement that states;

Requirement 1: Define algorithm for computing corner vertices \(P_n \in \mathbb{R}^2\) of square \(S\) having a center point at a given \(O\), when \(n\) is in the sequence \(\langle\) 1 \(,\) 2 \(,\) 3 \(,\) 4 \(\rangle\).

The requirement can be satisfied as follows.

Method 1: Applying the Rotation Matrix

Say that the side length of \(S\) is \(2l\). Then the first corner vertex, can be thought of as a vector

as is shown in the next figure.

Graph of the vector   \(P\)

Next coordinate, \(P\), can then be derived from \(P\) by multiplying it with the rotation matrix

angle of the two vertices \(P\) and \(P\) is \(\frac{\pi}{2}\). In particular, \(P\) \(=RP\), as is shown below

Graph of transformation   \(RP\)

The statement \(P\) \(=RP\) can be generalized as \(P_n = RP_{n-1}\) if positive integer \(n\) is in the range \[ \begin{matrix} \left| \left(\pi - \frac{\pi}{4} \right)-\pi\right| \cdot \frac{4}{\pi} \ = \\[1em] 1 \ \leq \ n \ \leq \ 4 \\[1em] = \ \left| \left(\pi - \frac{4\pi}{4} \right)-\pi\right| \cdot \frac{4}{\pi}, \end{matrix} \] \(\theta\) is as above, and the side length of the square is known.

Furthermore, equation \(P_n = RP_{n-1}\) can now be expressed as

that fulfills Requirement 1.

The algorithm initializes with rotation matrix \(R\), having the given right angle, and a array containing the first corner coordinate which can be determined from the argument \(S\), as is shown at line \(1\).

After this, a \(\text{\textbf{while}}\)-loop begins, that iterates from \[ 3\frac{\pi}{4} \ \rightarrow \ 2\frac{\pi}{4} \ \rightarrow \ \frac{\pi}{4} \ \rightarrow \ 0. \tag{1} \] This means that the examined loop iterates three times and stops when there is no more \(\pi\) left at the line \(2\). That is, the outermost loop invariant reaches its termination condition when \(i \rightarrow 0\) just after the first three values of sequence \((1)\). The values are then converted in to indexor value sequence \(\langle 1,2,3 \rangle\) in the first index at line \(5\), where these values designate any current index in \(A\) where to place any coordinate \(P_n\) of \(S\) under iteration \(n\).

Additionally, the \(A.\)Prev \(A[ \ \left( |i - \pi| \cdot \frac{4}{\pi} \right) - 1 \ ]\). Or in other words, a member function of \(S\) that returns \(P_{n-1}\) in \(A\).

Rest of the routine from lines \(3\) to \(5\) are for computing the matrix transform \(RP\). When the loop terminates, solution array \(A\) is returned with computed corner coordinates of \(S\).

Is recommended to replace the lines \(\footnotesize 2\) and \(\footnotesize 5\) of \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{1}}\) with lines

\( \quad \footnotesize \text{\scriptsize 2} \ \ \ \text{\textbf{for }} i = 1 \text{\textbf{ to }} A.\text{length} - 1 \text{\textbf{ do}} \)

and

\( \quad \footnotesize \text{\scriptsize 5} \ \quad \ \quad \ \quad \ A[i][k] \mathrel{+}= R[k][j] \cdot A[k-1][j] \)

in a production code in order to ease the readability and maintainability.

A possible implementation of \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{1}}\) is given in the following


\(\text{\small\textbf{Source}}\) \(\text{\small\textbf{Code}}\) \(\text{\small\textbf{Listing}}\) \(\text{\small\textbf{1}}\)

std::vector<Vertex> verticesOf(Square const& S)
{
    using namespace std;
    using namespace std::numbers;
    using rotation_matrix = vector<vector<double>>;

    double const theta = pi / 2;
    double const l     = S.side_len * 0.5;
    double i           = pi;

    auto const n = [](double const& i) {
        return static_cast <unsigned> (abs(i - pi) * (4 / pi));
    };

    rotation_matrix const R = {
        {  cos(theta), sin(theta) },
        { -sin(theta), cos(theta) }
    };

    vector<Vertex> A = {
        {l,l}, {0,0}, {0,0}, {0,0}
    };

    while (i -= pi / 4)
        for (unsigned j = 0; j < A[0].size() - 1; j++)
            for (unsigned k = 0; k < A[0].size() - 1; k++)
                A[n(i)][j] += R[j][k] * A[n(i) - 1][k];

    return A;
}

for some Vertex and Square. Notice that \(\text{\small \textbf{Source Code Listing 1}}\) uses of the new <numbers> [See cppreference] header from the C++20 language standard (ISO/IEC, 2020) for the value of \(\pi\).

Taking account of the recommendations from above, the line 24 could be written as

24 for (unsigned i = 0; i < A.size() - 1; i++)

and the line 27 as

27 A[i][j] += R[j][k] * A[i - 1][k];

Doing this allows discarding variable i and the lambda n of lines 9 and 11-13, and indeed eases the readability and maintainability by reducing the implementation into 25 LOC, if the \(\text{\textbf{\small Source Code Listing 1}}\) is the chosen implementation of \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{1}}\).


Method 2: Applying the Rotation Matrix with a twist

The matrix multiplication \(RP_{n-1}\) can be expressed differently.

If \[ P_n = R(\theta)P_{n-1} = \lambda P_{n-1}, \tag{2} \] for some \(\lambda\), is multiplied by \(P_{n-1}^{-1}\), then \[ R(\theta) = \lambda, \] and so \[ R(\theta) - \lambda = 0. \] Moreover, multiplying the rotation matrix \(R\) and the proposed \(\lambda\) with identity matrix \(\mathbf{I}\) implies that \[ R(\theta) - \lambda\mathbf{I} = 0. \] Then the \[ \begin{aligned} &\text{det}(R(\theta) - \lambda \mathbf{I}) \\[3pt] &= \text{det} \begin{bmatrix} \cos \theta - \lambda & - \sin \theta \\ \sin \theta & \sin \theta - \lambda \\ \end{bmatrix} \\[1em] &= (\cos \theta - \lambda)^2 + \sin^2 \theta \\[3pt] &= \theta^2 - 2 \lambda \cos \theta + 1 \\[3pt] &= 0, \end{aligned} \] where from it can be proven that \[ \begin{aligned} \lambda &= \cos \theta \pm \sqrt{\cos^2 \theta - 1} \\ &= \cos \theta \pm \sin \theta i \\ &= \exp \pm \theta i, \end{aligned} \] as is shown by Haber (2019).

Hence, continuing with the \(\theta\) of previous method, it follows from equation \((\)\(2\)\()\), that \[ P_n = \lambda P_{n-1} = \left( \exp \pm \frac{\pi}{4} i \right) P_{n-1}. \] Therefore it's possible to refactor the matrix multiplication logic of the \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{1}}\) by modifying it with the newfound information concerning the eigenvalue \(\lambda\).

But before this can be done, few more insights are required.

First off, observe the following \(\mathbf{while}\)-loop of \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{1}}\);

\( \footnotesize \quad {\scriptsize 2} \ \quad \text{\textbf{while}} \ (i = \frac{\pi}{4}) \leftarrow (\pi - i) \ \text{\textbf{do}} \) \((3)\)

where the indexor \(i\) can now be thought of as being a complex number due to the observed imaginary term in \(\lambda\). If the indexor is then initialized with \(-\sqrt{-1}\) and iterated with multiplication by \(\sqrt{-1}\) until \(i \rightarrow 1\), the loop \((3)\) can be expressed as

\( \footnotesize \quad \text{\scriptsize 2}' \quad \text{\textbf{while}} \ (i = -\sqrt{-1}) \leftarrow (i \cdot \sqrt{-1}) \neq 1 \ \text{\textbf{do}} \) \((4)\)

or even

\( \footnotesize \quad \text{\scriptsize 2}' \quad \text{\textbf{while not}} \ (i = -\sqrt{-1}) \leftarrow (i \cdot \sqrt{-1}) \ \text{\textbf{do}} \) \((5)\)

where the iteration goes from

\[ i \ \rightarrow \ -1 \ \rightarrow \ -i \ \rightarrow \ 1, \tag{6} \]

In addition, the given sequence can be used to drive the computation of corner coordinates of square \(S\) as is shown in the second observation.

By inspecting the evaluation of the positive conjugate of \(\lambda\), a constant \[ \exp \left( \frac{\pi\sqrt{-1}}{4}\right) = 0.7\ldots + 0.7\ldots i = c' \tag{7} \] is revealed, which can be multiplied by some complex number \(i\) of the sequence \((6)\) generated in the pair of loops above, in order to derive a signed complex scalar constant \(c'i\) for computing \(x\) and \(y\) values of a corner in \(S\).

Each multiplication by \(c'i\) results in a complex number that uncovers magnitudes of some coordinate in \(S\), when the imaginary part in \((7)\) is ignored. However, the constants \(0.7\ldots\) of \(c'\) are necessary to be set to \(1\) — this can be thought of as a squaring the circle step taken at design time — so that the product \(c'i\) does not undershoot by a fraction of \(0.3\ldots\) during the corner computations, as is demonstrated in the next figure

Squaring the circle step visualized

where the orange arrow indicates squaring of the circle for \(S\). Therefore the revised multiplier constant simplifies to \[ c = 1 + i. \tag{8} \]

Now, given the transformed \(\text{\textbf{while}}\)-loop \((4)\) or \((5)\), and the new multiplication constant \((8)\), the \(x\) and \(y\) components of a chosen corner coordinate \(P\) of \(S\) can be computed by including the following lines

\( \footnotesize \begin{aligned} \quad {\scriptsize 3}' \quad P.x \leftarrow l \cdot \operatorname{Re} \ ci\\ \quad {\scriptsize 4}' \quad P.y \leftarrow l \cdot \operatorname{Im} \ ci \end{aligned} \)

into the new reformumlation of the current algorithm under developmet, when \(l = S.\text{side\_len} \cdot \frac{1}{2}\).

This is shown in the redesign of \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{1}}\) into

Suffice to say, this algorithm can be developed further.

If the scaling and assignment lines \(3\) and \(4\) are the done in-place of the Insert method call of \(A\), the initialization of the vector \(P\) at line \(1\) becomes redundant, because the vector's purpose is simply to function as a temporary variable for holding data during the iteration of the algorithm.

Since it is good practice not to clutter code with unnecessary variables in general, \(P\) should be refactored away with a more sufficient design, such as the instruction descriped;

\( \footnotesize \quad \text{\scriptsize 5}' \quad \ A. \)Insert\( \footnotesize ( \ l \cdot (\operatorname{Re} c, \operatorname{Im} c ) \ ) \)

More over, it's evident that the \(\text{\textbf{while}}\)-loops \((4-5)\)discussed earlier are or can be ambiguous, but fortunately an inquiry into them can open a possibility for a new formalization.

Main logic of the \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{2}}\)'s loop is to multiply \(c\) with a updating \(i\) that iterates trough its characteristic sequence \((6)\) in order to generate corner coordinatorial complex numbers where from to extract the coordinate information throughout the iteration. On the other hand, it is equally valid to update the value of \(c\) and hold \(i\) constant during the evaluation of \(ci\).

When this idea is applied trough the usage of a \(\text{\textbf{for}}\)-loop construct, it is possible to evaluate \(ci\) as

\( \footnotesize \quad {\scriptsize n} \ \ \textbf{for} \ c \leftarrow ci \ \textbf{to} \ \bar{c} \\ \)

at some line \(n\).

Depending on aesthetic inclinations, the line \(n\) could also be expressed more compactly as

\( \footnotesize \quad {\scriptsize n} \ \ \textbf{for} \ c \mathrel{\cdot}= i \ \textbf{to} \ \bar{c} \)

where a the \(\mathrel{\cdot}=\) operator from \(\mathbb{C} \times \mathbb{C} \rightarrow \mathbb{C}\) can be thought of as having the meaning \[ \mathrel{\cdot}=(i, j) = (i \leftarrow i \cdot j). \]

In any case, the variable \(c\) is immediately multiplied with \(i\) to get the first value of the product \(ci\) where from to read the required information. Note that the iteration now just takes three steps, which is acceptable because array \(A\) can be initialized with the first corner value \((l,l)\) since the \(l\) is known at the beginning of \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{2}}\).

This loop design of lines \(n\) above reverts the loop logic such that the line \(5' = n + 1\), so that the evaluation of \(S\)' corners happen a posteriori with respect to the main loop driving the computation — which is on contrary to the previous algorithm, where the evaluation is done a priori.

Putting together the ideas considered, a consequent algorithm can be structured;

that considerably clarifies the previous \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{2}}\) and fulfills the Requirement 1.

Likewise this algorithm can be expressed in C++ as


\(\text{\small\textbf{Source}}\) \(\text{\small\textbf{Code}}\) \(\text{\small\textbf{Listing}}\) \(\text{\small\textbf{2}}\)

auto verticesOf_(Square const& S)
{
    std::complex<double> c {1.0, 1.0}, i {0.0, 1.0}, c_ = c;
    double const l = S.side_len * 0.5;
    auto A = std::make_unique <std::vector<Vertex>> (
        std::initializer_list <Vertex> ({{l , l}})
    );

    for ( c *= i; c != c_; c *= i )
        A->emplace_back( Vertex{ l * c.real(), l * c.imag() } );

    return A.release();
}


Conclusion

This essay began with a question on how to compute a corner vertices of a given square, and formalized it into a Requirement 1. From there on the problem was inspected on and two solutions were presented along with example C++ implementaitons.

Reasons for the study include, but are not limited to,

  • study the problem in order to utilize possible new ideas in the future,
  • interest for learning to develop software on a scientific basis,
  • writing,
  • develop the site,
  • and learn and use common web-development tools.
Future post will present at least one more solution that simplifies \(\small \text{\textbf{Algorithm}}\) \(\small \text{\textbf{3}}\) further.

-Matti


References:

Haber, H. (2019). Eigenvalues and eigenvectors of rotation matrices. [PDF]. Physics 116A. Retrieved: 17 February 2021, from http://scipp.ucsc.edu/~haber/ph116A/Rotation2.pdf

ISO/IEC. (2020). International Standard ISO/IEC 14882:2020(E) – Programming Language C++. Geneva, Switzerland: International Organization for Standardization (ISO). Retrieved: 23 February 2021, from https://isocpp.org/std/the-standard

Weisstein, E. W. (n.d.). Circle Squaring. Retrieved: 16 February 2021, from https://mathworld.wolfram.com/CircleSquaring.html