Search code examples
pythonnumpymappingcoordinate-transformation

3D-cartesian mesh to spherical section mesh conversion


I am modelling convection in spherical section using Finite element method. I am trying to deform 3D cartesian mesh into spherical section with 3480 <= radius <= 6371 km and the angular range in theta and phi is 80 degrees. Please check this image, left side figure is the 3D cartesian mesh which I have and right side is the required deformed geometry of the mesh. I have used some transformation function to deform the mesh but unable to reach final goal. I am looking for some suggestions on algorithm or transformation function to deform the mesh into a spherical section.


Solution

  • Not sure this is really the same, but by eye it looks close enough:

    import numpy as np
    
    phi, theta, R = np.ogrid[-40:40:21j, -40:40:21j, 3480:6371:21j]
    
    y, x = np.tan(phi*np.pi/180), np.tan(theta*np.pi/180)
    z = R / np.sqrt(x*x + y*y + 1)
    y, x = y*z, x*z
    
    from mpl_toolkits.mplot3d import Axes3D
    import pylab
    fig = pylab.figure(1)
    ax = fig.add_subplot(111, projection='3d')
    ax.plot_wireframe(x[..., 20], y[..., 20], z[..., 20])
    ax.plot_wireframe(0.01*x[..., 20], 0.01*y[..., 20], 0.01*z[..., 20]) # hack to scale z axis
    ax.plot_wireframe(x[:, 20, :], y[:, 20, :], z[:, 20, :])
    ax.plot_wireframe(x[0, ...], y[0, ...], z[0, ...])
    fig.show()
    

    enter image description here