From 0b77526af0588604cea15433e05be41fe592c7f3 Mon Sep 17 00:00:00 2001
From: "Leaf, Andrew T" <aleaf@usgs.gov>
Date: Fri, 3 Jan 2025 14:11:29 -0600
Subject: [PATCH] fix(mover.py::get_sfr_package_connections): when two outlets
 are connected across a LGR model interface and the connection threshold
 distance is large, a circular or double connection going both ways between
 the two models may result, which can prevent the Mover Package from
 converging in the simulation. This issue was fixed by adding a check for
 double connections, and then only retaining the connection with the smallest
 distance (which should always be correct if the flowline input to the SFR
 Package was digitized correctly).

---
 mfsetup/mover.py            | 32 ++++++++++++++++++++++++++++++--
 mfsetup/tests/test_lgr.py   |  8 ++++++++
 mfsetup/tests/test_mover.py |  4 ++++
 3 files changed, 42 insertions(+), 2 deletions(-)

diff --git a/mfsetup/mover.py b/mfsetup/mover.py
index 2c9e12b2..bb61182a 100644
--- a/mfsetup/mover.py
+++ b/mfsetup/mover.py
@@ -103,6 +103,10 @@ def neighbors(i, j):
     # make the connections
     parent_to_inset = dict()
     inset_to_parent = dict()
+    # record connection distances
+    # for breaking circular routing cases
+    parent_to_inset_distances = dict()
+    inset_to_parent_distances = dict()
     for _, r in reach_data1.iterrows():
         # check for connections to the other model
         # if the next reach is in another cell
@@ -120,9 +124,11 @@ def neighbors(i, j):
             # consider this each to be an outlet
             reach_end = Point(r['geometry'].coords[-1])
             distances = reach_end.distance(reach_data2['reach_start'])
-            if np.min(distances) < distance_threshold:
+            min_distance = np.min(distances)
+            if min_distance < distance_threshold:
                 next_reach = reach_data2.iloc[np.argmin(distances)]
                 parent_to_inset[r['rno']] = next_reach['rno']
+                parent_to_inset_distances[r['rno']] = min_distance
         # next reach is somewhere else in the parent model
         else:
             continue
@@ -144,12 +150,34 @@ def neighbors(i, j):
             # consider this each to be an outlet
             reach_end = Point(r['geometry'].coords[-1])
             distances = reach_end.distance(reach_data1['reach_start'])
-            if np.min(distances) < distance_threshold:
+            min_distance = np.min(distances)
+            if min_distance < distance_threshold:
                 next_reach = reach_data1.iloc[np.argmin(distances)]
                 inset_to_parent[r['rno']] = next_reach['rno']
+                inset_to_parent_distances[r['rno']] = min_distance
         # next reach is somewhere else in this model
         else:
             continue
+    # check for circular connections (going both ways between two models)
+    # retain connection with the smallest distance
+    delete_parent_to_inset_items = set()
+    for parent_reach, inset_reach in parent_to_inset.items():
+        if parent_reach in inset_to_parent.values():
+            parent_to_inset_distance = parent_to_inset_distances[parent_reach]
+            inset_to_parent_distance = inset_to_parent_distances[inset_reach]
+            if inset_to_parent_distance < parent_to_inset_distance:
+                delete_parent_to_inset_items.add(parent_reach)
+            elif parent_to_inset_distance < inset_to_parent_distance:
+                del inset_to_parent[inset_reach]
+            else:
+                raise ValueError("Circular connection between SFR Packages in the Mover Package input.\n"
+                                 f"Connection distance between the end of parent reach {parent_reach} "
+                                 f"in parent model and start of inset reach {inset_reach} in inset model "
+                                 f"is equal to\nthe distance between the end of inset reach {inset_reach} "
+                                 f"and start of parent reach {parent_reach}.\nCheck input linework."
+                                 )
+    for parent_reach in delete_parent_to_inset_items:
+        del parent_to_inset[parent_reach]
     return parent_to_inset, inset_to_parent
 
 
diff --git a/mfsetup/tests/test_lgr.py b/mfsetup/tests/test_lgr.py
index dddf829a..665ae7da 100644
--- a/mfsetup/tests/test_lgr.py
+++ b/mfsetup/tests/test_lgr.py
@@ -328,6 +328,14 @@ def test_mover_get_sfr_package_connections(pleasant_lgr_setup_from_yaml):
     # {inset_reach: parent_reach, ...}
     assert to_parent == {29: 13, 41: 1}
 
+    # test for no circular connections when two outlets are connected
+    # and distance_threshold is large
+    parent_reach_data.loc[parent_reach_data['rno'] == list(to_parent.values())[0], 'outreach'] = 0
+    to_inset, to_parent = get_sfr_package_connections(
+    gwfgwf_exchangedata,
+    parent_reach_data, inset_reach_data, distance_threshold=1e4)
+    assert not any(to_inset)
+
 
 def test_meandering_sfr_connections(shellmound_cfg, project_root_path, tmpdir):
     """Test for SFR routing continuity in LGR cases
diff --git a/mfsetup/tests/test_mover.py b/mfsetup/tests/test_mover.py
index 950a3787..26cb6b3d 100755
--- a/mfsetup/tests/test_mover.py
+++ b/mfsetup/tests/test_mover.py
@@ -1,3 +1,7 @@
+"""Test functions in the mover.py module. See test_lgr.py for
+a test of the test_mover_get_sfr_package_connections function.
+
+"""
 from copy import deepcopy
 
 import pytest