Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Diffusive Data Assimilation v2 with USGS streamflow data #552

Open
wants to merge 7 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
280 changes: 182 additions & 98 deletions src/kernel/diffusive/diffusive.f90

Large diffs are not rendered by default.

39 changes: 21 additions & 18 deletions src/kernel/diffusive/pydiffusive.f90
Original file line number Diff line number Diff line change
Expand Up @@ -5,43 +5,46 @@ module diffusive_interface

implicit none
contains
subroutine c_diffnw(timestep_ar_g, nts_ql_g, nts_ub_g, nts_db_g, ntss_ev_g, &
nts_qtrib_g, mxncomp_g, nrch_g, z_ar_g, bo_ar_g, traps_ar_g, &
tw_ar_g, twcc_ar_g, mann_ar_g, manncc_ar_g, so_ar_g, dx_ar_g, &
iniq, frnw_col, frnw_g, qlat_g, ubcd_g, dbcd_g, qtrib_g, &
paradim, para_ar_g, mxnbathy_g, x_bathy_g, z_bathy_g, &
mann_bathy_g, size_bathy_g, q_ev_g, elv_ev_g) bind(c)
subroutine c_diffnw(timestep_ar_g, nts_ql_g, nts_ub_g, nts_db_g, ntss_ev_g, nts_qtrib_g, nts_da_g, &
mxncomp_g, nrch_g, z_ar_g, bo_ar_g, traps_ar_g, tw_ar_g, twcc_ar_g, mann_ar_g, &
manncc_ar_g, so_ar_g, dx_ar_g, &
iniq, frnw_col, frnw_ar_g, qlat_g, ubcd_g, dbcd_g, qtrib_g, &
paradim, para_ar_g, mxnbathy_g, x_bathy_g, z_bathy_g, mann_bathy_g, size_bathy_g, &
usgs_da_g, usgs_da_reach_g, q_ev_g, elv_ev_g) bind(c)


integer(c_int), intent(in) :: nts_ql_g, nts_ub_g, nts_db_g, nts_qtrib_g
integer(c_int), intent(in) :: nts_ql_g, nts_ub_g, nts_db_g, nts_qtrib_g, nts_da_g
integer(c_int), intent(in) :: ntss_ev_g
integer(c_int), intent(in) :: mxncomp_g, nrch_g
integer(c_int), intent(in) :: frnw_col
integer(c_int), intent(in) :: paradim
integer(c_int), intent(in) :: mxnbathy_g
integer(c_int), dimension(nrch_g, frnw_col), intent(in) :: frnw_g
integer(c_int), intent(in) :: mxnbathy_g
integer(c_int), dimension(nrch_g), intent(in) :: usgs_da_reach_g
integer(c_int), dimension(nrch_g, frnw_col), intent(in) :: frnw_ar_g
integer(c_int), dimension(mxncomp_g, nrch_g), intent(in) :: size_bathy_g
real(c_double), dimension(nts_db_g), intent(in) :: dbcd_g
real(c_double), dimension(paradim), intent(in) :: para_ar_g
real(c_double), dimension(:), intent(in) :: timestep_ar_g(8)
real(c_double), dimension(:), intent(in) :: timestep_ar_g(9)
real(c_double), dimension(mxncomp_g, nrch_g), intent(in) :: z_ar_g, bo_ar_g, traps_ar_g
real(c_double), dimension(mxncomp_g, nrch_g), intent(in) :: tw_ar_g, twcc_ar_g
real(c_double), dimension(mxncomp_g, nrch_g), intent(in) :: mann_ar_g, manncc_ar_g
real(c_double), dimension(mxncomp_g, nrch_g), intent(inout) :: so_ar_g
real(c_double), dimension(mxncomp_g, nrch_g), intent(in) :: dx_ar_g, iniq
real(c_double), dimension(nts_ub_g, nrch_g), intent(in) :: ubcd_g
real(c_double), dimension(nts_qtrib_g, nrch_g), intent(in) :: qtrib_g
real(c_double), dimension(nts_da_g, nrch_g), intent(in) :: usgs_da_g
real(c_double), dimension(nts_ql_g, mxncomp_g, nrch_g), intent(in) :: qlat_g
real(c_double), dimension(mxnbathy_g, mxncomp_g, nrch_g), intent(in ) :: x_bathy_g
real(c_double), dimension(mxnbathy_g, mxncomp_g, nrch_g), intent(in ) :: z_bathy_g
real(c_double), dimension(mxnbathy_g, mxncomp_g, nrch_g), intent(in ) :: mann_bathy_g
real(c_double), dimension(ntss_ev_g, mxncomp_g, nrch_g), intent(out) :: q_ev_g, elv_ev_g

