Plotting streamlines of the velocity field within an ice sheet

[verhaegenyoni]

Hi all,

I have a flow field with horizontal (realxh) and vertical (realzh) coordinates and ice flow velocities (vertical (wi) and horizontal (ui)), and I want to plot 2D velocity vectors within the ice sheet. I can do it with quiver:

fig, ax = plt.subplots()
ax.quiver(realxh[0:nz,0:nx-1], realzh[0:nz,0:nx-1], ui[0:nz,0:nx-1], wi[0:nz,0:nx-1]*(domainsize/(np.max(sur)-np.min(bed))),color='b')
ax.set_title('Ice sheet 2D flow field')
plt.ylabel('Elevation (m)')
plt.xlabel('Distance (m)')
plt.show()

which results in the plot in attachment. I think it looks okay, but the vectors at the right are quite big (so big they go out the ice sheet outlines) and at the bottom they are quite small. I therefore wanted to try streamlines instead of quiver (which should plot continuous lines of the ice flow instead of vectors), to make it look more aesthetic/realistic:

fig, ax = plt.subplots()
ax.streamplot(realxh[0:nz,0:nx-1], realzh[0:nz,0:nx-1], ui[0:nz,0:nx-1], wi[0:nz,0:nx-1]*(domainsize/(np.max(sur)-np.min(bed))),color='b')
ax.set_title('Ice sheet 2D flow field')
plt.ylabel('Elevation (m)')
plt.xlabel('Distance (m)')
plt.show()

But it results in the following error:

Error:
Traceback (most recent call last):
  File "/Users/yv/PycharmProjects/icesheet/main.py", line 382, in <module>
    ax.streamplot(realxh[0:nz,0:nx-1], realzh[0:nz,0:nx-1], ui[0:nz,0:nx-1], wi[0:nz,0:nx-1]*(domainsize/(np.max(sur)-np.min(bed))),color='b')
  File "/Users/yv/PycharmProjects/icesheet/.venv/lib/python3.12/site-packages/matplotlib/streamplot.py", line 342, in __init__
    raise ValueError("The columns of 'y' must be equal")
ValueError: The columns of 'y' must be equal

Does anyone know how to resolve this? Apparently the shape of my coordinate files are not valid in this case. I have looked online for hours but did not find the answer. It works fine for quiver, but not for streamlines. I put the .mat files of the matrices in attachment (zip) so you can test it out (can be done by from scipy.io import loadmat). You can put nx= 31, nz = 15, and domainsize/(np.max(sur)-np.min(bed)) = 175. The latter scale factor is done to scale the y-axis with the size of the x-axis (the vertical ice thickness and its discretization scheme is lots thinner than the ice sheet horizontal extent).

Larz60+ write Mar-21-2024, 10:31 AM:
Please post all code, output and errors (it it's entirety) between their respective tags. Refer to BBCode help topic on how to post. Use the "Preview Post" button to make sure the code is presented as you expect before hitting the "Post Reply/Thread" button.
Tags have been added this time, please use BBCode tags on future posts.

[verhaegenyoni]

Nobody?:(