Skip to content

- Add option where IG source term uses mean wave energy Eprev_ig in s…#283

Draft
Leynse wants to merge 39 commits into
mainfrom
281-add-version-of-ig-source-term-where-dsxxdx-does-not-depend-on-the-bins
Draft

- Add option where IG source term uses mean wave energy Eprev_ig in s…#283
Leynse wants to merge 39 commits into
mainfrom
281-add-version-of-ig-source-term-where-dsxxdx-does-not-depend-on-the-bins

Conversation

@Leynse

@Leynse Leynse commented Feb 12, 2026

Copy link
Copy Markdown
Collaborator

…ource term calculation instead of bin dependent eeprev_ig, thanks to suggestion of Yasmine - results look promising

…ource term calculation instead of bin dependent eeprev_ig, thanks to suggestion of Yasmine - results look promising
@Leynse Leynse linked an issue Feb 12, 2026 that may be closed by this pull request
…ave energy, instead of the bin one

- Works, but not sure whether it is better than snapwave_ig_opt =3
…king point when IG waves should be released

- Determined based on delta Dw
…Zijlema (2022)

- Does at first glance not seem to influence much/enough
…h) where waves start breaking

- not sure this is the best place to do this...
- Add local wave height water depth ration 'gam' to his output
- For now with Qb determined based on upwind point
- For now with taking maxval over itheta for both qb and gam output
…king is started the velocity of cg_ig is changed to sqrt(g*h)

- In igopt 8 for now set only point 'k' to sqrt(g*h), not yet k1 and k2, to be discussed
- Does not work well in shadow zones currently, question is mostly how to determine Sxx_cons
… on integrated energies, and at the redistributing them over bins like for Dw

- Result pretty much the same as before with igopt 9
…cig as alrready mutliplied by ee_itheta,k) and in 'B' term

- ig opt 12 can be the working version in R
- Results don't look good yet, not sure whats wrong, it is like Snapwave trunk commit of Dano of 03-04-2025
…alue to multiply with snapwave_gamma factor for estimating the incident wave breaking point initiation, and setting the srcig to 0

- description: factor times gamma that is used to determine the maximum incident wave breaking point in the surf zone using local incident wave height over water depth ratio, among others used to set the IG source term to 0 shallower than this point
- Now include the /sqrt(2) so that value becomes gam (Hrms based) > 0.47 * gamma
- used for ig_opt >= 11
…lead to issues, therefore remove for ig_opt 12&12
- Only keep the (already precleaned) ig_opt 11&12 version of subroutine determine_infragravity_source_sink_term
- Add description qb_local, gam_local
…e_sink_term'

- Rename the main option 11 to ig_opt = 1 (=the default)
- For now keep also option 11 for ease, but exactly the same as 1
- Option 2 is name the same as ig_opt 1, but without using the conservative shoaling (as before)
- ig_opt = 3 is the same as option 1, but with including the cg_ig adjustment
- turn off ig_opt 12 - currently does not work, need to further check moving from R- to B-term
Leynse added 10 commits March 10, 2026 18:26
… is a more gradual decreasing of srcig to 0 .

- 'transition_factor' is controled around 'gamma_fac_br * gamma' as center
- And width controlled by transition_factor_width > 0.005 seems to work well
- With this, found best gamma_fac_br to be 0.45 now - adjusted default value
…with limiting beta_local to 0.07, as suggested by Maarten for steep slopes
…transitioned so that for beta 0.05 the factor is 0.45 (default snapwave_gamma_fac_br), but for steeper coasts (beta > 0.07) it is fully 1.0 (no reduction of srcig) because then driven by breakpoint forcing
- Remove some of the options
- let ig_opt =1 &2 be the simple one without smooth srcig to 0 functions
… for addition to alphaig

- Values set to initial estimate by Yasmine
- Tested to work for 1:3 test
- Todo: make factors user definable
…coded values of before. Result is therefore the same
Leynse added a commit that referenced this pull request Apr 2, 2026
Leynse added a commit that referenced this pull request Apr 29, 2026
@Leynse Leynse modified the milestone: 2026.02 release May 6, 2026
Leynse added a commit that referenced this pull request May 27, 2026
* some refactoring of snapwave solver

* deleted snapwave_solver.f90.ori

