OpenKnot

Class for working with open (linear) curves, that do not form closed loops. OpenKnot provides methods for visualising these curves and analysing their topology via different kinds of closures.

API documentation

class pyknotid.spacecurves.openknot.OpenKnot(*args, **kwargs)[source]

Bases: pyknotid.spacecurves.spacecurve.SpaceCurve

Class for holding the vertices of a single line that is assumed to be an open curve. This class inherits from SpaceCurve, replacing any default arguments that assume closed curves, and providing methods for statistical analysis of knot invariants on projection and closure.

All knot invariant methods return the results of a sampling over many projections of the knot, unless indicated otherwise.

alexander_fractions(number_of_samples=10, **kwargs)[source]

Returns each of the Alexander polynomials from self.alexander_polynomials, with the fraction of that type.

alexander_polynomials(number_of_samples=10, radius=None, recalculate=False, zero_centroid=False, optimise_closure=True)[source]

Returns a list of Alexander polynomials for the knot, closing on a sphere of the given radius, with the given number of sample points approximately evenly distributed on the sphere.

The results are cached by number of samples and radius.

Parameters:
  • number_of_samples (int) – The number of points on the sphere to sample. Defaults to 10.
  • optimise_closure (bool) – If True, doesn’t really close on a sphere but at infinity. This lets the calculation be optimised slightly, and so is the default.
  • radius (float) – The radius of the sphere on which to close the knot. Defaults to None, which picks 10 times the largest Cartesian deviation from 0. This is only used if optimise_closure=False.
  • zero_centroid (bool) – Whether to first move the average position of vertices to (0, 0, 0). Defaults to True.
Returns:

A number_of_samples by 3 array of angles and alexander polynomials.

Return type:

ndarray

alexander_polynomials_multiroots(number_of_samples=10, radius=None, zero_centroid=False)[source]

Returns a list of Alexander polynomials for the knot, closing on a sphere of the given radius, with the given number of sample points approximately evenly distributed on the sphere. The Alexander polynomials are found at three different roots (2, 3 and 4) and a the knot types corresponding to these roots are returned also.

The results are cached by number of samples and radius.

Parameters:
  • number_of_samples (int) – The number of points on the sphere to sample. Defaults to 10.
  • radius (float) – The radius of the sphere on which to close the knot. Defaults to None, which picks 10 times the largest Cartesian deviation from 0.
  • zero_centroid (bool) – Whether to first move the average position of vertices to (0, 0, 0). Defaults to True.
Returns:

A number_of_samples by 3 array of angles and alexander polynomials.

Return type:

ndarray

arclength()[source]

Calls pyknotid.spacecurves.spacecurve.SpaceCurve.arclength(), automatically not including the closure.

closing_distance()[source]

Returns the distance between the first and last points.

closure_alexander_polynomial(theta=0, phi=0)[source]

Returns the Alexander polynomial of the knot, when projected in the z plane after rotating the given theta and phi to the North pole.

Parameters:
  • theta (float) – The sphere angle theta
  • phi (float) – The sphere angle phi
generalised_alexander()[source]

Returns the generalised Alexander polynomial for the default projection of the open knot

multiroots_fractions(number_of_samples=10, **kwargs)[source]

Returns each of the knot types from self.alexander_polynomials_multiroots, with the fraction of that type.

plot_alexander_map(number_of_samples=10, scatter_points=False, mode='imshow', interpolation=100, **kwargs)[source]

Creates (and returns) a projective diagram showing each different Alexander polynomial in a different colour according to a closure on a far away point in this direction.

Parameters:
  • number_of_samples (int) – The number of points on the sphere to close at.
  • scatter_points (bool) – If True, plots a dot at each point on the map projection where a closure was made.
  • mode (str) – ‘imshow’ to plot the pixels of an image, otherwise plots filled contours. Defaults to ‘imshow’.
  • interpolation (int) – The (short) side length of the interpolation grid on which the map projection is made. Defaults to 100.
plot_alexander_shell(number_of_samples=100, mode='mesh', radius=None, **kwargs)[source]

Plots the curve in 3d via self.plot(), along with a translucent sphere coloured by the type of knot obtained by closing on each point.

Parameters are all passed to OpenKnot.alexander_polynomials(), except opacity and kwargs which are given to mayavi.mesh, and sphere_radius_factor which gives the radius of the enclosing sphere in terms of the maximum Cartesian distance of any point in the line from the origin.

plot_projections(number_of_samples)[source]

Plots the projection of the knot at each of the given number of samples squared, rotated such that the sample direction is vertical.

The output (and return) is a matplotlib plot with number_of_samples x number_of_samples axes.

plot_self_linking_map(number_of_samples=10, scatter_points=False, mode='imshow', **kwargs)[source]

Creates (and returns) a projective diagram showing each different self linking number in a different colour according to a projection in this direction.

plot_self_linking_shell(number_of_samples=100, **kwargs)[source]

Plots the curve in 3d via self.plot(), along with a translucent sphere coloured by the self linking number obtained by projecting from this point.

