PredefinedCurves

VortexPasta.PredefinedCurves Module

PredefinedCurves

Module defining some commonly used parametric curves.

source

Exported curves

VortexPasta.PredefinedCurves.Ring Type

Ring(; r = Returns(0), z = Returns(0))

Describes a circular ring in 3D space.

By default, the ring has radius $R=1$, it is defined on the XY plane, and its centre is located at (0, 0, 0).

The circular ring can be optionally perturbed in the radial and orthogonal (Z) directions using user-defined functions $r(t)$ and $z(t)$ respectively. These functions must be periodic with period 1. For example, one can pass r(t) = 0.1 * cos(2πmt) for a mode $m$ sinusoidal perturbation (associated to a wavelength $\lambda =2\pi R/m$) whose amplitude is 10% of the ring radius.

Moreover, the geometry can be modified via the optional transform and translate arguments of the define_curve function.

Definition

$$\begin{array}{rl}x(t)& =(1+r(t))\cos (2\pi t)\\ y(t)& =(1+r(t))\sin (2\pi t)\\ z(t)& =z(t)\end{array}$$

for $t\in [0,1]$.

Examples

Define a circular ring of radius 2:

julia> p = Ring();

julia> S = define_curve(p; scale = 2);  # scale the curve by 2 (to get radius = 2)

julia> ts = range(0, 1; length = 16 + 1);

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [2.0, 0.0, 0.0]
 [1.8477590650225735, 0.7653668647301796, 0.0]
 [1.4142135623730951, 1.414213562373095, 0.0]
 [0.7653668647301796, 1.8477590650225735, 0.0]
 [0.0, 2.0, 0.0]
 [-0.7653668647301796, 1.8477590650225735, 0.0]
 [-1.4142135623730951, 1.414213562373095, 0.0]
 [-1.8477590650225735, 0.7653668647301796, 0.0]
 [-2.0, 0.0, 0.0]
 [-1.8477590650225735, -0.7653668647301796, 0.0]
 [-1.4142135623730951, -1.414213562373095, 0.0]
 [-0.7653668647301796, -1.8477590650225735, 0.0]
 [0.0, -2.0, 0.0]
 [0.7653668647301796, -1.8477590650225735, 0.0]
 [1.4142135623730951, -1.414213562373095, 0.0]
 [1.8477590650225735, -0.7653668647301796, 0.0]
 [2.0, 0.0, 0.0]

Define a circular ring of radius 4 with a radial perturbation of 10%:

julia> rfun(t) = 0.1 * cospi(4t);  # this is a mode m = 2 perturbation with 10% amplitude

julia> p = Ring(r = rfun);

julia> S = define_curve(p; scale = 4);

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [4.4, 0.0, 0.0]
 [3.9568307230204227, 1.6389729494895986, 0.0]
 [2.8284271247461903, 2.82842712474619, 0.0]
 [1.4224945094311199, 3.4342055370698716, 0.0]
 [0.0, 3.6, 0.0]
 [-1.4224945094311199, 3.4342055370698716, 0.0]
 [-2.8284271247461903, 2.82842712474619, 0.0]
 [-3.9568307230204227, 1.6389729494895986, 0.0]
 [-4.4, 0.0, 0.0]
 [-3.9568307230204227, -1.6389729494895986, 0.0]
 [-2.8284271247461903, -2.82842712474619, 0.0]
 [-1.4224945094311199, -3.4342055370698716, 0.0]
 [0.0, -3.6, 0.0]
 [1.4224945094311199, -3.4342055370698716, 0.0]
 [2.8284271247461903, -2.82842712474619, 0.0]
 [3.9568307230204227, -1.6389729494895986, 0.0]
 [4.4, 0.0, 0.0]

source

VortexPasta.PredefinedCurves.TrefoilKnot Type

TrefoilKnot()

Describes a trefoil knot in 3D space.

Definition

Before transformations, the trefoil knot is defined by:

$$\begin{array}{rl}x(u)& =\sin u+2\sin 2u\\ y(u)& =\cos u-2\cos 2u\\ z(u)& =-\sin 3u\end{array}$$

for $u=2\pi t\in [0,2\pi ]$.

source

VortexPasta.PredefinedCurves.Lemniscate Type

Lemniscate(; Az = 0,)

Describes a lemniscate (or figure-eight curve) in 3D space.

