[SCM] Gerris Flow Solver branch, upstream, updated. b3aa46814a06c9cb2912790b23916ffb44f1f203
Stephane Popinet
s.popinet at niwa.co.nz
Fri May 15 02:51:36 UTC 2009
The following commit has been merged in the upstream branch:
commit 88364058e1eff9d6b82162aece5b3e69b45343a0
Author: Stephane Popinet <s.popinet at niwa.co.nz>
Date: Tue Feb 8 11:08:20 2005 +1100
Replaced gfs_youngs_gradient with gfs_youngs_normal
Which computes the three components in one operation and should be about
three times as fast as three calls to gfs_youngs_gradient.
darcs-hash:20050208000820-fbd8f-e24a546ee1e613c53723dd41f42ee73357e4ec3e.gz
diff --git a/src/tension.c b/src/tension.c
index e142dea..8a0736a 100644
--- a/src/tension.c
+++ b/src/tension.c
@@ -75,10 +75,9 @@ static void foreach_cell_normal (FttCell * cell, GfsSourceTension * s)
gdouble sigh = s->sigma/ftt_cell_size (cell);
FttComponent c;
- for (c = 0; c < FTT_DIMENSION; c++) {
- (&n.x)[c] = gfs_youngs_gradient (cell, c, s->c);
+ gfs_youngs_normal (cell, s->c, &n);
+ for (c = 0; c < FTT_DIMENSION; c++)
nn += (&n.x)[c]*(&n.x)[c];
- }
nn = sqrt (nn + 1e-50);
GFS_STATE (cell)->g[0] = sigh*n.x*n.x/nn;
GFS_STATE (cell)->g[1] = sigh*n.y*n.y/nn;
@@ -88,11 +87,14 @@ static void foreach_cell_normal (FttCell * cell, GfsSourceTension * s)
static void foreach_cell_tension (FttCell * cell, GfsSourceTension * s)
{
gdouble h = ftt_cell_size (cell);
+ FttVector nx, ny, nxy;
+
+ gfs_youngs_normal (cell, gfs_gx, &nx);
+ gfs_youngs_normal (cell, gfs_gy, &ny);
+ gfs_youngs_normal (cell, gfs_div, &nxy);
- GFS_VARIABLE (cell, s->t[0]->i) = (gfs_youngs_gradient (cell, FTT_X, gfs_gy) -
- gfs_youngs_gradient (cell, FTT_Y, gfs_div))/h;
- GFS_VARIABLE (cell, s->t[1]->i) = (gfs_youngs_gradient (cell, FTT_Y, gfs_gx) -
- gfs_youngs_gradient (cell, FTT_X, gfs_div))/h;
+ GFS_VARIABLE (cell, s->t[0]->i) = (ny.x - nxy.y)/h;
+ GFS_VARIABLE (cell, s->t[1]->i) = (nx.y - nxy.x)/h;
}
static void gfs_source_tension_event (GfsEvent * event,
diff --git a/src/vof.c b/src/vof.c
index 8a83100..81d6c74 100644
--- a/src/vof.c
+++ b/src/vof.c
@@ -360,47 +360,45 @@ gdouble gfs_plane_alpha (FttVector * m, gdouble c)
}
#endif /* 3D */
-static gdouble avg (FttCell * cell, FttDirection d, GfsVariable * v)
-{
- FttDirection d1[2*(FTT_DIMENSION - 1)][FTT_DIMENSION];
- FttComponent c = FTT_ORTHOGONAL_COMPONENT (d/2);
- gdouble v1 = 0.;
- guint i;
-
-#if FTT_2D
- d1[0][0] = d1[1][0] = d;
- d1[0][1] = 2*c;
- d1[1][1] = 2*c + 1;
-#else /* 3D */
- FttComponent c1 = FTT_ORTHOGONAL_COMPONENT (c);
-
- d1[0][0] = d1[1][0] = d1[2][0] = d1[3][0] = d;
- d1[0][1] = 2*c; d1[0][2] = 2*c1;
- d1[1][1] = 2*c; d1[1][2] = 2*c1 + 1;
- d1[2][1] = 2*c + 1; d1[2][2] = 2*c1;
- d1[3][1] = 2*c + 1; d1[3][2] = 2*c1 + 1;
-#endif /* 3D */
-
- for (i = 0; i < 2*(FTT_DIMENSION - 1); i++)
- v1 += gfs_cell_corner_value (cell, d1[i], v, -1);
- return v1/(2*(FTT_DIMENSION - 1));
-}
-
/**
- * gfs_youngs_gradient:
+ * gfs_youngs_normal:
* @cell: a #FttCell.
- * @c: a component.
* @v: a #GfsVariable.
+ * @n: a #FttVector.
*
- * Returns: the Youngs-averaged gradient of variable @v in direction
- * @c normalised by the size of @cell.
+ * Fills @n with the Youngs-averaged gradients of @v
+ * normalised by the size of @cell.
*/
-gdouble gfs_youngs_gradient (FttCell * cell, FttComponent c, GfsVariable * v)
+void gfs_youngs_normal (FttCell * cell, GfsVariable * v, FttVector * n)
{
- g_return_val_if_fail (cell != NULL, 0.);
- g_return_val_if_fail (v != NULL, 0.);
+ static FttDirection d[(FTT_DIMENSION - 1)*4][FTT_DIMENSION] = {
+#if FTT_2D
+ {FTT_RIGHT, FTT_TOP}, {FTT_LEFT, FTT_TOP}, {FTT_LEFT, FTT_BOTTOM}, {FTT_RIGHT, FTT_BOTTOM}
+#else /* 3D */
+ {FTT_RIGHT, FTT_TOP, FTT_FRONT}, {FTT_LEFT, FTT_TOP, FTT_FRONT},
+ {FTT_LEFT, FTT_BOTTOM, FTT_FRONT}, {FTT_RIGHT, FTT_BOTTOM, FTT_FRONT},
+ {FTT_RIGHT, FTT_TOP, FTT_BACK}, {FTT_LEFT, FTT_TOP, FTT_BACK},
+ {FTT_LEFT, FTT_BOTTOM, FTT_BACK}, {FTT_RIGHT, FTT_BOTTOM, FTT_BACK},
+#endif /* 3D */
+ };
+ gdouble u[(FTT_DIMENSION - 1)*4];
+ guint i;
+
+ g_return_if_fail (cell != NULL);
+ g_return_if_fail (v != NULL);
+ g_return_if_fail (n != NULL);
+
+ for (i = 0; i < (FTT_DIMENSION - 1)*4; i++)
+ u[i] = gfs_cell_corner_value (cell, d[i], v, -1);
- return avg (cell, 2*c, v) - avg (cell, 2*c + 1, v);
+#if FTT_2D
+ n->x = (u[0] + u[3] - u[1] - u[2])/2.;
+ n->y = (u[0] + u[1] - u[2] - u[3])/2.;
+#else /* 3D */
+ n->x = (u[0] + u[3] + u[4] + u[7] - u[1] - u[2] - u[5] - u[6])/4.;
+ n->y = (u[0] + u[1] + u[4] + u[5] - u[2] - u[3] - u[6] - u[7])/4.;
+ n->z = (u[0] + u[1] + u[2] + u[3] - u[4] - u[5] - u[6] - u[7])/4.;
+#endif /* 3D */
}
static void gfs_cell_vof_advected_face_values (FttCell * cell,
@@ -433,7 +431,7 @@ static void gfs_cell_vof_advected_face_values (FttCell * cell,
GFS_STATE (cell)->f[left].v = f;
}
else {
- FttVector m;
+ FttVector m1, m;
gdouble alpha, n;
static FttComponent d[FTT_DIMENSION][FTT_DIMENSION - 1] = {
#if FTT_2D
@@ -444,9 +442,10 @@ static void gfs_cell_vof_advected_face_values (FttCell * cell,
};
guint i;
- m.x = - gfs_youngs_gradient (cell, c, par->v);
+ gfs_youngs_normal (cell, par->v, &m1);
+ m.x = - (&m1.x)[c];
for (i = 0; i < FTT_DIMENSION - 1; i++)
- (&m.x)[i + 1] = - gfs_youngs_gradient (cell, d[c][i], par->v);
+ (&m.x)[i + 1] = - (&m1.x)[d[c][i]];
if (m.x < 0.) {
n = - u_left; u_left = - u_right; u_right = n;
diff --git a/src/vof.h b/src/vof.h
index 5d5d71e..041314e 100644
--- a/src/vof.h
+++ b/src/vof.h
@@ -49,9 +49,9 @@ void gfs_plane_center (FttVector * m,
gdouble a,
FttVector * p);
#endif /* 3D */
-gdouble gfs_youngs_gradient (FttCell * cell,
- FttComponent c,
- GfsVariable * v);
+void gfs_youngs_normal (FttCell * cell,
+ GfsVariable * v,
+ FttVector * n);
void gfs_cell_vof_advection (FttCell * cell,
FttComponent c,
GfsAdvectionParams * par);
--
Gerris Flow Solver
More information about the debian-science-commits
mailing list