Parameters are all passed to OpenKnot.virtual_checks(), except opacity and kwargs which are given to mayavi.mesh, and sphere_radius_factor which gives the radius of the enclosing sphere in terms of the maximum Cartesian distance of any point in the line from the origin.

plot_virtual_map(number_of_samples=10, scatter_points=False, mode='imshow', **kwargs)[source]

Creates (and returns) a projective diagram showing each different virtual Boolean in a different colour according to a projection in this direction.

plot_virtual_shell(number_of_samples=10, zero_centroid=False, sphere_radius_factor=2.0, opacity=0.3, **kwargs)[source]

Plots the curve in 3d via self.plot(), along with a translucent sphere coloured according to whether or not the projection from this point corresponds to a virtual knot or not.

Parameters are all passed to OpenKnot.virtual_checks(), except opacity and kwargs which are given to mayavi.mesh, and sphere_radius_factor which gives the radius of the enclosing sphere in terms of the maximum Cartesian distance of any point in the line from the origin.

points

The points of the spacecurve, as an Nx3 numpy array.

projection_invariant(**kwargs)[source]

First checks if the projection of an open curve is virtual or classical. If virtual, a virtual knot invariant is calculated. Otherwise a classical invariant is calculated.

raw_crossings(mode='use_max_jump', virtual_closure=False, recalculate=False, try_cython=False)[source]

Calls pyknotid.spacecurves.spacecurve.SpaceCurve.raw_crossings(), but without including the closing line between the last and first points (i.e. setting include_closure=False).

self_linking(theta=0, phi=0)[source]

Takes an open curve, finds its Gauss code (for the default projection) and calculates its self linking number, J(K). See Kauffman 2004 for more information.

Returns:The self linking number of the open curve
Return type:self_link_counter : int
self_linking_fractions(number_of_samples=10, **kwargs)[source]

Returns each of the self linking numbers from self.virtual.self_link.projections, with the fraction of each type.

self_linkings(number_of_samples=10, zero_centroid=False, **kwargs)[source]

Returns a list of self linking numbers for the curve with a given number of projections taken from directions approximately evenly distributed.

Parameters:
  • number_of_samples (int) – The number of points on the sphere to project from. Defaults to 10.
  • zero_centroid (bool) – Whether to first move the average position of vertices to (0, 0, 0). Defaults to False.
Returns:

A number_of_samples by 3 array of angles and self linking number

Return type:

ndarray

smooth(repeats=1, window_len=10, window='hanning')[source]

Calls pyknotid.spacecurves.spacecurve.SpaceCurve.smooth(), with the periodic argument automatically set to False.

vassiliev_degree_2_average(samples=10, recalculate=False, **kwargs)[source]

Returns the average Vassliev degree 2 invariant calculated by averaging its combinatorial value over many different projection directions.

Parameters:
  • samples (int) – The number of directions to average over. Defaults to 10.
  • recalculate (bool) – Whether to recalculate the writhe even if a cached result is available. Defaults to False.
  • **kwargs – These are passed directly to raw_crossings().
virtual_check()[source]

Takes an open curve and checks (for the default projection) if its Gauss code corresponds to a virtual knot or not. Returns a Boolean of this information.

Warning

This only checks the distance by which entries in the Gauss code are separated, it is not guaranteed to detect virtual knots.

Returns:virtual – True if the Gauss code corresponds to a virtual knot. False otherwise.
Return type:bool
virtual_checks(number_of_samples=10, zero_centroid=False)[source]

Returns a list of virtual Booleans for the curve with a given number if projections taken from directions approximately evenly distributed. A value of True corresponds to the projection giving a virtual knot, with False returned otherwise.

Parameters:
  • number_of_samples (int) – The number of points on the sphere to project from. Defaults to 10.
  • zero_centroid (bool) – Whether to first move the average position of vertices to (0, 0, 0). Defaults to False.
Returns:

A number_of_samples by 3 array of angles and virtual Booleans (True if virtual, False otherwise)

Return type:

ndarray

virtual_fractions(number_of_samples=10, **kwargs)[source]

Returns each of the virtual booleans from self.virtual.check.projections, with the fraction of each type.

pyknotid.spacecurves.openknot.gall_peters(theta, phi)[source]

Converts spherical coordinates to the Gall-Peters projection of the sphere, an area-preserving projection in the shape of a Rectangle.

Parameters:
  • theta (float) – The latitude, in radians.
  • phi (float) – The longitude, in radians.
pyknotid.spacecurves.openknot.knot_db_to_string(database_object)[source]

Takes output from from_invariants() and returns knot type as decimal. For example: <Knot 3_1> becomes 3.1 and <Knot K13n1496> becomes 13.1496

pyknotid.spacecurves.openknot.mollweide(phi, lambda_)[source]

Converts spherical coordinates to the Mollweide projection of the sphere, an area-preserving projection in the shape of an ellipse.

Parameters:
  • phi (float) – The latitude, in radians.
  • lambda (float) – The longitude, in radians.