2
\$\begingroup\$

In order to plot a squircle using its parametric equations:

enter image description here

The equations can be found in page 9 of this this article:

Fong, Chamberlain. "Squircular Calculations." arXiv preprint arXiv:1604.02174 (2016).

Unfortunately the following singularities are not handled by those equations:

  • For t = 0, t = pi/2, t = pi, t = 3pi/2, and t=2*pi

Here is the code to plot these equations in Python

"""Plot squircle shape using its parametric equations"""
import numpy as np
from numpy import pi, sqrt, cos, sin, sign
import matplotlib.pyplot as plt

# FG stands for Fernandez-Guasti who introduced the algebric equations.
def FG_param(r, s, t):
    x = np.ones_like(t)
    y = np.ones_like(t)
    for i, _t in enumerate(t):
        if np.isclose(_t, 0.0):
            x[i] = r
            y[i] = 0.0
        elif np.isclose(_t, pi/2.0):
            x[i] = 0.0
            y[i] = r
        elif np.isclose(_t, pi):
            x[i] = -r
            y[i] = 0.0
        elif np.isclose(_t, 3*pi/2.0):
            x[i] = 0.0
            y[i] = -r
        elif np.isclose(_t, 2*pi):
            x[i] = r
            y[i] = 0.0
        else:
            x[i] = r*sign(cos(t[i]))/(s*sqrt(2)*np.abs(sin(t[i])))*sqrt(1 - sqrt(1-s**2*sin(2*t[i])**2))
            y[i] = r*sign(sin(t[i]))/(s*sqrt(2)*np.abs(cos(t[i])))*sqrt(1 - sqrt(1-s**2*sin(2*t[i])**2))
    return x, y

s = 0.7
r = 1.0
NTHETA = 90
t = np.linspace(0, 2*pi, NTHETA)

x, y = FG_param(r, s, t)

plt.gca().set_aspect('equal')
plt.scatter(x, y, s=5)
plt.show()

The function FG_param looks very clumsy. I appreciate any help to make it more compact and neat.

Thank you

\$\endgroup\$
  • \$\begingroup\$ It is possible to simplify via knowing that boolean statements return 0 or 1 if they are coerced to; for example, np.isclose(_t, pi/2.0) returns 1 if true, so that can be inverted (display zero if matched) and then multiplied by the RHS to remove it; the equation can be restructured to be r + xsin(...), where x is the condition. This can be repeated for all remaining conditions and structured accordingly. I may provide a code sample \$\endgroup\$ – alexyorke Dec 6 '19 at 4:39
  • \$\begingroup\$ Something sort of like this, but this is not 100% correct: x[i] = r * \ (np.isclose(_t, 0.0) or np.isclose(_t, 2*pi)) * \ (not np.isclose(_t, 0.0) or np.isclose(_t, 2*pi)) * \ (1 - (not np.isclose(_t, pi))) * \ (not np.isclose(_t, 3*pi/2.0)) \$\endgroup\$ – alexyorke Dec 6 '19 at 4:40
3
\$\begingroup\$

Here's a different solution. The singularities occur when sin(_t) or cos(_t) is zero. So check for those conditions.

def FG_param(r, s, t):
    x = np.ones_like(t)
    y = np.ones_like(t)

    for i, _t in enumerate(t):
        sin_t = sin(_t)
        cos_t = cos(_t)

        if np.isclose(sin_t, 0.0):
            x[i] = -r if np.isclose(_t, pi) else r
            y[i] = 0.0

        elif np.isclose(cos_t, 0.0):
            x[i] = 0.0
            y[i] = r if _t < pi else -r

        else:
            radical = sqrt(1 - sqrt(1 - s**2 * sin(2*_t)**2))

            x[i] = r*sign(cos_t) / (s*sqrt(2)*np.abs(sin_t)) * radical
            y[i] = r*sign(sin_t) / (s*sqrt(2)*np.abs(cos_t)) * radical

    return x, y
|improve this answer|||||
\$\endgroup\$
  • 1
    \$\begingroup\$ The singularities occur also when s=0 \$\endgroup\$ – Navaro Dec 6 '19 at 12:13
  • \$\begingroup\$ @Navaro: Exactly. In that case (i.e: s=0) the map will be a circle with radius r. \$\endgroup\$ – IamNotaMathematician Dec 7 '19 at 0:31
  • \$\begingroup\$ @Navaro s==0 wasn't handled in the original code, so I didn't handle it either. It would be an easy case to handle. \$\endgroup\$ – RootTwo Dec 7 '19 at 1:08
  • \$\begingroup\$ @RootTwo: Yes, you could just add something like this at the top of the function: if s == 0: return r*cos(t), r*sin(t) \$\endgroup\$ – Navaro Dec 7 '19 at 1:12
2
\$\begingroup\$

This solution is an order of magnitude faster than the other ones. It uses numpy functions, so there aren't any explicit loops. np.where chooses values that keep the equations from blowing up where the denominator would go to zero. Case for s == 0 is handled separately.

def FG_param1(r, s, t):

    if s == 0:
        return r*cos(t), r*sin(t)


    sin_t = sin(t)
    cos_t = cos(t)

    radical = sqrt(1 - sqrt(1 - s**2 * sin(2*t)**2))

    x_denom = np.where(np.isclose(sin_t, 0.0), 1.0, s*sqrt(2)*np.abs(sin_t))
    x = r * sign(cos_t) * np.where(np.isclose(sin_t, 0.0, 1e-12), 1.0, radical / x_denom)

    y_denom = np.where(np.isclose(cos_t, 0.0), 1.0, s*sqrt(2)*np.abs(cos_t))
    y = r * sign(sin_t) * np.where(np.isclose(cos_t, 0.0, 1e-12), 1.0, radical / y_denom)

    return x, y
|improve this answer|||||
\$\endgroup\$
1
\$\begingroup\$

Off the top of my head, and untested, how about:

def FG_param(r, s, t):
    x = np.ones_like(t)
    y = np.ones_like(t)

    for i, _t in enumerate(t):
        radical = sqrt(1 - sqrt(1 - s**2 * sin(2*_t)**2))

        try:
            x[i] = r*sign(cos(_t)) / (s*sqrt(2)*np.abs(sin(_t))) * radical

        except ZeroDivisionError:
            x[i] = -r if np.isclose(_t, pi) else r

        try:
            y[i] = r*sign(sin(_t)) / (s*sqrt(2)*np.abs(cos(_t))) * radical

        except ZeroDivisionError:
            y[i] = r if _t < pi/2.0 else -r

    return x, y
|improve this answer|||||
\$\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.