BUG: streamlines only returns the dx_array, dy_array of last frame

[lilyminium]

Expected behavior

From a skim over the docs and code my impression was that streamlines.generate_streamlines would incorporate all frames between start_frame and end_frame in the dx_array, but looking more closely at the code I'm not sure that happens. Is it supposed to?

Actual behavior

I can get this test to pass with start_frame=0, start_frame=1, and start_frame=2.

mdanalysis/testsuite/MDAnalysisTests/visualization/test_streamlines.py

Lines 54 to 80 in 6ec838b

	def test_streamplot_2D(membrane_xtc, univ):
	# regression test the data structures
	# generated by the 2D streamplot code
	u1, v1, avg, std = streamlines.generate_streamlines(topology_file_path=Martini_membrane_gro,
	trajectory_file_path=membrane_xtc,
	grid_spacing=20,
	MDA_selection='name PO4',
	start_frame=1,
	end_frame=2,
	xmin=univ.atoms.positions[...,0].min(),
	xmax=univ.atoms.positions[...,0].max(),
	ymin=univ.atoms.positions[...,1].min(),
	ymax=univ.atoms.positions[...,1].max(),
	maximum_delta_magnitude=2.0,
	num_cores=1)
	assert_allclose(u1, np.array([[0.79999924, 0.79999924, 0.80000687, 0.79999542, 0.79998779],
	[0.80000019, 0.79999542, 0.79999924, 0.79999542, 0.80001068],
	[0.8000021, 0.79999924, 0.80001068, 0.80000305, 0.79999542],
	[0.80000019, 0.79999542, 0.80001068, 0.80000305, 0.80000305],
	[0.79999828, 0.80000305, 0.80000305, 0.80000305, 0.79999542]]))
	assert_allclose(v1, np.array([[0.53999901, 0.53999996, 0.53999996, 0.53999996, 0.54000092],
	[0.5399971, 0.54000092, 0.54000092, 0.54000092, 0.5399971 ],
	[0.54000473, 0.54000473, 0.54000092, 0.5399971, 0.54000473],
	[0.54000092, 0.53999329, 0.53999329, 0.53999329, 0.54000092],
	[0.54000092, 0.53999329, 0.53999329, 0.54000092, 0.53999329]]))
	assert avg == pytest.approx(0.965194167)
	assert std == pytest.approx(4.444808820e-06)

Because xy_deltas_to_write is returned from this function, it looks like the results get re-set for each frame:

mdanalysis/package/MDAnalysis/visualization/streamlines.py

Lines 201 to 214 in 6ec838b

	xy_deltas_to_write = []
	for square_1_centroid, square_2_centroid in zip(list_centroids_this_frame_using_indices_from_last_frame,
	list_previous_frame_centroids):
	if square_1_centroid is None or square_2_centroid is None:
	xy_deltas_to_write.append([0, 0])
	else:
	xy_deltas_to_write.append(np.subtract(square_1_centroid, square_2_centroid).tolist())
	#xy_deltas_to_write = np.subtract(np.array(
	# list_centroids_this_frame_using_indices_from_last_frame),np.array(list_previous_frame_centroids))
	xy_deltas_to_write = np.array(xy_deltas_to_write)
	#now filter the array to only contain distances in the range [-8,8] as a placeholder for dealing with PBC
	# issues (Matthieu seemed to use a limit of 8 as well);
	xy_deltas_to_write = np.clip(xy_deltas_to_write, -maximum_delta_magnitude, maximum_delta_magnitude)

Current version of MDAnalysis

• Which version are you using? (run python -c "import MDAnalysis as mda; print(mda.__version__)")
• Which version of Python (python -V)?
• Which operating system?

[orbeckst]

@tylerjereddy could you please spare a quick glance 👀 here?

[tylerjereddy]

I saw the ping, I'll try to check soon-ish. I'm under some pressure to get SciPy released to support NumPy 2 series, but will try to make some time.

[tylerjereddy]

I agree with Lily that there's a bug here.

The lack of input validation is also pretty jarring, I can give an obscure floating point value to the test start frame value and it still passes (pretty silly):

--- a/testsuite/MDAnalysisTests/visualization/test_streamlines.py
+++ b/testsuite/MDAnalysisTests/visualization/test_streamlines.py
@@ -88,7 +88,7 @@ def test_streamplot_2D_zero_return(membrane_xtc, univ, tmpdir):
                                                         trajectory_file_path=membrane_xtc,
                                                         grid_spacing=20,
                                                         MDA_selection='name POX',
-                                                        start_frame=1,
+                                                        start_frame=np.sqrt(17),
                                                         end_frame=2,
                                                         xmin=univ.atoms.positions[...,0].min(),
                                                         xmax=univ.atoms.positions[...,0].max(),

So yeah, definitely some issues here. Not sure that Lily should be burdened with fixing them all at this time, but definitely good to have issue(s) open so folks can be aware.

[orbeckst]

@tylerjereddy is this bug something you could fix?

[tylerjereddy]

Probably, it sounds like Lily may fix it during the refactor to AnalysisBase in any case though.