The specific definition used here corresponds to the lemniscate of Bernoulli on the XY plane. Under the normalised parametrisation $\mathit{x}(t)$ with $t\in [0,1]$, the curve crosses the origin at $t=1/4$ and $t=3/4$, intersecting itself at that point.

One can perturb the curve in the third direction by passing a non-zero value of Az so that the curve doesn't exactly intersect itself.

Definition

$$\begin{array}{rl}x(u)& =\frac{\cos u}{1+{\sin }^{2}u}\\ y(u)& =\frac{\sin u\cos u}{1+{\sin }^{2}u}\\ z(u)& =A_z\sin u\end{array}$$

for $u=2\pi t\in [0,2\pi ]$.

source

VortexPasta.PredefinedCurves.PeriodicLine Type

PeriodicLine(; x = Returns(0), y = Returns(0), r = Returns(0))

Describes an infinite (but unclosed) line which repeats itself periodically.

Definition

By default, this defines a straight infinite line passing through the origin and oriented in the Z direction.

The line can be perturbed along the X and/or Y directions using user-defined functions $x(t)$ and $y(t)$ which must be periodic with period $T=1$. Alternative, one can pass a complex function $r(t)=x(t)+im\ast y(t)$, which may be more convenient in some cases.

The line is defined by:

$$\begin{array}{rl}x(t)& =x(t)+\Re [r(t)]\\ y(t)& =y(t)+\Im [r(t)]\\ z(t)& =t-1/2\end{array}$$

for $t\in [0,1]$. Note that the line is "centred" at $z=0$, in the sense that $z(1/2)=0$.

Note that one can change the default period $T=1$ of the curve using the scale argument of define_curve. See below for some examples.

Examples

Define a $2\pi$-periodic curve with a sinusoidal perturbation along X:

julia> xfun(t) = 0.1 * sinpi(2t)  # this function satisfies having period T = 1
xfun (generic function with 1 method)

julia> p = PeriodicLine(x = xfun);

julia> S = define_curve(p; scale = 2π);  # we scale the curve by 2π (this also scales the perturbation amplitude!)

julia> ts = range(0, 1; length = 16 + 1);

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [0.0, 0.0, -3.141592653589793]
 [0.24044709195373853, 0.0, -2.748893571891069]
 [0.4442882938158367, 0.0, -2.356194490192345]
 [0.5804906304278862, 0.0, -1.9634954084936207]
 [0.6283185307179586, 0.0, -1.5707963267948966]
 [0.5804906304278862, 0.0, -1.1780972450961724]
 [0.4442882938158367, 0.0, -0.7853981633974483]
 [0.24044709195373853, 0.0, -0.39269908169872414]
 [0.0, 0.0, 0.0]
 [-0.24044709195373853, 0.0, 0.39269908169872414]
 [-0.4442882938158367, 0.0, 0.7853981633974483]
 [-0.5804906304278862, 0.0, 1.1780972450961724]
 [-0.6283185307179586, 0.0, 1.5707963267948966]
 [-0.5804906304278862, 0.0, 1.9634954084936207]
 [-0.4442882938158367, 0.0, 2.356194490192345]
 [-0.24044709195373853, 0.0, 2.748893571891069]
 [0.0, 0.0, 3.141592653589793]

Note that the amplitude of the resulting curve is not the original $0.1$ but instead $0.1\times 2\pi =0.6823\dots$. This is because passing scale = 2π scales the curve in all directions, including the direction of the perturbation.

Instead, if one wanted to keep the original perturbation while extending the curve period from $1$ to $2\pi$, one would need to apply scaling only along the Z direction:

julia> using StaticArrays: SDiagonal

julia> S = define_curve(p; scale = SDiagonal(1, 1, 2π));

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [0.0, 0.0, -3.141592653589793]
 [0.03826834323650898, 0.0, -2.748893571891069]
 [0.07071067811865477, 0.0, -2.356194490192345]
 [0.09238795325112868, 0.0, -1.9634954084936207]
 [0.1, 0.0, -1.5707963267948966]
 [0.09238795325112868, 0.0, -1.1780972450961724]
 [0.07071067811865477, 0.0, -0.7853981633974483]
 [0.03826834323650898, 0.0, -0.39269908169872414]
 [0.0, 0.0, 0.0]
 [-0.03826834323650898, 0.0, 0.39269908169872414]
 [-0.07071067811865477, 0.0, 0.7853981633974483]
 [-0.09238795325112868, 0.0, 1.1780972450961724]
 [-0.1, 0.0, 1.5707963267948966]
 [-0.09238795325112868, 0.0, 1.9634954084936207]
 [-0.07071067811865477, 0.0, 2.356194490192345]
 [-0.03826834323650898, 0.0, 2.748893571891069]
 [0.0, 0.0, 3.141592653589793]

