Simple wave sediment modelling¶
Note
Using a wave propagation model based on linear wave theory, we investigate the cumulative effect of a given wave forcing condition on erosional and depositional trend across a regional coastal system.
Loading required modules for this exercise. The main module here is wavesed a piece of code that was developed for this course and that can be found on GitHub
!cd wavesed/; f2py -c -m ocean ocean.f90
running build
running config_cc
unifing config_cc, config, build_clib, build_ext, build commands --compiler options
running config_fc
unifing config_fc, config, build_clib, build_ext, build commands --fcompiler options
running build_src
build_src
building extension "ocean" sources
f2py options: []
f2py:> /tmp/tmp09p_7588/src.linux-x86_64-3.6/oceanmodule.c
creating /tmp/tmp09p_7588/src.linux-x86_64-3.6
Reading fortran codes...
Reading file 'ocean.f90' (format:free)
Post-processing...
Block: ocean
Block: ocean
Block: airymodel
Block: transport
Block: diffusion
Post-processing (stage 2)...
Block: ocean
Block: unknown_interface
Block: ocean
Block: airymodel
Block: transport
Block: diffusion
Building modules...
Building module "ocean"...
Constructing F90 module support for "ocean"...
Variables: pi onpi pi2 pion2 onpi2 grav
Constructing wrapper function "ocean.airymodel"...
c,l,travel,waveh = airymodel(dx,dd,h0,depth,src,inland,shadow,[numrow,numcol])
Constructing wrapper function "ocean.transport"...
dz,dist = transport(iter,depth,hent,transx,transy,[numrow,numcol])
Constructing wrapper function "ocean.diffusion"...
depo = diffusion(oelev,dz,coeff,maxth,tstep,nstep,[numrow,numcol])
Wrote C/API module "ocean" to file "/tmp/tmp09p_7588/src.linux-x86_64-3.6/oceanmodule.c"
Fortran 90 wrappers are saved to "/tmp/tmp09p_7588/src.linux-x86_64-3.6/ocean-f2pywrappers2.f90"
adding '/tmp/tmp09p_7588/src.linux-x86_64-3.6/fortranobject.c' to sources.
adding '/tmp/tmp09p_7588/src.linux-x86_64-3.6' to include_dirs.
copying /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/f2py/src/fortranobject.c -> /tmp/tmp09p_7588/src.linux-x86_64-3.6
copying /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/f2py/src/fortranobject.h -> /tmp/tmp09p_7588/src.linux-x86_64-3.6
adding '/tmp/tmp09p_7588/src.linux-x86_64-3.6/ocean-f2pywrappers2.f90' to sources.
build_src: building npy-pkg config files
running build_ext
customize UnixCCompiler
C compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-cc -DNDEBUG -fwrapv -O2 -Wall -Wstrict-prototypes -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -DNDEBUG -D_FORTIFY_SOURCE=2 -O2 -isystem /usr/share/miniconda/envs/coast/include -fPIC
creating /tmp/tmpcnr870fu/tmp
creating /tmp/tmpcnr870fu/tmp/tmpcnr870fu
compile options: '-MMD -MF /tmp/tmpcnr870fu/file.c.d -c'
x86_64-conda-linux-gnu-cc: /tmp/tmpcnr870fu/file.c
customize UnixCCompiler using build_ext
get_default_fcompiler: matching types: '['gnu95', 'intel', 'lahey', 'pg', 'nv', 'absoft', 'nag', 'vast', 'compaq', 'intele', 'intelem', 'gnu', 'g95', 'pathf95', 'nagfor']'
customize Gnu95FCompiler
Found executable /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran
Found executable /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-ld
Found executable /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-ar
Found executable /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-ranlib
customize Gnu95FCompiler
customize Gnu95FCompiler using build_ext
building 'ocean' extension
compiling C sources
C compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-cc -DNDEBUG -fwrapv -O2 -Wall -Wstrict-prototypes -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -DNDEBUG -D_FORTIFY_SOURCE=2 -O2 -isystem /usr/share/miniconda/envs/coast/include -fPIC
creating /tmp/tmp09p_7588/tmp
creating /tmp/tmp09p_7588/tmp/tmp09p_7588
creating /tmp/tmp09p_7588/tmp/tmp09p_7588/src.linux-x86_64-3.6
compile options: '-I/tmp/tmp09p_7588/src.linux-x86_64-3.6 -I/usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include -I/usr/share/miniconda/envs/coast/include/python3.6m -c'
x86_64-conda-linux-gnu-cc: /tmp/tmp09p_7588/src.linux-x86_64-3.6/oceanmodule.c
x86_64-conda-linux-gnu-cc: /tmp/tmp09p_7588/src.linux-x86_64-3.6/fortranobject.c
In file included from /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/ndarraytypes.h:1822,
from /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/ndarrayobject.h:12,
from /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/arrayobject.h:4,
from /tmp/tmp09p_7588/src.linux-x86_64-3.6/fortranobject.h:13,
from /tmp/tmp09p_7588/src.linux-x86_64-3.6/fortranobject.c:2:
/usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/npy_1_7_deprecated_api.h:17:2: warning: #warning "Using deprecated NumPy API, disable it with " "#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION" [-Wcpp]
17 | #warning "Using deprecated NumPy API, disable it with " \
| ^~~~~~~
In file included from /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/ndarraytypes.h:1822,
from /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/ndarrayobject.h:12,
from /usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/arrayobject.h:4,
from /tmp/tmp09p_7588/src.linux-x86_64-3.6/fortranobject.h:13,
from /tmp/tmp09p_7588/src.linux-x86_64-3.6/oceanmodule.c:16:
/usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include/numpy/npy_1_7_deprecated_api.h:17:2: warning: #warning "Using deprecated NumPy API, disable it with " "#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION" [-Wcpp]
17 | #warning "Using deprecated NumPy API, disable it with " \
| ^~~~~~~
/tmp/tmp09p_7588/src.linux-x86_64-3.6/oceanmodule.c:109:12: warning: 'f2py_size' defined but not used [-Wunused-function]
109 | static int f2py_size(PyArrayObject* var, ...)
| ^~~~~~~~~
compiling Fortran 90 module sources
Fortran f77 compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -ffixed-form -fno-second-underscore -fPIC -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -O3 -funroll-loops
Fortran f90 compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -fno-second-underscore -fPIC -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -O3 -funroll-loops
Fortran fix compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -ffixed-form -fno-second-underscore -Wall -g -fno-second-underscore -fPIC -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -O3 -funroll-loops
compile options: '-I/tmp/tmp09p_7588/src.linux-x86_64-3.6 -I/usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include -I/usr/share/miniconda/envs/coast/include/python3.6m -c'
extra options: '-J/tmp/tmp09p_7588/ -I/tmp/tmp09p_7588/'
x86_64-conda-linux-gnu-gfortran:f90: ocean.f90
compiling Fortran sources
Fortran f77 compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -ffixed-form -fno-second-underscore -fPIC -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -O3 -funroll-loops
Fortran f90 compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -fno-second-underscore -fPIC -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -O3 -funroll-loops
Fortran fix compiler: /usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -ffixed-form -fno-second-underscore -Wall -g -fno-second-underscore -fPIC -fopenmp -march=nocona -mtune=haswell -ftree-vectorize -fPIC -fstack-protector-strong -fno-plt -O2 -ffunction-sections -pipe -isystem /usr/share/miniconda/envs/coast/include -O3 -funroll-loops
compile options: '-I/tmp/tmp09p_7588/src.linux-x86_64-3.6 -I/usr/share/miniconda/envs/coast/lib/python3.6/site-packages/numpy/core/include -I/usr/share/miniconda/envs/coast/include/python3.6m -c'
extra options: '-J/tmp/tmp09p_7588/ -I/tmp/tmp09p_7588/'
x86_64-conda-linux-gnu-gfortran:f90: /tmp/tmp09p_7588/src.linux-x86_64-3.6/ocean-f2pywrappers2.f90
/usr/share/miniconda/envs/coast/bin/x86_64-conda-linux-gnu-gfortran -Wall -g -Wall -g -shared -Wl,-O2 -Wl,--sort-common -Wl,--as-needed -Wl,-z,relro -Wl,-z,now -Wl,--disable-new-dtags -Wl,--gc-sections -Wl,-rpath,/usr/share/miniconda/envs/coast/lib -Wl,-rpath-link,/usr/share/miniconda/envs/coast/lib -L/usr/share/miniconda/envs/coast/lib /tmp/tmp09p_7588/tmp/tmp09p_7588/src.linux-x86_64-3.6/oceanmodule.o /tmp/tmp09p_7588/tmp/tmp09p_7588/src.linux-x86_64-3.6/fortranobject.o /tmp/tmp09p_7588/ocean.o /tmp/tmp09p_7588/tmp/tmp09p_7588/src.linux-x86_64-3.6/ocean-f2pywrappers2.o -L/usr/share/miniconda/envs/coast/bin/../x86_64-conda-linux-gnu/sysroot/lib/../lib -L/usr/share/miniconda/envs/coast/bin/../x86_64-conda-linux-gnu/sysroot/lib/../lib -lgfortran -o ./ocean.cpython-36m-x86_64-linux-gnu.so
Removing build directory /tmp/tmp09p_7588
import time
import numpy as np
import wavesed.wave as ocean
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
Model set-up¶
The model requires several parameters to be defined before being run. The most important one being the initial bathymetry. For this example, we will use the region around Fraser Island at a resolution of 1000 m. The digital elevation model (DEM) comes from the Project 3DGBR.
The file is called FraserDEM.csv and is located in the data folder. It is made of 3 columns containing the X, Y, and Z coordinates.
Model domain / grid parameters¶
It is possible to change the initial resolution by using the rfac
parameter, as an example:
rfac = 1
corresponds to a model at 1000 m resolution (similar to the initial DEM)rfac = 2
transforms the dataset to 2000 m resolution
# Bathymetric filename
bfile = '../pracenv/dataset/FraserDEM.csv'
# Resolution factor
rfac = 1
Definition of wave parameters¶
Now we need to define the forcing parameters for our simulation:
H0
the deep water wave height in [m],dir
the deep water wave source direction at the boundary in degrees (angle is defined counterclock wise from horizontal axis)wbase
the maximum depth for wave influence on seabed [m]slvl
the sea-level position in [m]
# Wave height
H0 = 2
# Wave direction at boundary
dir = 0
# Maximum depth for wave influence (m)
wbase = 10
# Sea level position (m)
slvl = 0.
# Output figure name
outname = 'wdir0'
Definition of sediment parameters¶
To estimate wave-induced sediment transport, this simple model only accounts for a single grain size that needs to be defined (d50
).
There are also additional parameters that can be set to perform longshore and diffusive transport but we will use the default values for this example…
# Mean grain size diameter in m
d50 = 0.0001
Model initialisation function¶
Now that all the parameters are defined in our python environment, we can initialise the simulation by calling the wavesed module:
#help(ocean.wave.__init__)
wavesed = ocean.wave(filename=bfile, wavebase=wbase, resfac=rfac, dia=d50)
As you can see from the command above there are still some parameters that have not been passed to the module:
slvl
dir
H0
These parameters are sent separately so that multiple climatic conditions could be simulate more easily.
First we define the sea-level position slvl
.
t0 = time.clock()
wavesed.findland(slvl)
print('Definition of shoreline position took (s):',time.clock()-t0)
Definition of shoreline position took (s): 0.0030769999999999964
Initialising wave boundary direction dir
is then done with the wavesource
function:
wdir = wavesed.wavesource(dir)
Waves evolution¶
The waves are then transformed using the linear wave theory from deep (set-up conditions) to shallow water assuming shore-parallel depth contours.
The orientation of wave fronts is determined by wave celerity and refraction due to depth variations.
Travel time in the domain is calculated from Huygen’s principle (using an order \(\sqrt{5}\) approximation).
Assuming no refraction or loss of energy due to bottom friction, wave power \(P\) is conserved from deep to shallow water.
#help(wavesed.cmptwaves)
t0 = time.clock()
wavesed.cmptwaves(src=wdir, h0=H0, sigma=1.)
print('Wave evolution computation took (s): ',time.clock()-t0)
Wave evolution computation took (s): 0.6260160000000001
Plotting simulated wave characteristics¶
A series of graph can be visualised to show the following characteristics:
bathy: coarse bathymetric map
wlength: coarse wave lenght map
travel: coarse wave travel time map
wcelerity: coarse wave celerity map
wheight: coarse wave height map
wpower: coarse wave power map
ubot: coarse wave induced bottom velocity map
shear: coarse wave induced shear stress map
#help(wavesed.plotData)
size = (6,6)
# i1 = 0
# i2 = -1
# j1 = 0
# j2 = -1
# Zooming to a specific region
i1 = 80
i2 = 250
j1 = 0
j2 = 170
wavesed.plotData(data='bathy', figsize=size, vmin=-100, vmax=100,
fontsize=10, imin=i1, imax=i2, jmin=j1, jmax=j2)
wavesed.plotData(data='travel', figsize=size, tstep=100, vmin=0, vmax=5.e4,
fontsize=10, imin=i1, imax=i2, jmin=j1, jmax=j2)
wavesed.plotData(data='wcelerity', figsize=size, vmin=0, vmax=15,
fontsize=10, stream=3, imin=i1, imax=i2, jmin=j1, jmax=j2)
wavesed.plotData(data='ubot', figsize=size, vmin=0, vmax=2,
fontsize=10, imin=i1, imax=i2, jmin=j1, jmax=j2)
wavesed.plotData(data='shear', figsize=size, vmin=-10, vmax=10,
fontsize=10, imin=i1, imax=i2, jmin=j1, jmax=j2)
Sediment entrainment, transport & deposition¶
Sediment entrainment relates to wave induced shear stress.
The transport is computed according to both wave direction and longshore transport. Deposition is dependent of shear stress and diffusion.
#help(wavesed.cmptsed)
t0 = time.clock()
wavesed.cmptsed(sigma=1.)
print('Sediment erosion/deposition computation took (s): ',time.clock()-t0)
Sediment erosion/deposition computation took (s): 0.5779190000000014
Plotting morphological changes¶
Using the same plotting function (plotData
), we can now plot the following metrics:
ent: coarse sediment entrainment map
dz: coarse erosion/deposition map
#help(wavesed.plotData)
size = (6,6)
# i1 = 0
# i2 = -1
# j1 = 0
# j2 = -1
# Zooming to a specific region
i1 = 80
i2 = 250
j1 = 0
j2 = 170
wavesed.plotData(data='bathy', figsize=size, vmin=-100, vmax=100,
fontsize=10, stream=0, imin=i1, imax=i2, jmin=j1, jmax=j2)
wavesed.plotData(data='ent', figsize=size, vmin=0., vmax=5.,
fontsize=10, stream=0, imin=i1, imax=i2, jmin=j1, jmax=j2)
wavesed.plotData(data='dz', figsize=size, vmin=-1., vmax=1.,
fontsize=10, stream=0, imin=i1, imax=i2, jmin=j1, jmax=j2) #, save=outname)
Saving wave/sedimentation data¶
#waveparams.outputCSV(filename='erodep.csv', seddata=1)
Note
Change the initial wave direction conditions to 90 and then to 270. Evaluate the impact on the morphological patterns.
Change the resolution factor to see how it impacts on erosion/deposition resolution.