[SCM] Gerris Flow Solver branch, upstream, updated. b3aa46814a06c9cb2912790b23916ffb44f1f203
Stephane Popinet
popinet at users.sourceforge.net
Fri May 15 02:51:20 UTC 2009
The following commit has been merged in the upstream branch:
commit 5da804f5d3303fc4f6582edc68a7f3fffbf04841
Author: Stephane Popinet <popinet at users.sourceforge.net>
Date: Tue Nov 2 11:49:45 2004 +1100
Fractions initialisation uses new solid fraction algorithm (gerris--mainline--0.7--patch-20)
gerris--mainline--0.7--patch-20
Keywords:
The gfs_cell_init_fraction() function has been replaced by
gfs_domain_init_fraction() which now calls
gfs_domain_init_solid_fractions().
The GfsInitFraction event has been changed accordingly.
darcs-hash:20041102004945-aabb8-44d6609bfeed3ff82382f86474af2c77f446fb78.gz
diff --git a/doc/gfs-sections.txt b/doc/gfs-sections.txt
index e2e7c8a..1d7173b 100644
--- a/doc/gfs-sections.txt
+++ b/doc/gfs-sections.txt
@@ -167,7 +167,6 @@ gfs_get_from_below_extensive
gfs_get_from_above
gfs_cell_coarse_init
gfs_cell_fine_init
-gfs_cell_init_fraction
gfs_velocity_norm
gfs_velocity_norm2
gfs_divergence
@@ -319,6 +318,7 @@ gfs_domain_norm_velocity
gfs_domain_stats_balance
gfs_domain_stats_solid
gfs_domain_solid_force
+gfs_domain_init_fraction
gfs_domain_read
gfs_domain_split
gfs_domain_variable_sources
diff --git a/doc/tmpl/init_fraction.sgml b/doc/tmpl/init_fraction.sgml
index 9490487..60fd6dc 100644
--- a/doc/tmpl/init_fraction.sgml
+++ b/doc/tmpl/init_fraction.sgml
@@ -15,7 +15,7 @@ The syntax in parameter files is as follows:
[ #GfsGenericInit ] C surface.gts
</programlisting>
</informalexample>
-where C is a variable name and surface.gts is a closed, manifold surface.
+where C is a variable name and surface.gts is an orientable surface.
</para>
<!-- ##### SECTION See_Also ##### -->
diff --git a/src/event.c b/src/event.c
index fcf87ab..c3aea2d 100644
--- a/src/event.c
+++ b/src/event.c
@@ -1087,8 +1087,6 @@ static void gfs_init_fraction_destroy (GtsObject * object)
if (init->surface)
gts_object_destroy (GTS_OBJECT (init->surface));
- if (init->stree)
- gts_bb_tree_destroy (init->stree, TRUE);
(* GTS_OBJECT_CLASS (gfs_init_fraction_class ())->parent_class->destroy)
(object);
@@ -1098,7 +1096,6 @@ static void gfs_init_fraction_read (GtsObject ** o, GtsFile * fp)
{
GfsInitFraction * init;
GfsDomain * domain;
- GtsSurface * self;
if (GTS_OBJECT_CLASS (gfs_init_fraction_class ())->parent_class->read)
(* GTS_OBJECT_CLASS (gfs_init_fraction_class ())->parent_class->read)
@@ -1162,20 +1159,7 @@ static void gfs_init_fraction_read (GtsObject ** o, GtsFile * fp)
gts_file_error (fp, "surface is not orientable");
return;
}
- if (!gts_surface_is_closed (init->surface)) {
- gts_file_error (fp, "surface is not closed");
- return;
- }
- if ((self = gts_surface_is_self_intersecting (init->surface))) {
- gts_object_destroy (GTS_OBJECT (self));
- gts_file_error (fp, "surface is self-intersecting");
- return;
- }
- init->stree = gts_bb_tree_surface (init->surface);
- if (gts_surface_volume (init->surface) < 0.)
- init->is_open = TRUE;
-
gts_file_next_token (fp);
}
@@ -1191,18 +1175,13 @@ static void gfs_init_fraction_write (GtsObject * o, FILE * fp)
fputs ("}\n", fp);
}
-static void box_init_fraction (GfsBox * box, GfsInitFraction * init)
-{
- gfs_cell_init_fraction (box->root,
- init->surface, init->stree, init->is_open,
- init->c);
-}
-
static gboolean gfs_init_fraction_event (GfsEvent * event, GfsSimulation * sim)
{
- if ((* GFS_EVENT_CLASS (GTS_OBJECT_CLASS (gfs_init_fraction_class ())->parent_class)->event) (event, sim)) {
- gts_container_foreach (GTS_CONTAINER (sim),
- (GtsFunc) box_init_fraction, event);
+ if ((* GFS_EVENT_CLASS (GTS_OBJECT_CLASS (gfs_init_fraction_class ())->parent_class)->event)
+ (event, sim)) {
+ gfs_domain_init_fraction (GFS_DOMAIN (sim),
+ GFS_INIT_FRACTION (event)->surface,
+ GFS_INIT_FRACTION (event)->c);
return TRUE;
}
return FALSE;
diff --git a/src/event.h b/src/event.h
index d48e0b6..690a208 100644
--- a/src/event.h
+++ b/src/event.h
@@ -199,8 +199,6 @@ struct _GfsInitFraction {
GfsVariable * c;
GtsSurface * surface;
- GNode * stree;
- gboolean is_open;
};
typedef struct _GfsInitFractionClass GfsInitFractionClass;
diff --git a/src/solid.c b/src/solid.c
index 4ace14f..a6a25ad 100644
--- a/src/solid.c
+++ b/src/solid.c
@@ -224,7 +224,7 @@ static void set_solid_fractions_from_surface (FttCell * cell,
}
default:
g_log (G_LOG_DOMAIN, G_LOG_LEVEL_ERROR,
- "the solid surface is probably not closed (n1 = %d)", n1);
+ "the surface is probably not closed (n1 = %d)", n1);
}
}
#else /* 3D */
@@ -374,7 +374,7 @@ static void set_solid_fractions_from_surface (FttCell * cell, GtsSurface * s)
}
default:
g_log (G_LOG_DOMAIN, G_LOG_LEVEL_ERROR,
- "the solid surface is probably not closed (n2 = %d)", n2);
+ "the surface is probably not closed (n2 = %d)", n2);
}
}
@@ -603,10 +603,18 @@ static void solid_fractions_from_children (FttCell * cell, InitSolidParams * p)
GFS_STATE (cell)->div = 1.;
else {
gfs_cell_init_solid_fractions_from_children (cell);
- GFS_STATE (cell)->div = 0.;
+ if (p->destroy_solid)
+ GFS_STATE (cell)->div = 0.;
+ else if (!GFS_IS_MIXED (cell)) {
+ ftt_cell_children (cell, &child);
+ GFS_STATE (cell)->div = 1.;
+ for (i = 0; i < FTT_CELLS; i++)
+ if (child.c[i] && GFS_STATE (child.c[i])->div == 2.)
+ GFS_STATE (cell)->div = 2.;
+ }
}
}
- if (GFS_STATE (cell)->div == 1.) {
+ if (p->destroy_solid && GFS_STATE (cell)->div == 1.) {
if (FTT_CELL_IS_ROOT (cell))
g_log (G_LOG_DOMAIN, G_LOG_LEVEL_ERROR,
"root cell is entirely outside of the fluid domain\n"
@@ -812,294 +820,62 @@ gboolean gfs_refine_mixed (const FttCell * cell)
return FALSE;
}
-/* GfsFace: Header */
-
-typedef struct _GfsFace GfsFace;
-typedef struct _GfsFaceClass GfsFaceClass;
-
-struct _GfsFace {
- GtsFace parent;
-
- FttDirection dir;
-};
-
-struct _GfsFaceClass {
- GtsFaceClass parent_class;
-};
-
-#define GFS_FACE(obj) GTS_OBJECT_CAST (obj,\
- GfsFace,\
- gfs_face_class ())
-#define GFS_FACE_CLASS(klass) GTS_OBJECT_CLASS_CAST (klass,\
- GfsFaceClass,\
- gfs_face_class())
-#define GFS_IS_FACE(obj) (gts_object_is_from_class (obj,\
- gfs_face_class ()))
-
-static GfsFaceClass * gfs_face_class (void);
-static GfsFace * gfs_face_new (GfsFaceClass * klass,
- GtsEdge * e1,
- GtsEdge * e2,
- GtsEdge * e3,
- guint dir);
-
-/* GfsFace: Object */
-
-static void gfs_face_link (GtsObject * object, GtsObject * with)
+static void save_solid (FttCell * cell, GfsVariable * c)
{
- GFS_FACE (object)->dir = GFS_FACE (with)->dir;
+ GFS_DOUBLE_TO_POINTER (GFS_VARIABLE (cell, c->i)) = GFS_STATE (cell)->solid;
+ GFS_STATE (cell)->solid = NULL;
}
-static void gfs_face_class_init (GfsFaceClass * klass)
+static void restore_solid (FttCell * cell, gpointer * data)
{
- GTS_OBJECT_CLASS (klass)->attributes = gfs_face_link;
-}
-
-static void gfs_face_init (GfsFace * f)
-{
- f->dir = 0;
-}
+ GfsVariable * c = data[0];
+ gboolean * not_cut = data[1];
+ GfsSolidVector * solid = GFS_STATE (cell)->solid;
-static GfsFaceClass * gfs_face_class (void)
-{
- static GfsFaceClass * klass = NULL;
-
- if (klass == NULL) {
- GtsObjectClassInfo gfs_face_info = {
- "GfsFace",
- sizeof (GfsFace),
- sizeof (GfsFaceClass),
- (GtsObjectClassInitFunc) gfs_face_class_init,
- (GtsObjectInitFunc) gfs_face_init,
- (GtsArgSetFunc) NULL,
- (GtsArgGetFunc) NULL
- };
- klass = gts_object_class_new (GTS_OBJECT_CLASS (gts_face_class ()),
- &gfs_face_info);
+ GFS_STATE (cell)->solid = GFS_DOUBLE_TO_POINTER (GFS_VARIABLE (cell, c->i));
+ if (solid) {
+ GFS_VARIABLE (cell, c->i) = solid->a;
+ g_free (solid);
+ *not_cut = FALSE;
}
-
- return klass;
-}
-
-static GfsFace * gfs_face_new (GfsFaceClass * klass,
- GtsEdge * e1,
- GtsEdge * e2,
- GtsEdge * e3,
- guint dir)
-{
- GfsFace * f = GFS_FACE (gts_face_new (GTS_FACE_CLASS (klass), e1, e2, e3));
-
- f->dir = dir;
- return f;
-}
-
-static void surface_add_box (GtsSurface * s,
- gdouble x1, gdouble y1, gdouble z1,
- gdouble x2, gdouble y2, gdouble z2)
-{
- GtsVertex * v0 = gts_vertex_new (s->vertex_class, x1, y1, z1);
- GtsVertex * v1 = gts_vertex_new (s->vertex_class, x1, y1, z2);
- GtsVertex * v2 = gts_vertex_new (s->vertex_class, x1, y2, z2);
- GtsVertex * v3 = gts_vertex_new (s->vertex_class, x1, y2, z1);
- GtsVertex * v4 = gts_vertex_new (s->vertex_class, x2, y1, z1);
- GtsVertex * v5 = gts_vertex_new (s->vertex_class, x2, y1, z2);
- GtsVertex * v6 = gts_vertex_new (s->vertex_class, x2, y2, z2);
- GtsVertex * v7 = gts_vertex_new (s->vertex_class, x2, y2, z1);
-
- GtsEdge * e1 = gts_edge_new (s->edge_class, v0, v1);
- GtsEdge * e2 = gts_edge_new (s->edge_class, v1, v2);
- GtsEdge * e3 = gts_edge_new (s->edge_class, v2, v3);
- GtsEdge * e4 = gts_edge_new (s->edge_class, v3, v0);
- GtsEdge * e5 = gts_edge_new (s->edge_class, v0, v2);
-
- GtsEdge * e6 = gts_edge_new (s->edge_class, v4, v5);
- GtsEdge * e7 = gts_edge_new (s->edge_class, v5, v6);
- GtsEdge * e8 = gts_edge_new (s->edge_class, v6, v7);
- GtsEdge * e9 = gts_edge_new (s->edge_class, v7, v4);
- GtsEdge * e10 = gts_edge_new (s->edge_class, v4, v6);
-
- GtsEdge * e11 = gts_edge_new (s->edge_class, v3, v7);
- GtsEdge * e12 = gts_edge_new (s->edge_class, v2, v6);
- GtsEdge * e13 = gts_edge_new (s->edge_class, v1, v5);
- GtsEdge * e14 = gts_edge_new (s->edge_class, v0, v4);
-
- GtsEdge * e15 = gts_edge_new (s->edge_class, v1, v6);
- GtsEdge * e16 = gts_edge_new (s->edge_class, v2, v7);
- GtsEdge * e17 = gts_edge_new (s->edge_class, v3, v4);
- GtsEdge * e18 = gts_edge_new (s->edge_class, v0, v5);
-
- GfsFaceClass * klass = gfs_face_class ();
-
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e1, e2, e5, 1)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e5, e3, e4, 1)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e6, e10, e7, 0)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e10, e9, e8, 0)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e2, e15, e12, 4)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e15, e13, e7, 4)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e3, e16, e11, 2)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e16, e12, e8, 2)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e17, e14, e4, 5)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e17, e11, e9, 5)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e18, e13, e1, 3)));
- gts_surface_add_face (s, GTS_FACE (gfs_face_new (klass, e18, e14, e6, 3)));
-}
-
-static void cell_size (FttCell * cell, GtsVector s)
-{
-#if FTT_2D3
- s[0] = s[1] = ftt_cell_size (cell);
- s[2] = 1.;
-#else /* 2D or 3D */
- s[0] = s[1] = s[2] = ftt_cell_size (cell);
-#endif /* 2D or 3D */
-}
-
-static GtsBBox * bbox_cell (GtsBBoxClass * klass, FttCell * cell)
-{
- FttVector p;
- GtsVector size;
-
- ftt_cell_pos (cell, &p);
- cell_size (cell, size);
- return gts_bbox_new (klass, cell,
- p.x - size[0]/2., p.y - size[1]/2., p.z - size[2]/2.,
- p.x + size[0]/2., p.y + size[1]/2., p.z + size[2]/2.);
-}
-
-static void gfs_cell_set_fraction (FttCell * cell,
- GfsVariable * c,
- gdouble val)
-{
- g_return_if_fail (cell != NULL);
-
- GFS_VARIABLE (cell, c->i) = val;
- if (!FTT_CELL_IS_LEAF (cell)) {
- FttCellChildren child;
- guint i;
-
- ftt_cell_children (cell, &child);
- for (i = 0; i < FTT_CELLS; i++)
- if (child.c[i])
- gfs_cell_set_fraction (child.c[i], c, val);
+ else if (GFS_STATE (cell)->div == 0.) {
+ g_assert (*not_cut);
+ GFS_VARIABLE (cell, c->i) = 0.;
}
-}
-
-static void set_full_or_empty (FttCell * cell,
- GNode * tree,
- gboolean is_open,
- GfsVariable * c)
-{
- FttVector pos;
- GtsPoint * p;
-
- ftt_cell_pos (cell, &pos);
- p = gts_point_new (gts_point_class (), pos.x, pos.y, pos.z);
- if (gts_point_is_inside_surface (p, tree, is_open))
- gfs_cell_set_fraction (cell, c, 1.);
- else
- gfs_cell_set_fraction (cell, c, 0.);
- gts_object_destroy (GTS_OBJECT (p));
-}
-
-static void bbox_size (GtsBBox * bbox, GtsVector s)
-{
- s[0] = bbox->x2 - bbox->x1;
- s[1] = bbox->y2 - bbox->y1;
- s[2] = bbox->z2 - bbox->z1;
-}
-
-static void set_fraction_from_surface (FttCell * cell,
- GtsBBox * bbox,
- GtsSurface * s,
- GNode * stree,
- gboolean is_open,
- GfsVariable * c)
-{
- GtsSurface * cs;
- GNode * ctree;
- GtsSurfaceInter * si;
- gboolean closed = TRUE;
-
- cs = gts_surface_new (gts_surface_class (),
- GTS_FACE_CLASS (gfs_face_class ()),
- gts_edge_class (),
- gts_vertex_class ());
- surface_add_box (cs,
- bbox->x1, bbox->y1, bbox->z1,
- bbox->x2, bbox->y2, bbox->z2);
- ctree = gts_bb_tree_surface (cs);
- si = gts_surface_inter_new (gts_surface_inter_class (),
- cs, s, ctree, stree, FALSE, is_open);
- g_assert (gts_surface_inter_check (si, &closed));
- if (si->edges == NULL)
- set_full_or_empty (cell, stree, is_open, c);
else {
- GtsSurface * inter = gts_surface_new (gts_surface_class (),
- gts_face_class (),
- gts_edge_class (),
- gts_vertex_class ());
- GtsVector size;
-
- g_assert (closed);
- gts_surface_inter_boolean (si, inter, GTS_1_IN_2);
- gts_surface_inter_boolean (si, inter, GTS_2_IN_1);
- bbox_size (bbox, size);
- GFS_VARIABLE (cell, c->i) = gts_surface_volume (inter)/(size[0]*size[1]*size[2]);
- g_assert (GFS_VARIABLE (cell, c->i) > -1e-9 &&
- GFS_VARIABLE (cell, c->i) < 1. + 1e-9);
- gts_object_destroy (GTS_OBJECT (inter));
+ g_assert (GFS_STATE (cell)->div == 1. || GFS_STATE (cell)->div == 2.);
+ GFS_VARIABLE (cell, c->i) = GFS_STATE (cell)->div - 1.;
}
-
- gts_object_destroy (GTS_OBJECT (si));
- gts_bb_tree_destroy (ctree, TRUE);
- gts_object_destroy (GTS_OBJECT (cs));
}
/**
- * gfs_cell_init_fraction:
- * @root: the root #FttCell of the cell tree.
- * @s: a closed, orientable surface defining the solid boundary.
- * @stree: a bounding box tree of the faces of @s.
- * @is_open: %TRUE if @s is an "open" boundary i.e. the signed volume
- * enclosed by @s is negative, %FALSE otherwise.
+ * gfs_domain_init_fraction:
+ * @domain: a #GfsDomain.
+ * @s: an orientable surface defining the interface boundary.
* @c: a #GfsVariable.
*
- * Initializes the fraction @c of all the cells of the cell tree
- * starting at @root.
+ * Initializes the fraction @c of the interface @s contained in all
+ * the cells of @domain.
*/
-void gfs_cell_init_fraction (FttCell * root,
- GtsSurface * s,
- GNode * stree,
- gboolean is_open,
- GfsVariable * c)
+void gfs_domain_init_fraction (GfsDomain * domain,
+ GtsSurface * s,
+ GfsVariable * c)
{
- GtsBBox * bbox;
+ gboolean not_cut = TRUE;
+ gpointer data[2];
- g_return_if_fail (root != NULL);
+ g_return_if_fail (domain != NULL);
g_return_if_fail (s != NULL);
- g_return_if_fail (stree != NULL);
g_return_if_fail (c != NULL);
- bbox = bbox_cell (gts_bbox_class (), root);
- if (gts_bb_tree_is_overlapping (stree, bbox)) {
- if (FTT_CELL_IS_LEAF (root))
- set_fraction_from_surface (root, bbox, s, stree, is_open, c);
- else {
- FttCellChildren child;
- guint i;
-
- ftt_cell_children (root, &child);
- for (i = 0; i < FTT_CELLS; i++)
- if (child.c[i])
- gfs_cell_init_fraction (child.c[i], s, stree, is_open, c);
- gfs_get_from_below_extensive (root, c);
- GFS_VARIABLE (root, c->i) /= 4.;
- }
- }
- else /* !gts_bb_tree_is_overlapping (stree, bbox) */
- set_full_or_empty (root, stree, is_open, c);
-
- gts_object_destroy (GTS_OBJECT (bbox));
+ gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_ALL, -1,
+ (FttCellTraverseFunc) save_solid, c);
+ gfs_domain_init_solid_fractions (domain, s, FALSE, NULL, NULL);
+ data[0] = c;
+ data[1] = ¬_cut;
+ gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_ALL, -1,
+ (FttCellTraverseFunc) restore_solid, data);
+ gfs_domain_bc (domain, FTT_TRAVERSE_LEAFS, -1, c);
}
/**
diff --git a/src/solid.h b/src/solid.h
index 51b0492..8907a35 100644
--- a/src/solid.h
+++ b/src/solid.h
@@ -37,10 +37,8 @@ void gfs_domain_init_solid_fractions (GfsDomain * domain,
void gfs_cell_init_solid_fractions_from_children (FttCell * cell);
gboolean gfs_cell_check_solid_fractions (FttCell * root);
gboolean gfs_refine_mixed (const FttCell * cell);
-void gfs_cell_init_fraction (FttCell * root,
+void gfs_domain_init_fraction (GfsDomain * domain,
GtsSurface * s,
- GNode * stree,
- gboolean is_open,
GfsVariable * c);
void gfs_cell_cm (const FttCell * cell,
FttVector * cm);
--
Gerris Flow Solver
More information about the debian-science-commits
mailing list