* - Make exponent for multiplying for extra dissipation user defineable with 'snapwave_baldock_exponent', with 0 meaning unused as default.
- can be set to 1 or 2 for more dissipation for steep coasts
- baldock_opt not used anymore I see

* - Bump version for clarity
- nc output fix

* - As done in PR of jre/snapwave_sfincs, use the SFINCS input routines directly in SnapWave
- Made clearer by putting the read functions in new filed 'sfincs_read', and only call those in sfincs_snapwave.f90

* - Also in sfincs_spiderweb.f90 replace local read_real_input and read_int_input by use sfincs_read

* - add posibility to inspect snapwave grid as in PR jre/snapwave_sfincs

* - restart already fixed to true in sfincs_snapwave

* - Remove unused subroutine 'neuboundaries' that is replaced by neuboundaries_light

* - Add back some description

* - Bugfix, actually update Hk so it is used in IG bottom friction determination

* - Preprocess snapwave_ncoutput

* - Keep snapwave_crit= 0.001 as default as in main

* - make (:,k) to (:, k) consistently, everywhere I hope

* - Try adding back missing stuff

* making it work for wind

* removed second call to celerities

* - Bump version for clarity when testing

* - Make Herbers IG boundary warning dimensionless, using hsinc / depth
- Check further for best ratio later, start at 0.5

* - Add user definable snapwave_relax_factor_DoverE and snapwave_relax_factor_DoverA keywords to sfincs_snapwave
- By default set to 0.25 for both for now

* - For now do simple Hmx(k)    = gamma * depth(k) and IG as in trunk, TBD

* - As discussed, remove fdrspr for now. Could be added back later on (with actual working function)

* - Add back comment

* - Remove subroutine bj78

* - Added calculation and write_log message of error_ig and %ok_ig with thanks to Mr Claude
- For now criterion 'ok' of incident waves for determining iteration ended is unaltered

* - Add omp parallellisation

* - Already do the bugfix of PR #283 by switching to srcig based on E and redistribute by /E * ee

* - Add back Neumann message

* - Bugfix in calculating kwav_ig!

* - Suggestion of Marlies to keep the original formulation for Hmx_ig

* Reinit Tp for inner/Neumann; ignore docs/_build

Reinitialize the per-node Tp only for internal and Neumann-connected cells during the wind iteration so boundary cells keep their prescribed Tp (as set by update_boundaries). Implemented a loop that sets Tp(k)=Tpini for inner(k) and Tp(neumannconnected(k)) where present. Also added docs/_build to .gitignore and cleaned up minor whitespace.

* - Add gamma_fac_br as in PR #281

