sandbox/popinet/steel.c
How to impose boundary conditions on a masked box
#include "grid/octree.h"
#include "navier-stokes/centered.h"
// NOTE: this uses global coordinates for u, i.e. for any condition below
// u.n -> u.x
// u.t -> u.y
// u.r -> u.z
// top is an outflow
u.t[top] = neumann(0);
p[top] = dirichlet(0.);
pf[top] = dirichlet(0.);
u.n[top] = neumann(0);
u.r[top] = neumann(0);
// bottom is a slip wall
u.t[bottom] = dirichlet(0);
u.n[bottom] = neumann(0);
u.r[bottom] = neumann(0);
// right and left are outflows
u.n[left] = neumann(0);
p[left] = dirichlet(0.);
pf[left] = dirichlet(0.);
u.t[left] = neumann(0.);
u.r[left] = neumann(0.);
u.n[right] = neumann(0.);
p[right] = dirichlet(0.);
pf[right] = dirichlet(0.);
u.t[right] = neumann(0.);
u.r[right] = neumann(0.);
// front and back are slip walls
u.r[front] = dirichlet(0.);
u.t[front] = neumann(0.);
u.n[front] = neumann(0.);
u.r[back] = dirichlet(0.);
u.t[back] = neumann(0.);
u.n[back] = neumann(0.);
// the steel band u.y = 1, u.x = u.z = 0
bid steel;
u.n[steel] = dirichlet(0.);
u.r[steel] = dirichlet(0.);
u.t[steel] = dirichlet(1.);
int main()
{
init_grid (16);
origin (-L0/2.,-L0/2.,-L0/2.);
DT = 0.01;
mu = unityf;
run ();
}
event init (i = 0) {
mask (fabs(x) < 0.25 && fabs(z) < 0.25 ? steel : none);
foreach()
u.x[] = u.y[] = u.z[] = -1;
for (scalar s in {u})
foreach_dimension()
s.v.x.i = -1; // u is now not a vector
boundary ((scalar *){u});
foreach_boundary (steel)
printf ("steel %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
foreach_boundary (left)
printf ("left %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
foreach_boundary (right)
printf ("right %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
foreach_boundary (top)
printf ("top %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
foreach_boundary (bottom)
printf ("bottom %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
foreach_boundary (front)
printf ("front %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
foreach_boundary (back)
printf ("back %g %g %g %g %g %g\n", x, y, z,
(u.x[] + u.x[ghost])/2.,
(u.y[] + u.y[ghost])/2.,
(u.z[] + u.z[ghost])/2.);
}
The values on the boundaries look OK.
set xlabel 'x'
set xlabel 'y'
set xlabel 'z'
scale = 0.1
set view 120, 326, 1, 1
splot '< grep steel out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'steel', '< grep left out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'left', '< grep right out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'right', '< grep top out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'top', '< grep bottom out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'bottom', '< grep front out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'front', '< grep back out' u 2:3:4:($5*scale):($6*scale):($7*scale) w vector t 'back'