Skip to content

Commit

Permalink
fix: water and sediment routing bugs that caused inconsistent behavior
Browse files Browse the repository at this point in the history
There were two issues with the refactoring. One was that water routing 
developed a flat flow field because I wasn't adding each passing water 
parcel to the Q field. The other was related to sediment parcel routing, 
where I had transcribed the i and j of stepping incorrectly in the 
jitted functions.
  • Loading branch information
ericbarefoot committed May 26, 2020
1 parent 10025fa commit fd238c6
Show file tree
Hide file tree
Showing 3 changed files with 15 additions and 10 deletions.
6 changes: 3 additions & 3 deletions pyDeltaRCM/sed_tools.py
Original file line number Diff line number Diff line change
Expand Up @@ -207,14 +207,14 @@ def sed_parcel(self, theta_sed, sed, px, py):

depoPart = self.Vp_res / 2 / self.dt / self.dx

px, py, self.qs[px, py] = shared_tools.partition_sand(self.qs[px, py], depoPart, py, px, dist, istep, jstep)
px, py, self.qs = shared_tools.partition_sand(self.qs, depoPart, py, px, dist, istep, jstep)

self.sand_dep_ero(px, py)

if sed == 'mud': # mud

px = px + istep
py = py + jstep
px = px + jstep
py = py + istep

self.mud_dep_ero(px, py)

Expand Down
17 changes: 11 additions & 6 deletions pyDeltaRCM/shared_tools.py
Original file line number Diff line number Diff line change
Expand Up @@ -62,13 +62,14 @@ def get_filtered_weight(weight, px, depth_ind, cell_type_ind, dry_depth):
def partition_sand(qs, depoPart, py, px, dist, istep, jstep):
if dist > 0:
# deposition in current cell
qs = qs + depoPart
px = px + istep
py = py + jstep
qs[px, py] += depoPart

px = px + jstep
py = py + istep

if dist > 0:
# deposition in downstream cell
qs = qs + depoPart
qs[px, py] += depoPart
return px, py, qs


Expand All @@ -86,13 +87,17 @@ def get_steps(new_cells, iwalk, jwalk):

@njit
def update_dirQfield(qfield, dist, inds, astep, dirstep):
qfield[inds[astep]] += dirstep[astep] / dist[astep]
for i, ii in enumerate(inds):
if astep[i]:
qfield[ii] += dirstep[i] / dist[i]
return qfield


@njit
def update_absQfield(qfield, dist, inds, astep, Qp_water, dx):
qfield[inds[astep]] += Qp_water / dx / 2
for i, ii in enumerate(inds):
if astep[i]:
qfield[ii] += Qp_water / dx / 2
return qfield


Expand Down
2 changes: 1 addition & 1 deletion pyDeltaRCM/water_tools.py
Original file line number Diff line number Diff line change
Expand Up @@ -89,7 +89,7 @@ def run_water_iteration(self):

current_inds, self.looped, self.free_surf_flag = shared_tools.check_for_loops(
self.indices,
current_inds,
next_index,
iter,
self.L0,
self.looped,
Expand Down

0 comments on commit fd238c6

Please sign in to comment.