Now the perturbation amplitude is $0.1$ as we wanted.

source

Functions

VortexPasta.PredefinedCurves.define_curve Function

define_curve(
    p::ParametricCurve;
    scale = 1, rotate = LinearAlgebra.I, translate = 0,
    orientation::Int = 1,
) -> Function

Return the definition of the parametric curve p as a function.

The returned function S(t) returns a coordinate $\mathit{x}$ for any given value of the scalar parameter $t\in [0,1]$. In particular, closed curves satisfy S(0) == S(1).

Optional arguments

Coordinate transformations

The original curve can be transformed by (1) scaling, (2) rotation and (3) translation operations. Note that transformations are applied in that order.

• scale: scales the curve by a given factor. The argument can be a scalar value for isotropic scaling (same scaling in all directions) or a Tuple of values for anisotropic scaling (see below for some examples).

• rotate: in 3D, this should be a 3×3 orthogonal matrix describing pure rotation. For convenience, one can use the Rotations.jl package for defining such rotations using different parametrisations (see below for some examples).

• translate: a scalar or a vector describing a translation operation. Note that translations are performed after scaling and rotation.

Curve orientation

• orientation: allows to set the curve orientation. In particular, this determines the orientation of the tangent vector along the curve. Should be either 1 or -1.

Extended help

Examples

Translated and scaled circular ring

Define a circular ring of radius $R=2$ centred at $\mathit{x}₀=(0,0,1)$ and evaluate its coordinates over equispaced points.

julia> S = define_curve(Ring(); translate = (0, 0, 1), scale = 2);

julia> ts = range(0, 1; length = 16 + 1)
0.0:0.0625:1.0

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [2.0, 0.0, 1.0]
 [1.8477590650225735, 0.7653668647301796, 1.0]
 [1.4142135623730951, 1.4142135623730951, 1.0]
 [0.7653668647301796, 1.8477590650225735, 1.0]
 [0.0, 2.0, 1.0]
 [-0.7653668647301796, 1.8477590650225735, 1.0]
 [-1.4142135623730951, 1.4142135623730951, 1.0]
 [-1.8477590650225735, 0.7653668647301796, 1.0]
 [-2.0, 0.0, 1.0]
 [-1.8477590650225735, -0.7653668647301796, 1.0]
 [-1.4142135623730951, -1.4142135623730951, 1.0]
 [-0.7653668647301796, -1.8477590650225735, 1.0]
 [0.0, -2.0, 1.0]
 [0.7653668647301796, -1.8477590650225735, 1.0]
 [1.4142135623730951, -1.4142135623730951, 1.0]
 [1.8477590650225735, -0.7653668647301796, 1.0]
 [2.0, 0.0, 1.0]

Anisotropic scaling: ellipse from circular ring

If one wants an ellipse instead of a circle, one can simply apply an anisotropic scaling transformation:

julia> using LinearAlgebra

julia> S = define_curve(Ring(); scale = (2, 1, 1));

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [2.0, 0.0, 0.0]
 [1.8477590650225735, 0.3826834323650898, 0.0]
 [1.4142135623730951, 0.7071067811865476, 0.0]
 [0.7653668647301796, 0.9238795325112867, 0.0]
 [0.0, 1.0, 0.0]
 [-0.7653668647301796, 0.9238795325112867, 0.0]
 [-1.4142135623730951, 0.7071067811865476, 0.0]
 [-1.8477590650225735, 0.3826834323650898, 0.0]
 [-2.0, 0.0, 0.0]
 [-1.8477590650225735, -0.3826834323650898, 0.0]
 [-1.4142135623730951, -0.7071067811865476, 0.0]
 [-0.7653668647301796, -0.9238795325112867, 0.0]
 [0.0, -1.0, 0.0]
 [0.7653668647301796, -0.9238795325112867, 0.0]
 [1.4142135623730951, -0.7071067811865476, 0.0]
 [1.8477590650225735, -0.3826834323650898, 0.0]
 [2.0, 0.0, 0.0]

