Wrapping a given coordinates (x,y,z) in a twisted structure

Question

I am trying to wrap (x,y,z) coordinates on a Mobius strip (a topological structure obtained by twisting once and connecting the ends of a strip).

The structure (for which I can create xyz coordinates) is as follow (from RG) using following code.

import numpy as np
import matplotlib.pyplot as plt

bLength=1.6
numPoints=10
radius = bLength*numPoints / (2 * np.pi)
theta = np.linspace(0,2*np.pi,numPoints,endpoint=False)
dtheta=theta[1]-theta[0]
x0,y0=(radius * np.cos(theta)), (radius * np.sin(theta))
x1,y1=(radius * np.cos(theta+dtheta/2)) , (radius * 
np.sin(theta+dtheta/2))
#plt.plot(x0,y0)
#plt.show()
#plt.plot(x1,y1)
#plt.show()
cons0=np.ones(x0.shape)*0
cons1=np.ones(x1.shape)*2
np.savetxt('cooRing1.csv',np.c_[x0,y0,cons0],delimiter=' ')
np.savetxt('cooRing2.csv',np.c_[x1,y1,cons1],delimiter=' ')
cat cooRing1.csv cooRing2.csv > coordinates.csv

[![enter image description here][3]][3]

I want to map above xyz co-ordinates on a Mobius strip.

Following example has been given in an website for a different strip (also the code uses a module which is not publicly available)

import Twister as Twister    
import math

def displacement(x, width, wrapping_angle):
    """
       Function for converting a nanosheet coordinate into a partly wrapped nanotube
       @param x : Coordinates of nanosheet atom
       @param width : Width of the nano-sheet
       @param wrapping_angle : maximum wrapping angle of the nanotube in radians
    """
    # calculate the average radius of the incomplete wrapped tube
    radius = width/wrapping_angle
    # find the angle of the current atom
    angle = (x[2]-width/2.)/radius
    # calculate the radius of the current atom
    atom_radius = radius+x[1]

    # return atom position of the wrapped atom
    return numpy.array([x[0], atom_radius*math.cos(angle),atom_radius*math.sin(angle)])

def configuration(n, m, repetition):
    """
       Function for generating a moebius molecule
       @param n : Chiral vector index
       @param m : Chiral vector index
       @param repetition : Repetition along z
    """

    # build n,m ribbon
    ribbon = NanoRibbon(n,m)
    ribbon = ribbon.repeat(1,1,repetition)

    # get properties of the ribbon
    lattice = ribbon.bravaisLattice()
    elements = ribbon.elements()
    cartesian_coordinates=ribbon.cartesianCoordinates().inUnitsOf(Angstrom)

    # calculate the length of the 1-d structure
    z_length = numpy.linalg.norm(lattice.primitiveVectors()[2].inUnitsOf(Angstrom)) 

    # calculate twist parameters
    rotation_angle_per_z =  math.pi /z_length
    rotation_axis = numpy.array([0,0,1])
    rotation_axis_center = numpy.sum(cartesian_coordinates,axis=0)/len(cartesian_coordinates)

    # define a function of one variable, f(c),  for displacing the atoms
    f = lambda c : Twister.displacement(c, rotation_angle_per_z, rotation_axis, 
                                rotation_axis_center, 0.,z_length)
    # apply the function to find new displaced coordinates
    cartesian_coordinates = numpy.apply_along_axis(f, 1, cartesian_coordinates)
    cartesian_center = numpy.sum(cartesian_coordinates,axis=0)/len(cartesian_coordinates)
    cartesian_coordinates = cartesian_coordinates - cartesian_center


    # define a function of one variable, f(c),  for displacing the atoms
    f = lambda c : displacement(c, z_length,2.0*math.pi) 
    # apply the function to find new displaced coordinates
    cartesian_coordinates = numpy.apply_along_axis(f, 1, cartesian_coordinates)

    return MoleculeConfiguration(
            elements=elements,
            cartesian_coordinates=cartesian_coordinates * Angstrom
            )

# Instantiate the builder object and choose our title
builder = Builder()
builder.title('Moebius ribbon')

# Set the configuration generator
builder.setConfigurationGenerator(configuration)

# Tube properties group
builder.newGroup('Tube parameters')
builder.integer( 'n', 6, 'n', min=1, max=1000)
builder.integer( 'm', 0, 'm', min=0, max=1000)
builder.integer( 'repetition', 40, 'C-direction', min=1, max=1000)

Is there any similar module in Python, so that I can build the structure I want and also create the xyz coordinates? Or How I can proceed further to complete the code?

Answer 1

There is no possible way to transform the coordinates for one of your 2xN rings into a Möbius strip. That's because an alternating strip with a half twist in it must have an odd number of atoms, while your current rings always have an even number. You'd need to cut one atom out or add one, in order to make the twist work. Maybe you could make a transform that could work by stacking the first and last atoms on top of each other, but I think it would be very ugly (both mathematically and perhaps also in a plot).

While that probably makes it impractical to what you're asking with a 3D-transform, you can instead create coordinates from scratch with the desired twist. Just produce a single array of points for the single edge of the strip (having only one edge is one of the things that makes a Möbius strip strange). The edge makes two loops around the circle, one on each side of the strip, and the twist happens at half the speed of the main rotation.

Here's my stab at it:

bLength=1.6
n = 10
numPoints = 2*n + 1
radius = bLength*numPoints / (4*np.pi)   # 4*pi because the edge goes around the circle twice
theta = np.linspace(0, 4*np.pi, numPoints, endpoint=False)   # here too
x = (radius + np.cos(theta/2)) * np.cos(theta)
y = (radius + np.cos(theta/2)) * np.sin(theta)
z = np.sin(theta/2)

plt.plot(x, y) # I don't know how to graph in 3d, so only x and y here
plt.show()

This makes a strip of width of 2, just like your current code does. If that's not the correct width, you should multiply the cos terms added to the radius in the calculation of x and y and the sin value used for z by half the desired width.

Note that this while this code may produce what you want, the coordinates are probably nonsense when it comes to physically describing an actual ring of atoms. Some of the distances between the coordinates will be very different than others (e.g on the inside of the strip versus the outside when it's nearly flat), and a real Möbius strip made of atoms (if such a thing can be made at all) would probably fold in some way rather than twisting at a uniform rate (a Möbius strip made of paper doesn't twist this way either). Finding exactly how real atoms would be arranged would be a much more difficult problem (and a problem for a physicist, not a programmer).