Modpath Timeseries Class way of work

[javgs-bd]

From a steady state simulation (like that one of Example 01a in the Modpath documentation) I want to plot particles positions at different time instants. For that simulation 21 particles are released and then tracked in the time series simulation at time intervals of 1000 days.

The ex01a_mf6.timeseries7 files looks like this:

It is sorted following a time instant grouping philosophy. However, when loaded with flopy using the TimeSeriesFile Class it gets reordered in a recarray containing one entry per particle id. That recarray entry contains all the times in the timeseries for that particle.

So, flopy is giving me a hard time to plot those 21 particles at the different time instants stored in the timeseries files. First, because the timeseries array read with 'flopy.utils.TimeseriesFile' is structuredthat way (kind of thinking on per particle IDs rather than in per time instant). Second, because (from inspecting the source code) all that plot_timeseries does is to parse a travel_time argument to plot_pathlines. Then the latter plots all particles at all time instants in the timeseries that are lower equal than the travel_time. Thus I always get past times plotted together with the one I want.

I'll try to build a minimum working example (the week run out on me), but the following code snippet could explain a bit of the problem. It already has a modflow model object (gwf), modpath object and modpath ouput files have been created.

ts_file     = os.path.join(mp_folder, mp_modelname + ".timeseries")
ts           = flopy.utils.TimeseriesFile(ts_file)
ts_array = ts.get_alldata()
       
times =  np.unique(ts._data.time)
      
norm_ts  = plt.Normalize()
ts_colors = plt.cm.jet(norm_ts(times))
        
fig, ax = plt.subplots(
                        nrows       = 1,
                        ncols       = 1,
                        figsize     = (5, 9),
                        constrained_layout = True
                       )
        
mview  = flopy.plot.PlotMapView(ax = ax, model = gwf, layer = 0) 
mview.plot_grid(color = (0.4,0.4,0.4,0.2), lw = 1)
mview.plot_timeseries(ts_array, travel_time=times[5], layer="all", marker='o', markercolor= ts_colors[5], lw=0)

I had to use ts._data to identify the unique time instants in the timeseries file. That way, the code snippet would result in all particles being plotted for all times in the timeseries before the fifth time instant in the timeseries file, instead of only those particles present in the model at that time instant.

Am I missing some obvious way to plot particles for a specific time instant in the timeseries file?

I'll post the example as soon as I can, but I wonder if the TimeSeriesFile class and plot_timeseries should not be revised, to be aimed to time instants instead of individual particles.

[jlarsen-usgs]

@javgs-bd

I just looked into this and the travel_time argument does not currently support a logical equal. I've put in a PR to fix this issue and add support for this; PR #1736.

Once this PR is accepted, you'll be able to filter particles by a specific time instance:

example usage:

mview.plot_timeseries(ts_array, travel_time=f"== {times[5]}", layer="all", marker='o', markercolor= ts_colors[5], lw=0)

[javgs-bd]

Thanks @jlarsen-usgs. That'll be just great. Since we are here, I'll dare to ask. I guess a method to get back the timeseries time instants would be the only thing missing.

[jlarsen-usgs]

Yes, it must've been overlooked when we updated the plotting methods. A simple lamba function was the only thing that needed to be added. The PR has been approved and merged, so if you update your flopy with the development branch you'll have access to this functionality.

[javgs-bd]

Many thanks again!