pygridgen & pygridtools

A successful collaboration with the private sector, academia, and the open source community

Transcript

1. pygridgen & pygridtools A successful collaboration with the private sector, academia, and the open source community Paul Hobson Geosyntec Consultants Scipy 2017 Austin, TX
2. Overview • Acknowledgements • Personal and topical background • Basics of curvilinear-orthogonal grids • The Foundation (suite of C libraries) • The Abstraction (basic python tools) • The Community (conda-forge) • The Enhancements (interactivity) • Questions
3. Acknowledgements • Pavel Sakov at the Australian BoM – wrote the code • Robert Hetland at Texas A&M – wrote the python bindings • Filipe Fernandes, Phil Elson, Rich Signell (USGS) • Lucas Nguyen & Rica Enriquez at Geosyntec (our employer) • Scipy conference organizers and reviewers
5. Grid Vocabulary • Boundary or Domain– Polygon enclosing the grid • Turning Point – Point on the boundary that is defined to either close the polygon (β=+1) or open a new side channel (β = -1) • Node or Vertex – a single point in the generated grid • Shape – number of nodes in each direction • Cell – The quadrilateral formed by 4 adjacent nodes • Centroid – The geometric center of a cell
6. Utilizing negative turning points for side channels
7. Grid generation from 100k ft up • “Numerical Conformal Mapping” • Find the boundary’s regional corners • Delaunay triangulation • Schwarz-Christoffel transforms • Computing the σ’s, β’s, & ζ’s • Partial differential equations (Broyden’s Method) • ??? • Grids!
8. The Foundation • Four foundational C-libraries by Pavel Sakov (Australian Bureau of Meteorology) • Hard to compile from source (*nix systems, very difficult) • Still can’t compile on Windows • Breadth of examples serving as user docs • Central control file points to boundary file, defined other parameters (e.g., number of nodes) • Called executable from terminal
9. Param file input xy.0 output grid.0 nx 120 ny 90 nnodes 10 precision 1.0e-12 newton 1 rectangle rect.0 sigmas sigmas.0
12. Grid Generation – The Abstraction • pygridgen – python bindings for the C libraries using ctypes • Vendored C sources • Authored by Dr. Robert Hetland at Texas A&M with contributions from Rich Signell at USGS. • Allows users to use familiar numpy arrays • Interactive grid boundary creation via MPL widgets • Exposed and API to utilize Focus…
13. What is focus? Uses the error function to refine/rarify the unit cartesian grid onto which the curvilinear grid is eventually mapped.
14. Adding detail with focus
15. Using pygridgen x = numpy.array([0.50, 2.00, 2.00, ...]) y = numpy.array([0.50, 0.50, 1.75, ...]) beta = numpy.array([1, 1, -1, 0, ...]) focus = pygridgen.Focus() focus.add_focus(0.50, 'y', factor=5, extent=0.25) focus.add_focus(0.50, 'x', factor=3, extent=0.50) grid = pygridgen.Gridgen(x, y, beta, shape=(10, 20), focus=focus)
16. The community • The C libraries were hard to compile and link properly • Equally difficult to link the Python bindings • Growing demand from oceanographers, hydrodynamic modelers • conda-forge saved us all
17. conda-forge • Community-driven conda packages • Package recipes hosted as “feedstocks” in individual Github repositories • Commits (changes to the repo) trigger continuous integration that build packages for all platforms when possible) conda install <pkg_name> --channel=conda-forge
18. Benefits of conda-forge • Dead simple install • conda install pygridgen --channel=conda-forge • Grabs built C libraries • Windows :( • Recipes serve as canonical build instructions for the C libraries
19. The Enhancements • Needed to make a grid for ~60 miles of braided river for a sediment and pollutant fate and transport model • Fixed a few bugs • Added and clarified some docstrings • Cached slow-to-compute grid domain properties (the σ’s, β’s, & ζ’s) • Sphinx docs at http://phobson.github.io/pygridgen/
21. pygridtools • Higher-level API wrapper around pygridgen • Visualization tools via matplotlib • Spatial data I/O via fiona • “Pandas-like” operations for manipulating the grids • Merging • Spliting • Refining • Masking • Etc • Docs at http://geosyntec.github.io/pygridtools/
23. Merge/refine/etc example
24. Merge/refine/etc example fig = ( grid2.merge(grid3, how='vert', where='-', shift=2) .merge(grid1, how='horiz', where='-', shift=11) .refine(10, axis=1, n_points=4) .refine(13, axis=0, n_points=2) .transform(lambda x: x*5 + 2) .update_cell_mask() .plotCells(ax=ax) )
25. Masking
26. Masking
27. Sullivan Creek, WA Grid boundary River bank Island
28. Sullivan Creek, WA nx, ny = 100, 20 focus = pgg.Focus() # tighten the grid in the channels # around the big island & coarsen upstream focus.add_focus(5./ny, 'y', 4., extent=8./ny) focus.add_focus(14.5/ny, 'y', 4., extent=4./ny) focus.add_focus(52./nx, 'x', 4., extent=20./nx) focus.add_focus(98./nx, 'x', 0.25, extent=4./nx) # generate the main grid grid = pgt.makeGrid(domain=gridbounds, nx=nx, ny=nx, focus=focus)
29. Sullivan Creek, WA # outside the banks grid.mask_cells_with_polygon(river, inside=False) # inside the multiple islands for island in island_arrays: grid.mask_cells_with_polygon(island, inside=True)
35. The Future • Clean up and pull in existing plotting code • Interactive visualization • Jupyter widgets to update grid parameters • Better docs • Cython rewrite? • Scipy rewrite?
36. Where to find me • Twitter @pmhobson • Github @phobson and @geosyntec • Email phobson@geosyntec.com