import numpy as np
import matplotlib.pyplot as plt
import mikeio
Let's consider a simple mesh consisting of 2 triangular elements.
fn = "data/two_elements.mesh"
with open(fn, "r") as f:
print(f.read())
msh = mikeio.open(fn)
msh
msh.plot(show_mesh=True);
msh.node_coordinates
msh.element_table
msh.element_coordinates
msh.get_element_area()
Let's plot the node and element coordinates:
xn, yn = msh.node_coordinates[:,0], msh.node_coordinates[:,1]
xe, ye = msh.element_coordinates[:,0], msh.element_coordinates[:,1]
ax = msh.plot(show_mesh=True)
ax.plot(xn, yn, 'ro', markersize=10)
ax.plot(xe, ye, 'bx', markersize=10)
It can sometimes be convenient to have mesh boundary as a polyline (or multiple in case of more complex meshes).
bxy = msh.boundary_polylines.exteriors[0].xy
plt.plot(bxy[:,0], bxy[:,1])
plt.axis("equal");
MIKE IO has a method for determining if a point (or a list of points) is inside the domain:
pt_1 = [2.0, 1.2]
msh.contains(pt_1)[0]
# or multiple points at the same time
pt_2 = [4.0, 1.2]
pts = np.array([pt_1, pt_2])
msh.contains(pts)
plt.plot(bxy[:,0], bxy[:,1], label='boundary')
plt.plot(xe[0], ye[0], 'b*', markersize=10, label="center, elem 0")
plt.plot(xe[1], ye[1], 'c*', markersize=10, label="center, elem 1")
plt.plot(*pt_1, 'go', markersize=10, label="pt_1")
plt.plot(*pt_2, 'rs', markersize=10, label="pt_2")
plt.axis("equal")
plt.legend(loc="upper right");
MIKE IO has a method for obtaining the index of the element containing a point:
g = msh.geometry
g.find_index(coords=pt_1)[0]
MIKE IO also has a method for obtaining a list of the n closest element centers:
g.find_nearest_elements(pt_1)
g.find_nearest_elements(pt_1, return_distances=True)
g.find_nearest_elements(pt_1, n_nearest=2)
# for multiple points
g.find_nearest_elements(pts, return_distances=True)
dfs = mikeio.open("data/FakeLake.dfsu")
g = dfs.geometry
g
g.plot();
# insert code here
g.max_nodes_per_element
msh = mikeio.open("data/FakeLake.dfsu").geometry
msh.plot();
msh.node_coordinates[:,2] = np.clip(msh.node_coordinates[:,2], -15, 0) # clip depth to interval [-15,0]
msh.plot(title="No change??")
del msh.element_coordinates # remove cached element coords calculated based on original node coords)
msh.plot(title="Updated")
msh.to_mesh('Fake_lake_clip15.mesh') # save to a new file
msh = mikeio.open("data/southern_north_sea.mesh")
msh
The default is to plot the elements and color them according to the bathymetry.
msh.plot();
msh.plot.outline();
msh.plot.mesh();
Maybe we would like to higlight the bathymetric variations in some range, in this case in the -40, -20m range.
msh.plot(vmin=-40, vmax=-20);
There are other options as well, such as explicit specification of which contour lines to show or choosing a specific colormap (matplotlib colormaps)
msh.plot.contour(show_mesh=True,
levels=[-50,-30,-20,-10,-5], cmap="tab10",
figsize=(12,12), title="Coarse North Sea model");
See the MIKE IO Mesh Example notebook for more Mesh operations (including shapely operations).