* 331 adapt wave force limiter for lab tests (#333)

* - issue for Buckley solved if no limiter at all
- Later improve limiter

* - Temp save, need to still make it work here/in sfincs_snapwave

* - Option to determine fwmaxfac once every 'update_wave_field' call in sfincs_snapwave

* 292 renewed branch vegetation effects on short waves and wave forces (#296)

* merged veggie shortwave drag branch143 into current main

* added missing vfproj update

* - First step to making mean flow correct on the uv points

* - clean repo

* - Clean up this PR, so that only the pure snapwave veggie changes remain, the mean flow effect is still kept in copy of this branch: 315-veggie-effect-on-mean-flow

* - Cleanup and switch to logicals

* - Start of sfincs_vegetation

* - Only for quadtree

* - Move veggie reference out of sfincs_quadtree

* - Add new generic subroutine 'read_netcdf_quadtree_get_dimension'
- cleanup

* - Neat way to import only needed veggie variables from sfincs_data into snapwave_domain

* - Cleanup code

* - Manual merge back swvegnonlin, momeqveg, update F(k)

* - For now limit vegetation_vertical_segments to max 4

* - move check from sfincs_snapwave to wavemaker (not veggie related), if snapwave yes but wavemaker no, the message before still appeared

* Add toml-f and TOML src_structure parsing

Add the toml-f third-party TOML parser (sources and licenses) and wire it into the build (Makefile.am and Visual Studio .vfproj). Implement TOML-based parsing for src_structure inputs in sfincs_src_structures.f90: new types, constants, read_toml_src_structures routine, validation helpers and a dispatcher in initialize_src_structures() that probes for TOML and falls back to the legacy reader. Also include minor non-functional whitespace/formatting cleanups in sfincs_discharges.f90.

* - Add TOML-based vegetation type lookup table with CF flag NetCDF input

Vegetation parameters can now be specified via a two-file approach:
- The NetCDF file (vegetationfile) stores only an integer vegetation_type
  per grid cell, following CF flag_values / flag_meanings conventions so
  the file is self-documenting in standard tools (ncview, Panoply, xarray).
- A companion TOML file (vegfile_toml) contains the parameter lookup table
  keyed by type name, supporting any number of vertical layers per type
  (mixed layer counts across types are handled by zero-padding).

The old single-NetCDF code path is fully preserved for backward
compatibility; the new path is activated only when vegfile_toml is set
in sfincs.inp.

Changes:
- sfincs_vegetation.f90: add read_vegetation_toml(); rewrite
  initialize_vegetation() with old/new path split; lift 4-layer cap to 64
- sfincs_ncinput.F90: add read_netcdf_quadtree_integer() and
  read_netcdf_flag_meanings() for the new input format
- sfincs_data.f90: add veggiefile_toml, vegetation_type_index; rename
  vegetation_cd -> vegetation_stems_cd and
  vegetation_stems_width -> vegetation_stems_diameter for consistency
- sfincs_input.f90: add vegfile_toml keyword
- snapwave_domain.f90: update to renamed variable names

* - vegetationtype_toml as input name - to discuss

* - Visual cleanup

* - Bump version

* - missing stop_sfincs reference

---------

Co-authored-by: Leynse <Tim.Leijnse@deltares.nl>
Co-authored-by: Tim Leijnse <timleynse94@gmail.com>
Co-authored-by: Maarten van Ormondt <maarten.vanormondt@deltares.nl>

---------

Co-authored-by: Leynse <Tim.Leijnse@deltares.nl>
Co-authored-by: Tim Leijnse <timleynse94@gmail.com>
Co-authored-by: Marlies <marlies.vanderlugt@deltares.nl>
Leynse added a commit that referenced this pull request Jun 15, 2026
* Snapwave domain updates (#300)

* some refactoring of snapwave solver

* deleted snapwave_solver.f90.ori

* - Make exponent for multiplying for extra dissipation user defineable with 'snapwave_baldock_exponent', with 0 meaning unused as default.
- can be set to 1 or 2 for more dissipation for steep coasts
- baldock_opt not used anymore I see

* - Bump version for clarity
- nc output fix

* - As done in PR of jre/snapwave_sfincs, use the SFINCS input routines directly in SnapWave
- Made clearer by putting the read functions in new filed 'sfincs_read', and only call those in sfincs_snapwave.f90

* - Also in sfincs_spiderweb.f90 replace local read_real_input and read_int_input by use sfincs_read

* - add posibility to inspect snapwave grid as in PR jre/snapwave_sfincs

* - restart already fixed to true in sfincs_snapwave

* - Remove unused subroutine 'neuboundaries' that is replaced by neuboundaries_light

* - Add back some description

* - Bugfix, actually update Hk so it is used in IG bottom friction determination

* - Preprocess snapwave_ncoutput

* - Keep snapwave_crit= 0.001 as default as in main

* - make (:,k) to (:, k) consistently, everywhere I hope

* - Try adding back missing stuff

* making it work for wind

* removed second call to celerities

* - Bump version for clarity when testing

* - Make Herbers IG boundary warning dimensionless, using hsinc / depth
- Check further for best ratio later, start at 0.5

* - Add user definable snapwave_relax_factor_DoverE and snapwave_relax_factor_DoverA keywords to sfincs_snapwave
- By default set to 0.25 for both for now

* - For now do simple Hmx(k)    = gamma * depth(k) and IG as in trunk, TBD

* - As discussed, remove fdrspr for now. Could be added back later on (with actual working function)

* - Add back comment

* - Remove subroutine bj78

* - Added calculation and write_log message of error_ig and %ok_ig with thanks to Mr Claude
- For now criterion 'ok' of incident waves for determining iteration ended is unaltered

* - Add omp parallellisation

* - Already do the bugfix of PR #283 by switching to srcig based on E and redistribute by /E * ee

* - Add back Neumann message

* - Bugfix in calculating kwav_ig!

* - Suggestion of Marlies to keep the original formulation for Hmx_ig

* Reinit Tp for inner/Neumann; ignore docs/_build

Reinitialize the per-node Tp only for internal and Neumann-connected cells during the wind iteration so boundary cells keep their prescribed Tp (as set by update_boundaries). Implemented a loop that sets Tp(k)=Tpini for inner(k) and Tp(neumannconnected(k)) where present. Also added docs/_build to .gitignore and cleaned up minor whitespace.

* - Add gamma_fac_br as in PR #281

* 331 adapt wave force limiter for lab tests (#333)

* - issue for Buckley solved if no limiter at all
- Later improve limiter

* - Temp save, need to still make it work here/in sfincs_snapwave

* - Option to determine fwmaxfac once every 'update_wave_field' call in sfincs_snapwave

* 292 renewed branch vegetation effects on short waves and wave forces (#296)

* merged veggie shortwave drag branch143 into current main

* added missing vfproj update

* - First step to making mean flow correct on the uv points

* - clean repo

* - Clean up this PR, so that only the pure snapwave veggie changes remain, the mean flow effect is still kept in copy of this branch: 315-veggie-effect-on-mean-flow

* - Cleanup and switch to logicals

* - Start of sfincs_vegetation

* - Only for quadtree

* - Move veggie reference out of sfincs_quadtree

* - Add new generic subroutine 'read_netcdf_quadtree_get_dimension'
- cleanup

* - Neat way to import only needed veggie variables from sfincs_data into snapwave_domain

* - Cleanup code

* - Manual merge back swvegnonlin, momeqveg, update F(k)

* - For now limit vegetation_vertical_segments to max 4

* - move check from sfincs_snapwave to wavemaker (not veggie related), if snapwave yes but wavemaker no, the message before still appeared

* Add toml-f and TOML src_structure parsing

Add the toml-f third-party TOML parser (sources and licenses) and wire it into the build (Makefile.am and Visual Studio .vfproj). Implement TOML-based parsing for src_structure inputs in sfincs_src_structures.f90: new types, constants, read_toml_src_structures routine, validation helpers and a dispatcher in initialize_src_structures() that probes for TOML and falls back to the legacy reader. Also include minor non-functional whitespace/formatting cleanups in sfincs_discharges.f90.

* - Add TOML-based vegetation type lookup table with CF flag NetCDF input

Vegetation parameters can now be specified via a two-file approach:
- The NetCDF file (vegetationfile) stores only an integer vegetation_type
  per grid cell, following CF flag_values / flag_meanings conventions so
  the file is self-documenting in standard tools (ncview, Panoply, xarray).
- A companion TOML file (vegfile_toml) contains the parameter lookup table
  keyed by type name, supporting any number of vertical layers per type
  (mixed layer counts across types are handled by zero-padding).

The old single-NetCDF code path is fully preserved for backward
compatibility; the new path is activated only when vegfile_toml is set
in sfincs.inp.

Changes:
- sfincs_vegetation.f90: add read_vegetation_toml(); rewrite
  initialize_vegetation() with old/new path split; lift 4-layer cap to 64
- sfincs_ncinput.F90: add read_netcdf_quadtree_integer() and
  read_netcdf_flag_meanings() for the new input format
- sfincs_data.f90: add veggiefile_toml, vegetation_type_index; rename
  vegetation_cd -> vegetation_stems_cd and
  vegetation_stems_width -> vegetation_stems_diameter for consistency
- sfincs_input.f90: add vegfile_toml keyword
- snapwave_domain.f90: update to renamed variable names

* - vegetationtype_toml as input name - to discuss

* - Visual cleanup

* - Bump version

* - missing stop_sfincs reference

---------

Co-authored-by: Leynse <Tim.Leijnse@deltares.nl>
Co-authored-by: Tim Leijnse <timleynse94@gmail.com>
Co-authored-by: Maarten van Ormondt <maarten.vanormondt@deltares.nl>

---------

Co-authored-by: Leynse <Tim.Leijnse@deltares.nl>
Co-authored-by: Tim Leijnse <timleynse94@gmail.com>
Co-authored-by: Marlies <marlies.vanderlugt@deltares.nl>

* - Bump version to release candidate

* Fix directional asymmetry in wavemaker initialization for south-facing wavemakers

Two bugs in initialize_wavemakers() that caused incorrect or zero forcing for wavemakers with phi in the south half-plane (Q3/Q4):

Q3 validity check (phi ∈ [π, 1.5π)): wrong index passed to uv_index_z_nm — used Z-cell index ip instead of UV index nmu for the first left-neighbor (md1) check. This caused nearly all Q3 candidates to pass the validity filter, producing spurious wavemaker points for south-west-facing lines.

"Below" direction angfac (Q3 nd1/nd2, Q4 nd1/nd2): sin(π − φ) ≡ sin(φ) and sin(φ) ≤ 0 for φ ∈ [π, 2π), so max(sin(π − φ), 0) always returned zero, silencing all south-pointing wavemaker UV faces. Fixed to max(−sin(φ), 0).

* - Please double check this Claude bugfix #3 @maarten !

Fix wavemaker side-pointer not set for refined (*2) UV half-edges

On quadtree-refined grids, each coarse cell side can have two UV faces (*1 and *2). Only the *1 face was stored in wavemaker_nmu/nmd/num/ndm, so the *2 face injected flux without the continuity solver tracking it, causing mass drift near the wavemaker. No effect on uniform grids.

* - only write hm0_ig min/mean/max for wavemaker points

* - Last small fixes

* - Remove write statements

* Warn on wavemaker points at refinement boundary

Add a logical flag to detect wavemaker points adjacent to quadtree refinement and emit a log warning. Declares refinement_warning, initializes it to .false., sets it to .true. in the wavemaker point checks (where kcs(iz)==1), and writes a WARNING message to the log if any such cases are found. This helps notify users that wavemaker points along refinement boundaries are not recommended.

* - Cleanup message

* - Add netcdf output of the veggie variables (change names later)

* - Add input variable 'snapwave_waveforces_ratio' which you can set to 0 to turn off wave forces and thus incident wave setup

* - Update names

* - Fix updated variable names

* - Not used anymore

* - snapwave_gamma_fac_br

* - remove unused RFveg.inc

* - cleanup

* - consistent naming

* - consistently as diameter

* - Set snapwave_baldock_exponent = 2 as the new default. Generally, only active for steep coastlines, where Baldock dissipation can be too low in the surf zone.

* Change default snapwave_gammax to 999 now, which works beceause of snapwave_baldock_exponent = 2

* - Bump version number license

* -Wrong order of checking
-Bump version

* Fix/bugs 11 12 13 16 (#325)

* Fix #11: guard tsrc index to prevent out-of-bounds at t < tsrc(1)

* Fix #12: skip drainage points outside grid before any index use

* Fix #13: clamp spiderweb bilinear weight dj1 >= 0 near cyclone eye

* Fix #16: prevent Horton rain_T1 going negative at storm onset

* Revert "Fix #11: guard tsrc index to prevent out-of-bounds at t < tsrc(1)"
* - Think better solution to check times
- Now it applies also to netcdf srcfile input

* - Update docs developments figure, release description

* - Bump version

* - Add slogan "Generating Accurate Large-scale Inundation: Better Insights for Emergency Response"

* - Change version name

* - Update for local build

* - Fix some linux compiler issues

* - Attempt to fix nvfortran compilation errors of tomlf library

* - More fixes for nvfortran

* - Last needed edits - Claude

* with some extra info on the epsg_code, renaming of variables and long_names, qgis can immediatly read the mesh2d output! (#345)

Co-authored-by: Roel de Goede <roel.degoede@deltares.nl>

* - Bump version

* - Also bump SnapWave version

* - Fix compilation on Windows

* - Update figure including the Qgis netcdf output for quadtree

* - Remove the overwriting values to 0.0 in case input values are below -99 or above 999.
- The warning is kept in tact as is

* - Bump version

* - Add SFINCS development status table for the first time, not release dependent

* - Add QGIS ready netcdf output file milestone

* - Update readme

* - Fix date

* - Add to docs

---------

Co-authored-by: Maarten van Ormondt <82946105+maartenvanormondt@users.noreply.github.com>
Co-authored-by: Marlies <marlies.vanderlugt@deltares.nl>
Co-authored-by: Maarten van Ormondt <maarten.vanormondt@deltares.nl>
Co-authored-by: Roel de Goede <roel.degoede@deltares.nl>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Add version of IG source term where dSxxdx does not depend on the bins

1 participant