Display the area in common between two surfaces

Common area between two surfaces? In this blog post, we will explore how to achieve this using Python, NumPy, and Matplotlib. By following the provided code and steps, you’ll be able to visualize the shared region between two surfaces efficiently.

To start, let’s examine the code snippet provided. The code generates a paraboloid and extracts the portion of the area above a lower limit or below an upper limit. However, you want to determine the common area between these two plots. Here’s how you can modify the code to achieve that:

pythonCopy code















































import numpy as np
import matplotlib.pyplot as plt

if __name__ == "__main__":
    r: float = 1.5
    r_min: float = 0.5 * r
    r_max: float = 0.75 * r
    n_points: int = 100

    r: np.ndarray = np.linspace(0, r, n_points)
    r_pellet: np.ndarray = np.linspace(r_min, r_max, n_points)
    theta: np.ndarray = np.linspace(0, 2 * np.pi, n_points)
    R, THETA = np.meshgrid(r, theta)

    # Deal with figures
    fig = plt.figure()
    ax = fig.add_subplot(projection="3d")

    # Create projectile surface
    X, Y = R * np.cos(THETA), R * np.sin(THETA)
    Z = (X ** 2 + Y ** 2)
    ax.plot_surface(X, Y, Z, cmap="cool", linewidth=0, antialiased=True)

    # Create pellet surface
    # Above lower limit
    X_SUP = np.where(r_min <= X ** 2 + Y ** 2, X, np.NAN)
    Y_SUP = np.where(r_min <= X ** 2 + Y ** 2, Y, np.NAN)
    Z_SUP = X_SUP ** 2 + Y_SUP ** 2

    # Under upper limit
    X_INF = np.where(X ** 2 + Y ** 2 <= r_max, X, np.NAN)
    Y_INF = np.where(X ** 2 + Y ** 2 <= r_max, Y, np.NAN)
    Z_INF = X_INF ** 2 + Y_INF ** 2

    ax.plot_surface(X, Y, Z_INF, cmap="hot", linewidth=0, antialiased=True)

    # Common area between the two surfaces
    X_INF_SUP = np.where((X ** 2 + Y ** 2 >= r_min) & (X ** 2 + Y ** 2 <= r_max), X, np.NAN)
    Y_INF_SUP = np.where((X ** 2 + Y ** 2 >= r_min) & (X ** 2 + Y ** 2 <= r_max), Y, np.NAN)
    Z_INF_SUP = X_INF_SUP ** 2 + Y_INF_SUP ** 2

    ax.plot_surface(X_INF_SUP, Y_INF_SUP, Z_INF_SUP, cmap="spring", linewidth=0, antialiased=True)

    ax.set_xlabel("X")
    ax.set_ylabel("Y")
    ax.set_zlabel("Z")

    plt.show()

In the modified code, we added the section for the common area between the two surfaces. Using NumPy’s logical AND operation (&), we ensure that both conditions, (X ** 2 + Y ** 2 >= r_min) and `(X ** 2 + Y ** 2 <= r_max