call diffnw(timestep_ar_g, nts_ql_g, nts_ub_g, nts_db_g, ntss_ev_g, &
nts_qtrib_g, mxncomp_g, nrch_g, z_ar_g, bo_ar_g, traps_ar_g, &
tw_ar_g, twcc_ar_g, mann_ar_g, manncc_ar_g, so_ar_g, dx_ar_g, &
iniq, frnw_col, frnw_g, qlat_g, ubcd_g, dbcd_g, qtrib_g, &
paradim, para_ar_g, mxnbathy_g, x_bathy_g, z_bathy_g, &
mann_bathy_g, size_bathy_g, q_ev_g, elv_ev_g)
real(c_double), dimension(ntss_ev_g, mxncomp_g, nrch_g), intent(out) :: q_ev_g, elv_ev_g
call diffnw(timestep_ar_g, nts_ql_g, nts_ub_g, nts_db_g, ntss_ev_g, nts_qtrib_g, nts_da_g, &
mxncomp_g, nrch_g, z_ar_g, bo_ar_g, traps_ar_g, tw_ar_g, twcc_ar_g, mann_ar_g, &
manncc_ar_g, so_ar_g, dx_ar_g, &
iniq, frnw_col, frnw_ar_g, qlat_g, ubcd_g, dbcd_g, qtrib_g, &
paradim, para_ar_g, mxnbathy_g, x_bathy_g, z_bathy_g, mann_bathy_g, size_bathy_g, &
usgs_da_g, usgs_da_reach_g, q_ev_g, elv_ev_g)

end subroutine c_diffnw
end module diffusive_interface
10 changes: 6 additions & 4 deletions src/troute-network/troute/nhd_network_utilities_v02.py
Original file line number Diff line number Diff line change
Expand Up @@ -511,7 +511,7 @@ def build_da_sets(da_params, run_sets, t0):
"usace_timeslices_folder",
None
)

# User-specified DA ON/OFF preferences
usace_da = False
usgs_da = False
Expand Down Expand Up @@ -552,7 +552,7 @@ def build_da_sets(da_params, run_sets, t0):
dt_timeslice = timedelta(minutes = 15)

da_sets = [] # initialize list to store TimeSlice set lists

# Loop through each run set and build lists of available TimeSlice files
for (i, set_dict) in enumerate(run_sets):

Expand Down Expand Up @@ -753,8 +753,10 @@ def build_data_assimilation_lastobs(data_assimilation_parameters):
da_parameter_dict = {}
da_parameter_dict["da_decay_coefficient"] = data_assimilation_parameters.get(
"da_decay_coefficient",
120
)
120)
da_parameter_dict["diffusive_streamflow_nudging"] = streamflow_da_parameters.get(
"diffusive_streamflow_nudging", False)

# TODO: Add parameters here for interpolation length (14/59), QC threshold (1.0)

return lastobs_df, da_parameter_dict
Expand Down
3 changes: 2 additions & 1 deletion src/troute-nwm/src/nwm_routing/__main__.py
Original file line number Diff line number Diff line change
Expand Up @@ -127,13 +127,14 @@ def nwm_route(

LOG.debug("MC computation complete in %s seconds." % (time.time() - start_time_mc))
start_time_diff = time.time()

# call diffusive wave simulation and append results to MC results
results.extend(
compute_diffusive_routing(
results,
diffusive_network_data,
cpu_pool,
t0,
dt,
nts,
q0,
Expand Down
2 changes: 1 addition & 1 deletion src/troute-nwm/src/nwm_routing/preprocess.py
Original file line number Diff line number Diff line change
Expand Up @@ -152,7 +152,7 @@ def nwm_network_preprocess(
connections, param_df, wbody_conn, gages = nnu.build_connections(
supernetwork_parameters,
)

link_gage_df = pd.DataFrame.from_dict(gages)
link_gage_df.index.name = 'link'

Expand Down
14 changes: 11 additions & 3 deletions src/troute-routing/troute/routing/compute.py
Original file line number Diff line number Diff line change
Expand Up @@ -1096,6 +1096,7 @@ def compute_diffusive_routing(
results,
diffusive_network_data,
cpu_pool,
t0,
dt,
nts,
q0,
Expand All @@ -1107,8 +1108,7 @@ def compute_diffusive_routing(
diffusive_parameters,
waterbodies_df,
topobathy_data,
):

):
results_diffusive = []
for tw in diffusive_network_data: # <------- TODO - by-network parallel loop, here.

Expand All @@ -1131,6 +1131,12 @@ def compute_diffusive_routing(
topobathy_data_bytw = topobathy_data.loc[diffusive_network_data[tw]['mainstem_segs']]
else:
topobathy_data_bytw = pd.DataFrame()

# diffusive streamflow DA activation switch
if da_parameter_dict['diffusive_streamflow_nudging']==True:
diff_usgs_df = usgs_df
else:
diff_usgs_df = pd.DataFrame()

# build diffusive inputs
diffusive_inputs = diff_utils.diffusive_input_data_v02(
Expand All @@ -1145,15 +1151,17 @@ def compute_diffusive_routing(
q0,
junction_inflows,
qts_subdivisions,
t0,
nts,
dt,
waterbodies_df,
topobathy_data_bytw,
diff_usgs_df,
)

# run the simulation
out_q, out_elv = diffusive.compute_diffusive(diffusive_inputs)

# unpack results
rch_list, dat_all = diff_utils.unpack_output(
diffusive_inputs['pynw'],
Expand Down
Loading