chlnddev / oceanmesh Goto Github PK
View Code? Open in Web Editor NEWAutomatic coastal ocean mesh generation in Python and C++. **under development**
License: GNU General Public License v3.0
Automatic coastal ocean mesh generation in Python and C++. **under development**
License: GNU General Public License v3.0
om
will re-project this data into Mercator to perform the meshing.for applications with phase resolving wave models.
Hi,
I found an off-by-one error while trying to read Etopo1 data with DEM
class. Well, first I need to convert Etopo data from 0-360 to the 180-180 referential, then I get this error upon reading:
bbox, min_edge_length = (-75.000, -70.001, 40.0001, 41.9000), 1e3 # from the example
dem = DEM("etopo1.nc", bbox)
Trace-back:
ValueError Traceback (most recent call last)
~/Github/ww3-brazil-wc/create_mesh.py in <module>
22
23 fname = "/home/stringari/GIS/etopo1fix.nc"
---> 24 dem = DEM(fname, bbox)
25
26 edge_length = distance_sizing_function(shore, max_edge_length=5e3)
~/Github/oceanmesh/oceanmesh/geodata.py in __init__(self, dem, bbox, verbose)
501 dy = numpy.abs(lats[1] - lats[0])
502 dx = numpy.abs(lons[1] - lons[0])
--> 503 super().__init__(
504 bbox=bbox,
505 dx=dx,
~/Github/oceanmesh/oceanmesh/grid.py in __init__(self, bbox, dx, dy, hmin, values, fill)
66 self.nx = int(floor((self.bbox[1] - self.bbox[0]) / self.dx))
67 self.ny = int(floor((self.bbox[3] - self.bbox[2]) / self.dy))
---> 68 self.values = values
69 self.eval = None
70 self.fill = fill
~/Github/oceanmesh/oceanmesh/grid.py in values(self, data)
118 pass
119 elif data.shape != (self.nx, self.ny):
--> 120 raise ValueError(
121 f"Shape of values {data.shape} does not match grid size {self.nx, self.ny}"
122 )
ValueError: Shape of values (300, 114) does not match grid size (299, 113)
If I modify lines 66 and 67 of grid.py
to
self.nx = int((self.bbox[1] - self.bbox[0]) // self.dx)
self.ny = int((self.bbox[3] - self.bbox[2]) // self.dy)
I can read both my file and the tests pass. Should I open a PR?
Caio
Hi there. Looking to download and use oceanmesh for generating a triangular mesh for my DG-FEM code for modelling flow in and around submarine canyons, given by some fluid-depth function varying in space. I am however having issues with getting oceanmesh to install and work on my windows machine, and could benefit from help. Unfortunately, the slack channel link you posted seems to no longer be active.
The @setter
for grid.values
is misleading during the calculation of the sizing functions that require a dem.
The values grid.dx
and grid.dy
change here for the bathy gradient and here for the wavelength
Here an example for the world model where both the wavelength and bathy gradient have an offset:
The quick fix (without touching the Grid
class) is to force in both functions the dx
and dy
in the sizing functions.
grid.dx = dem.dx
grid.dy = dem.dy
Line 22 in c59f039
def _cell_to_cell(t):
"""Cell to cell connectivity table.
Cell `i` is connected to cells `ctoc[ix[i]:ix[i+1]]`
By connected, I mean shares a mutual edge.
Parameters
----------
t: array-like
Mesh connectivity table.
Returns
-------
ctoc: array-like
Cell numbers connected to cells.
ix: array-like
indices into `ctoc`
"""
nt = len(t)
t = np.sort(t, axis=1)
# NB: we use order="F" to reshape because the `np.tile` command below
e = t[:, [[0, 1], [0, 2], [1, 2]]].reshape((nt * 3, 2), order="F")
trinum = np.tile(np.arange(nt), 3)
j = _arg_sortrows(e)
e = e[j, :]
trinum = trinum[j]
k = np.argwhere(~np.diff(e, axis=0).any(axis=1))
ctoc = np.concatenate((trinum[k], trinum[k + 1]), axis=1)
ctoc = np.append(ctoc, np.fliplr(ctoc), axis=0)
ctoc = ctoc[np.argsort(ctoc[:, 0]), :]
idx = np.argwhere(np.diff(ctoc[:, 0])) + 1
idx = np.insert(idx, 0, 0)
idx = np.append(idx, len(ctoc))
return ctoc, idx
Hello,
Thank you for your work on the development and improvement of OceanMesh.
I am currently studying OceanMesh and am curious to know if the functions of the Matlab-based version are fully implemented in the Python version.
a msh class to provide basic meta data and storage capabilities
a port of the topological cleaning routine in the original version.
We need Grid addition and min operations overloaded.
Mesh size grading
In the README.md, it is suggested that to ensure that the vertices of each element are arranged in a counter-clockwise fashion, one simply uses the command points, cells = fix_mesh(points, cells)
. I had the error of two many values to unpack. Upon investigation, I could see that fix_mesh
returns three outputs. To fix this in my own code, I simply changed it to
points, cells, jx = fix_mesh(points, cells)
Sorry, just thought others would find that helpful. (I have changed that in mine, and will submit a pull request once I am satisfied with other edits I have made.
This issue is related to what was mentioned in the PR for global meshes about the stereographic distortion.
Close to the poles, the triangles get distorted and don't have the right equilateralness.
As mentioned,
it has to do with the forces being too powerful
To solve this issue, I have implemented a parametric law in the mesh_generator()
:
def _parametric(lat):
ones = np.ones(lat.shape)
res = ((90 - lat) * 2 + 18) / 180 * np.pi
return np.minimum(res, ones)
and illustrated it in the following python notebook
I will propose a PR with this parametric law. But I am still not convinced this is the best way to go..
The meshes look almost perfect closer to the poles now:
Keen to get another opinion on this.
def dfs(graph):
seen = set()
lifo = []
start = random_node(graph)
lifo.append(start)
while lifo:
n = lifo.pop()
seen.add(n)
for nei in neighbours(n, graph):
if nei in seen:
continue
else:
lifo.append(n)
nflag = numpy.zeros(nnodes)
nflag[list(seen)] = 1
and growing outward in the radial direction.
Hi!
I am trying to figure out the different uses of DEM() since I have a netcdf that is not working. So I tried both the string and array options with an example that I know it works.
I have a netcdf, that if I load and plot with xarray:
plt.figure()
plt.contourf(ds.lon.data, ds.lat.data, ds.z.data)
plt.colorbar()
If I then use DEM() it all works correctly:
dem_etopo_built_nc = f"/home/filippog/Mesh/dem_etopo_built_{name}.nc"
dem_etopo_fine = om.DEM(dem_etopo_built_nc)
dem_etopo_fine.plot()
I wanted to get the same result giving an array as input.
1st ISSUE: you need to explicitly give bbox to DEM() since it fails to create it inside the function if not given (bbox=None) and then it crashes in lines 677-679 of geodata.py
bbox = (8, 10, 44, 45) # Liguria
dem_etopo_fine = om.DEM(ds.z.values, bbox=(8, 10, 44, 45))
dem_etopo_fine.plot()
After playing with the array I figured out that you need to give it as:
dem_etopo_fine = om.DEM(ds.z.values[::-1, :], bbox=(8, 10, 44, 45))
But cannot figure out why I need to invert the lon?
Any thoughts or ideas?
Thanks a lot!!
feature size
topographic length scale.
mesh size gradation.
cfl limiting
Hi guy, thanks for your work, please help me with this issue.
I try to run the examples you give in "README.md" or "test cases". Every time it runs into the code line with "distance_sizing_function" it has an error.
shoreline = om.Shoreline(fname, extent.bbox, min_edge_length)
edge_length = om.distance_sizing_function(shoreline, #rate=0.15)
Error:
raise ValueError("The points in dimension %d must be strictly "
ValueError: The points in dimension 1 must be strictly ascending
bbox
is a 2-tuple.@WPringle, do any points get generated in the mesh gen code other than the ones initially created?
oceanmesh/oceanmesh/mesh_generator.py
Line 146 in b763b1a
@jcharris @krober10nd @sgriffithjones I see you had the same linting issue in the CI lint with black :
would reformat /home/runner/work/oceanmesh/oceanmesh/oceanmesh/idw.py
would reformat /home/runner/work/oceanmesh/oceanmesh/oceanmesh/geodata.py
Oh no! 💥 💔 💥
2 files would be reformatted, 29 files would be left unchanged.
I see 2 options here:
I guess 1 is the best one in the long run. I think we should also update black locally for python 3.11
The grid.values
are given in SI units of metres in wavelength_sizing_function
of oceanmesh/edgefx.py
. These need to be translated to the grid coordinate system. Should be an easy fix?
In order to get a better representation of rivers and watersheds in the coastal ocean model, I'm trying to generate a mesh for extracted river polygons:
These polygons are extracted from DEM and are not necessarily connected. Apart from that I'd like to use all the vertices and edges that the extraction algorithm has outputed to guide the meshing. I can use pfix
for nodes, but is there anything similar for edges?
I know that probably I cannot rely on Shoreline
and mesh sizing functions provided for "normal" coastal meshing, so I tried to hack it by defining mine (without worrying about performance for now). This is how I define the signed distance function
def get_dist(x):
pts = gpd.points_from_xy(x[:, 0], x[:, 0])
sign = pts.within(shp_geom) * -2 + 1
dist = sign * pts.distance(shp_geom.boundary)
return dist
sdf = om.Domain(extent.bbox, get_dist)
And this is how I define the edge length.
edge_length = om.Grid(
bbox=extent.bbox,
dx=hmin,
hmin=hmin,
extrapolate=True,
values=hmax,
crs=gdf_river_polys.crs,
)
edge_length.build_interpolant()
points, cells = om.generate_mesh(
sdf,
edge_length,
pfix=df_river_coord[['lon', 'lat']].values,
max_iter=20
)
I take the hmin
to be the smallest of the initial edges and hmin
to be the size of the largest edge. I'm still playing with what to specify at as edge_length
function, but my ideal solution is to actually not have any size specified and the mesher to just triangulate the nodes specified as initial points pfix
.
When installing oceanmesh personalized python library with pip install -U -e .
in the oceanmeshmater directorym I get a long big thread of errors due to pybind11
.
I already installed the requirements sudo apt install libcgal-dev python3-pybind11
.
The log with the errors from the installation of oceanmesh with pip is bellow. Did you get the same kind of errors? How could I solve this?
Installing collected packages: oceanmesh
Running setup.py develop for oceanmesh
ERROR: Command errored out with exit status 1:
command: /home/nuno/resources/miniconda3/bin/python -c 'import sys, setuptools, tokenize; sys.argv[0] = '"'"'/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py'"'"'; __file__='"'"'/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py'"'"';f=getattr(tokenize, '"'"'open'"'"', open)(__file__);code=f.read().replace('"'"'\r\n'"'"', '"'"'\n'"'"');f.close();exec(compile(code, __file__, '"'"'exec'"'"'))' develop --no-deps
cwd: /home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/
Complete output (107 lines):
running develop
running egg_info
writing oceanmesh.egg-info/PKG-INFO
writing dependency_links to oceanmesh.egg-info/dependency_links.txt
writing requirements to oceanmesh.egg-info/requires.txt
writing top-level names to oceanmesh.egg-info/top_level.txt
reading manifest file 'oceanmesh.egg-info/SOURCES.txt'
reading manifest template 'MANIFEST.in'
writing manifest file 'oceanmesh.egg-info/SOURCES.txt'
running build_ext
-- OS detected: Linux
-- CXX Compiler detected: GNU
-- CMake additional search path for libraries:
-- Using header-only CGAL
-- Targetting Unix Makefiles
-- Using /usr/bin/c++ compiler.
-- Boost version: 1.65.1
-- Boost include dirs: /usr/include
-- Boost libraries:
-- Using gcc version 4 or later. Adding -frounding-math
-- CGAL version: 5.00.0.100
-- Build type: Release
-- USING CXXFLAGS = '-DVERSION_INFO=\"0.0.1\" -O3 -DNDEBUG'
-- USING EXEFLAGS = ' '
-- Requested component: MPFR
-- Requested component: GMP
-- Found pybind11 v2.0.1: /usr/include;/home/nuno/resources/miniconda3/include/python3.8
-- Configuring done
-- Generating done
-- Build files have been written to: /home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/build/temp.linux-x86_64-3.8
[ 25%] Building CXX object CMakeFiles/fast_geometry.dir/oceanmesh/cpp/fast_geometry.cpp.o
[ 50%] Building CXX object CMakeFiles/delaunay_class.dir/oceanmesh/cpp/delaunay_class.cpp.o
/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/oceanmesh/cpp/fast_geometry.cpp: In function ‘pybind11::array unique_edges(pybind11::array_t<int, 17>)’:
/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/oceanmesh/cpp/fast_geometry.cpp:76:23: error: no matching function for call to ‘pybind11::buffer_info::buffer_info(int*, long unsigned int, std::__cxx11::string, int, std::vector<long int>&, std::vector<long int>&)’
));
^
In file included from /usr/include/pybind11/pytypes.h:12:0,
from /usr/include/pybind11/cast.h:13,
from /usr/include/pybind11/attr.h:13,
from /usr/include/pybind11/pybind11.h:36,
from /usr/include/pybind11/complex.h:12,
from /home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/oceanmesh/cpp/fast_geometry.cpp:2:
/usr/include/pybind11/common.h:287:5: note: candidate: pybind11::buffer_info::buffer_info(pybind11::buffer_info&&)
buffer_info(buffer_info &&other) {
^~~~~~~~~~~
/usr/include/pybind11/common.h:287:5: note: candidate expects 1 argument, 6 provided
/usr/include/pybind11/common.h:274:14: note: candidate: pybind11::buffer_info::buffer_info(Py_buffer*, bool)
explicit buffer_info(Py_buffer *view, bool ownview = true)
^~~~~~~~~~~
/usr/include/pybind11/common.h:274:14: note: candidate expects 2 arguments, 6 provided
/usr/include/pybind11/common.h:270:5: note: candidate: pybind11::buffer_info::buffer_info(void*, pybind11::size_t, const string&, pybind11::size_t)
buffer_info(void *ptr, size_t itemsize, const std::string &format, size_t size)
^~~~~~~~~~~
/usr/include/pybind11/common.h:270:5: note: candidate expects 4 arguments, 6 provided
/usr/include/pybind11/common.h:262:5: note: candidate: pybind11::buffer_info::buffer_info(void*, pybind11::size_t, const string&, pybind11::size_t, const std::vector<long unsigned int>&, const std::vector<long unsigned int>&)
buffer_info(void *ptr, size_t itemsize, const std::string &format, size_t ndim,
^~~~~~~~~~~
/usr/include/pybind11/common.h:262:5: note: no known conversion for argument 5 from ‘std::vector<long int>’ to ‘const std::vector<long unsigned int>&’
/usr/include/pybind11/common.h:260:5: note: candidate: pybind11::buffer_info::buffer_info()
buffer_info() { }
^~~~~~~~~~~
/usr/include/pybind11/common.h:260:5: note: candidate expects 0 arguments, 6 provided
/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/oceanmesh/cpp/fast_geometry.cpp: At global scope:
/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/oceanmesh/cpp/fast_geometry.cpp:79:16: error: expected constructor, destructor, or type conversion before ‘(’ token
PYBIND11_MODULE(fast_geometry, m) {
^
CMakeFiles/fast_geometry.dir/build.make:62: recipe for target 'CMakeFiles/fast_geometry.dir/oceanmesh/cpp/fast_geometry.cpp.o' failed
make[2]: *** [CMakeFiles/fast_geometry.dir/oceanmesh/cpp/fast_geometry.cpp.o] Error 1
CMakeFiles/Makefile2:67: recipe for target 'CMakeFiles/fast_geometry.dir/all' failed
make[1]: *** [CMakeFiles/fast_geometry.dir/all] Error 2
make[1]: *** Waiting for unfinished jobs....
/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/oceanmesh/cpp/delaunay_class.cpp:66:16: error: expected constructor, destructor, or type conversion before ‘(’ token
PYBIND11_MODULE(delaunay_class, m) {
^
CMakeFiles/delaunay_class.dir/build.make:62: recipe for target 'CMakeFiles/delaunay_class.dir/oceanmesh/cpp/delaunay_class.cpp.o' failed
make[2]: *** [CMakeFiles/delaunay_class.dir/oceanmesh/cpp/delaunay_class.cpp.o] Error 1
CMakeFiles/Makefile2:104: recipe for target 'CMakeFiles/delaunay_class.dir/all' failed
make[1]: *** [CMakeFiles/delaunay_class.dir/all] Error 2
Makefile:83: recipe for target 'all' failed
make: *** [all] Error 2
Traceback (most recent call last):
File "<string>", line 1, in <module>
File "/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py", line 81, in <module>
setup(
File "/home/nuno/resources/miniconda3/lib/python3.8/site-packages/setuptools/__init__.py", line 153, in setup
return distutils.core.setup(**attrs)
File "/home/nuno/resources/miniconda3/lib/python3.8/distutils/core.py", line 148, in setup
dist.run_commands()
File "/home/nuno/resources/miniconda3/lib/python3.8/distutils/dist.py", line 966, in run_commands
self.run_command(cmd)
File "/home/nuno/resources/miniconda3/lib/python3.8/distutils/dist.py", line 985, in run_command
cmd_obj.run()
File "/home/nuno/resources/miniconda3/lib/python3.8/site-packages/setuptools/command/develop.py", line 34, in run
self.install_for_development()
File "/home/nuno/resources/miniconda3/lib/python3.8/site-packages/setuptools/command/develop.py", line 136, in install_for_development
self.run_command('build_ext')
File "/home/nuno/resources/miniconda3/lib/python3.8/distutils/cmd.py", line 313, in run_command
self.distribution.run_command(command)
File "/home/nuno/resources/miniconda3/lib/python3.8/distutils/dist.py", line 985, in run_command
cmd_obj.run()
File "/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py", line 44, in run
self.build_extension(ext)
File "/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py", line 76, in build_extension
subprocess.check_call(
File "/home/nuno/resources/miniconda3/lib/python3.8/subprocess.py", line 364, in check_call
raise CalledProcessError(retcode, cmd)
subprocess.CalledProcessError: Command '['cmake', '--build', '.', '--config', 'Release', '--', '-j2']' returned non-zero exit status 2.
----------------------------------------
ERROR: Command errored out with exit status 1: /home/nuno/resources/miniconda3/bin/python -c 'import sys, setuptools, tokenize; sys.argv[0] = '"'"'/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py'"'"'; __file__='"'"'/home/nuno/universidade/tese/unst_grid_generation/oceanmesh-master/setup.py'"'"';f=getattr(tokenize, '"'"'open'"'"', open)(__file__);code=f.read().replace('"'"'\r\n'"'"', '"'"'\n'"'"');f.close();exec(compile(code, __file__, '"'"'exec'"'"'))' develop --no-deps Check the logs for full command output.
This is one of the reasons that the CI actions don't pass:
Some of the examples in the README use EPSG = 4326 as an integer for the coordinate system.
However in edgefx.py, only strings are tested i.e. :
if crs == "EPSG:4326":
like in bathymetric_gradient_sizing_function()
and wavelength_sizing_function()
if crs == "EPSG:4326" or crs == 4326:
seems to be a good fix
numpy just upgraded to version 2.0.0 and now make the tests in CI fail:
A module that was compiled using NumPy 1.x cannot be run in
NumPy 2.0.0 as it may crash. To support both 1.x and 2.x
versions of NumPy, modules must be compiled with NumPy 2.0.
Some module may need to rebuild instead e.g. with 'pybind11>=2.12'.
If you are a user of the module, the easiest solution will be to
downgrade to 'numpy<2' or try to upgrade the affected module.
We expect that some modules will need time to support NumPy 2.
Not only CI fails, but also new installs.
For new comers to be able to use this package, they would just have to
pip install numpy==1.26
at the end of the install.
I think we should downgrade to numpy<2
- we wouldn't recommend using it after a few releases anyway.
At the same occasion, I think it might be good to set up a proper pyproject.toml
for new people to be able to install oceanmesh
properly.
And also update the README
accordlingly
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.