[SCM] Gerris Flow Solver branch, upstream, updated. b3aa46814a06c9cb2912790b23916ffb44f1f203

[Stephane Popinet]

The following commit has been merged in the upstream branch:
commit 7cae0f7dd19fcf6e3cde55a86bcc852d4ff12bc0
Author: Stephane Popinet <popinet at users.sourceforge.net>
Date:   Wed Nov 3 16:10:36 2004 +1100

    New function gfs_domain_remove_specks() (gerris--mainline--0.7--patch-26)
    
    gerris--mainline--0.7--patch-26
    Keywords:
    
    "Specks" are solid islands not containing any empty cell.
    
    They cause stability problems for the free-surface solver of the ocean
    model.
    
    They cannot be considered resolved anyway (by definition their maximum
    thickness is smaller than one grid cell).
    
    darcs-hash:20041103051036-aabb8-4a76a36de12d15be9d098d1ae116312de92bd522.gz

diff --git a/doc/gfs-sections.txt b/doc/gfs-sections.txt
index cf97fff..4caaf5e 100644
--- a/doc/gfs-sections.txt
+++ b/doc/gfs-sections.txt
@@ -324,6 +324,7 @@ gfs_domain_split
 gfs_domain_variable_sources
 gfs_domain_remove_droplets
 gfs_domain_remove_ponds
+gfs_domain_remove_specks
 gfs_domain_advect_point
 gfs_domain_cfl
 gfs_domain_size
diff --git a/src/domain.c b/src/domain.c
index c61dfc6..6f2f4cf 100644
--- a/src/domain.c
+++ b/src/domain.c
@@ -2569,6 +2569,86 @@ void gfs_domain_remove_ponds (GfsDomain * domain,
   g_array_free (sizes, TRUE);
 }
 
+static gboolean tag_speck (FttCell * cell)
+{
+  if (GFS_STATE (cell)->div == 0.) {
+    FttDirection d;
+    FttCellNeighbors n;
+    GfsSolidVector * solid = GFS_STATE (cell)->solid;
+    
+    g_assert (FTT_CELL_IS_LEAF (cell));
+    ftt_cell_neighbors (cell, &n);
+    for (d = 0; d < FTT_NEIGHBORS; d++)
+      if (!n.c[d])
+	return FALSE;
+    GFS_STATE (cell)->div = 1.;
+    for (d = 0; d < FTT_NEIGHBORS; d++)
+      if (GFS_STATE (n.c[d])->div == 0. && 
+	  !GFS_CELL_IS_BOUNDARY (n.c[d]) &&
+	  solid->s[d] > 0. && solid->s[d] < 1.) {
+	g_assert (GFS_IS_MIXED (n.c[d]));
+	if (FTT_CELL_IS_LEAF (n.c[d])) {
+	  if (!tag_speck (n.c[d])) {
+	    GFS_STATE (cell)->div = 0.;
+	    return FALSE;
+	  }
+	}
+	else {
+	  FttCellChildren child;
+	  FttDirection od = FTT_OPPOSITE_DIRECTION (d);
+	  guint i;
+	  
+#if FTT_2D3
+	  g_assert_not_implemented ();
+#endif	
+	  ftt_cell_children_direction (n.c[d], od, &child);
+	  for (i = 0; i < FTT_CELLS/2; i++)
+	    if (!child.c[i] || (GFS_STATE (child.c[i])->div == 0 && 
+				GFS_IS_MIXED (child.c[i]) &&
+				!tag_speck (child.c[i]))) {
+	      GFS_STATE (cell)->div = 0.;
+	      return FALSE;
+	    }
+	}
+      }
+  }
+  return TRUE;
+}
+
+static void fill_speck (FttCell * cell, gboolean * changed)
+{
+  if (GFS_STATE (cell)->div == 1.) {
+    g_free (GFS_STATE (cell)->solid);
+    GFS_STATE (cell)->solid = NULL;
+    *changed = TRUE;
+  }
+}
+
+/**
+ * gfs_domain_remove_specks:
+ * @domain: a #GfsDomain.
+ *
+ * Removes all the solid "specks" of @domain. Solid specks are islands
+ * which do not contain any empty cell.
+ */
+void gfs_domain_remove_specks (GfsDomain * domain)
+{
+  gboolean changed = FALSE;
+
+  g_return_if_fail (domain != NULL);
+
+  gfs_domain_traverse_mixed (domain, FTT_PRE_ORDER, FTT_TRAVERSE_ALL,
+			     (FttCellTraverseFunc) gfs_cell_reset, gfs_div);
+  gfs_domain_traverse_mixed (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS,
+			     (FttCellTraverseFunc) tag_speck, NULL);
+  gfs_domain_traverse_mixed (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS,
+			     (FttCellTraverseFunc) fill_speck, &changed);
+  if (changed)
+    gfs_domain_cell_traverse (domain, FTT_POST_ORDER, FTT_TRAVERSE_NON_LEAFS, -1,
+			      (FttCellTraverseFunc) gfs_cell_init_solid_fractions_from_children, 
+			      NULL);
+}
+
 /**
  * gfs_domain_timer_start:
  * @domain: a #GfsDomain.
diff --git a/src/domain.h b/src/domain.h
index 99465d3..388cafd 100644
--- a/src/domain.h
+++ b/src/domain.h
@@ -208,6 +208,7 @@ void         gfs_domain_remove_ponds          (GfsDomain * domain,
 					       gint min,
 					       FttCellCleanupFunc cleanup,
 					       gpointer data);
+void         gfs_domain_remove_specks         (GfsDomain * domain);
 void         gfs_domain_timer_start           (GfsDomain * domain, 
 					       const gchar * name);
 void         gfs_domain_timer_stop            (GfsDomain * domain, 

-- 
Gerris Flow Solver