Rotated circular ring

If we wanted a circular ring defined on the YZ plane instead of the default XZ plane, we can achieve this by rotating the original curve by 90° about the Y axis.

julia> using Rotations: RotY  # there's also RotX and RotZ

julia> rot = RotY(π / 2)  # rotation of 90° about the Y axis
3×3 RotY{Float64} with indices SOneTo(3)×SOneTo(3)(1.5708):
  6.12323e-17  0.0  1.0
  0.0          1.0  0.0
 -1.0          0.0  6.12323e-17

julia> rot = SMatrix(replace(x -> abs(x) < 1e-16 ? zero(x) : x, rot))  # remove spurious near-zero values
3×3 SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
  0.0  0.0  1.0
  0.0  1.0  0.0
 -1.0  0.0  0.0

julia> S = define_curve(Ring(); scale = 2, rotate = rot);

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [0.0, 0.0, -2.0]
 [0.0, 0.7653668647301796, -1.8477590650225735]
 [0.0, 1.4142135623730951, -1.4142135623730951]
 [0.0, 1.8477590650225735, -0.7653668647301796]
 [0.0, 2.0, 0.0]
 [0.0, 1.8477590650225735, 0.7653668647301796]
 [0.0, 1.4142135623730951, 1.4142135623730951]
 [0.0, 0.7653668647301796, 1.8477590650225735]
 [0.0, 0.0, 2.0]
 [0.0, -0.7653668647301796, 1.8477590650225735]
 [0.0, -1.4142135623730951, 1.4142135623730951]
 [0.0, -1.8477590650225735, 0.7653668647301796]
 [0.0, -2.0, 0.0]
 [0.0, -1.8477590650225735, -0.7653668647301796]
 [0.0, -1.4142135623730951, -1.4142135623730951]
 [0.0, -0.7653668647301796, -1.8477590650225735]
 [0.0, 0.0, -2.0]

More generally, to rotate about an arbitrary axis ê = [ex, ey, ez] by an angle θ (in radians), one can use rot = AngleAxis(θ, ex, ey, ez) from the Rotations.jl package.

Random rotations

In addition, the Rotations.jl package allows to easily generate random and uniformly distributed rotations:

julia> using Rotations: QuatRotation  # parametrise rotations using quaternions

julia> using Random: MersenneTwister

julia> rng = MersenneTwister(42);

julia> rot = rand(rng, QuatRotation)  # uniformly distributed random rotation
3×3 QuatRotation{Float64} with indices SOneTo(3)×SOneTo(3)(QuaternionF64(0.923391, -0.0602724, 0.307861, 0.22122)):
  0.712566  -0.445655   0.541886
  0.371433   0.894858   0.24752
 -0.59522    0.0249001  0.803177

julia> S = define_curve(Ring(); scale = 2, rotate = rot);

julia> S.(ts)
17-element Vector{SVector{3, Float64}}:
 [1.4251329706813418, 0.7428666122301657, -1.190439083830269]
 [0.9755612628402096, 1.3712140745918031, -1.0807646293652284]
 [0.37746919621652514, 1.790806624183378, -0.806553557235094]
 [-0.27808913376434097, 1.9377650989455064, -0.4095520174421194]
 [-0.8913109140138611, 1.7897164032775457, 0.04980010440813476]
 [-1.3688386873583265, 1.3691996090301746, 0.501570611801321]
 [-1.637973179106087, 0.7402345861333226, 0.8769815402966743]
 [-1.6577411025987892, -0.00142444225909442, 1.1188799791393182]
 [-1.4251329706813418, -0.7428666122301657, 1.190439083830269]
 [-0.9755612628402096, -1.3712140745918031, 1.0807646293652284]
 [-0.37746919621652514, -1.790806624183378, 0.806553557235094]
 [0.27808913376434097, -1.9377650989455064, 0.4095520174421194]
 [0.8913109140138611, -1.7897164032775457, -0.04980010440813476]
 [1.3688386873583265, -1.3691996090301746, -0.501570611801321]
 [1.637973179106087, -0.7402345861333226, -0.8769815402966743]
 [1.6577411025987892, 0.00142444225909442, -1.1188799791393182]
 [1.4251329706813418, 0.7428666122301657, -1.190439083830269]

source