u/No_Fun_3602

I am trying to refine clinker of cement and i dont have some CIF files i cant access from the free database

thanks

yes i tried it but coulding get it working

thank you

# I tried to go about it like this but the graph i get is different from what in the paper attached in the original post

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp

# Dome parameters
R = 385 
lambda_val = 0.1 
z_range = np.linspace(0, R, 100)  

# Radial coordinate r(z)
def r(z):
    return np.sqrt(R**2 - z**2)

def r_prime(z):
    return -z / np.sqrt(R**2 - z**2)

def r_double_prime(z):
    return -R**2 / (R**2 - z**2)**(3/2)
def dalpha_dz(z, alpha):
    r_val = r(z)
    r_prime_val = r_prime(z)
    r_double_prime_val = r_double_prime(z)
    
    
    term1 = np.sin(alpha) * np.tan(alpha) / r_val
    term2 = r_double_prime_val / (1 + r_prime_val**2) * np.cos(alpha)
    term3 = r_prime_val * np.tan(alpha) / r_val
    
    d_alpha_dz = lambda_val * (term1 - term2 - term3)
    return d_alpha_dz

# Initial condition
alpha_0 = 0
sol = solve_ivp(dalpha_dz, [z_range[0], z_range[-1]], [alpha_0], t_eval=z_range)

# Plot the results
plt.plot(sol.t, sol.y[0], label=r'$\alpha(z)$')
plt.xlabel('Axial Coordinate (z)')
plt.ylabel('Winding Angle (α)')
plt.title('Winding Angle α as a Function of z')
plt.legend()
plt.grid(True)
plt.show()

I cant seem to plot the figure shown in the research paper

i am trying to make a python app that i can give my input variables and it plot the trajectory based on those inputs

i tried this but it only plot the cylinder without the correct path
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

c = 1.5  # minimum radius at the ends in cm
R = 4.0  # equatorial (middle section) radius in cm
L = 8.0  # total length of the cylinder in cm
transition_length = L * 0.2  # length of the transition sections at the ends
winding_angle = np.radians(85)  # winding angle in degrees

z = np.linspace(-L/2, L/2, 100)

radii = np.piecewise(
    z,
    [z < -L/2 + transition_length,
     (z >= -L/2 + transition_length) & (z <= L/2 - transition_length),
     z > L/2 - transition_length],
    [lambda z: c + (R - c) * np.sqrt(1 - ((z + L/2 - transition_length) / transition_length) ** 2),
     R,
     lambda z: c + (R - c) * np.sqrt(1 - ((L/2 - transition_length - z) / transition_length) ** 2)]
)

theta = np.linspace(0, 2 * np.pi, 100)
theta_grid, z_grid = np.meshgrid(theta, z)
x_grid = np.outer(radii, np.cos(theta))
y_grid = np.outer(radii, np.sin(theta))

n_turns = 5  # Number of turns for the winding
z_winding = np.linspace(-L/2, L/2, 500)
theta_winding = z_winding * np.tan(winding_angle) / R  # Geodesic winding equation
x_winding = R * np.cos(theta_winding)
y_winding = R * np.sin(theta_winding)

fig = plt.figure(figsize=(8, 6))
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(x_grid, y_grid, z_grid, color='cyan', alpha=0.7, edgecolor='k')
ax.plot(x_winding, y_winding, z_winding, color='red', linewidth=2, label="Geodesic Helical Winding")

ax.set_xlabel('X (cm)')
ax.set_ylabel('Y (cm)')
ax.set_zlabel('Z (cm)')
ax.set_title('3D Model of Cylinder with Geodesic Helical Winding')
ax.legend()
plt.show()

the main part i need help with is coming up with the logic to show the trajectory on the cylinder

Yes that's essentially the main obstacle I came across

ok got it

will i be able to get a similar outcome as in the video?