diff --git a/Examples/FlowAroundSphere.py b/Examples/FlowAroundSphere.py
index e769394f1f3566c48883faac1f90db10ef77b70c..7ff6bf3f2dc1da533b0a63562458edd927f74f5e 100644
--- a/Examples/FlowAroundSphere.py
+++ b/Examples/FlowAroundSphere.py
@@ -49,11 +49,15 @@ VISCOSITY = 1. / 300.
# ======= Domain =======
dim = 3
-Nx = 257
-Ny = Nz = 129
+Nx = 513
+Ny = Nz = 257
+#Nx = 129
+#Ny = Nz = 129
g = 2
boxlength = npw.asrealarray([10.24, 5.12, 5.12])
boxorigin = npw.asrealarray([-2.0, -2.56, -2.56])
+#boxlength = npw.asrealarray([5.12, 5.12, 5.12])
+#boxorigin = npw.asrealarray([-2.56, -2.56, -2.56])
box = Box(length=boxlength, origin=boxorigin)
# A global discretization with ghost points
@@ -148,10 +152,6 @@ op['diffusion'] = Diffusion(viscosity=VISCOSITY, vorticity=vorti,
discretization=d3d)
rate = VariableParameter(formula=computeFlowrate)
-op['vort_absorption'] = AbsorptionBC(velo, vorti, discretization=d3d,
- req_flowrate=rate,
- x_coords_absorp=[7.24, 8.24])
-
op['poisson'] = Poisson(velo, vorti, discretization=d3d, flowrate=rate)
# ===== Discretization of computational operators ======
@@ -161,6 +161,13 @@ for ope in op.values():
topofft = op['poisson'].discreteFields[vorti].topology
topoadvec = op['advection'].discreteFields[vorti].topology
+# ===== Smooth vorticity absorption at the outlet =====
+op['vort_absorption'] = AbsorptionBC(velo, vorti, discretization=topofft,
+ req_flowrate=rate,
+ x_coords_absorp=[7.24, 8.24])
+# x_coords_absorp=[1.56, 2.56])
+op['vort_absorption'].discretize()
+
# ===== Penalization of the vorticity on a sphere inside the domain =====
from hysop.operator.penalize_vorticity import PenalizeVorticity
op['penalVort'] = PenalizeVorticity(velocity=velo, vorticity=vorti,
@@ -173,6 +180,9 @@ op['penalVort'].discretize()
# ==== Operators to map data between the different computational operators ===
# (i.e. between topologies)
distr = {}
+distr['fft2str'] = RedistributeIntra(source=op['poisson'],
+ target=op['stretching'],
+ variables=[velo, vorti])
distr['fft2str'] = RedistributeIntra(source=op['poisson'],
target=op['stretching'],
variables=[velo, vorti])
@@ -231,9 +241,11 @@ monitors['forcesMom'] = MomentumForces(velocity=velo,
# io_params=io_forcesPenal)
step_dir = ref_step[0]
-io_sliceXY = IO_params('sliceXY', frequency=10)
+io_sliceXY = IO_params('sliceXY', frequency=300)
thickSliceXY = ControlBox(parent=box, origin=[-2.0, -2.56, -2.0 * step_dir],
- length=[10.24, 5.12, 4.0 * step_dir])
+ length=[10.24- step_dir, 5.12- step_dir, 4.0 * step_dir])
+#thickSliceXY = ControlBox(parent=box, origin=[-2.56, -2.56, -2.0 * step_dir],
+# length=[5.12 - step_dir, 5.12 - step_dir, 4.0 * step_dir])
monitors['writerSliceXY'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
io_params=io_sliceXY, subset=thickSliceXY,
xmfalways=True)
@@ -290,11 +302,11 @@ def run(sequence):
distr['fft2str'].apply(simu)
distr['fft2str'].wait()
op['penalVort'].apply(simu) # Vorticity penalization
- distr['str2fft'].apply(simu)
- distr['str2fft'].wait()
- op['poisson'].apply(simu)
- distr['fft2str'].apply(simu)
- distr['fft2str'].wait()
+# distr['str2fft'].apply(simu)
+# distr['str2fft'].wait()
+# op['poisson'].apply(simu)
+# distr['fft2str'].apply(simu)
+# distr['fft2str'].wait()
op['stretching'].apply(simu) # Stretching
# monitors['forcesNoca'].apply(simu) # Forces Noca
distr['str2fft'].apply(simu)
diff --git a/HySoP/hysop/operator/discrete/absorption_BC.py b/HySoP/hysop/operator/discrete/absorption_BC.py
index 5ca65422f3f9d0e9e146ddb748605318b906e43a..33223c43e993c4c6e680d932eeb7222fe4fa35fc 100755
--- a/HySoP/hysop/operator/discrete/absorption_BC.py
+++ b/HySoP/hysop/operator/discrete/absorption_BC.py
@@ -20,8 +20,8 @@ class AbsorptionBC_D(DiscreteOperator):
"""
The periodic boundary condition is modified at the outlet
in the flow direction in order to discard
- in the dowstream region the eddies coming
- periodically from the oulet.
+ in the downstream region the eddies coming
+ periodically from the outlet.
The vorticity absorption conserves div(omega)=0.
The far field velocity is set to u_inf at the inlet.
"""
@@ -80,7 +80,7 @@ class AbsorptionBC_D(DiscreteOperator):
self.xb = self.x_coords_absorp[0]
self.xe = self.x_coords_absorp[1]
self.xc = self.xb + (self.xe - self.xb) / 2.0
- self.eps = 10
+ self.eps = 10.0
self.coeff = 1.0 / (np.tanh(self.eps * (self.xb - self.xc)) -
np.tanh(self.eps * (self.xe - self.xc)))
self.coords = self.topo.mesh.coords
@@ -98,18 +98,18 @@ class AbsorptionBC_D(DiscreteOperator):
self._filter = npw.ones_like(self.vorticity.data[0])
indFilter = np.where(np.logical_and(self.coords[0][:,0,0] >= self.xb,
self.coords[0][:,0,0] <= self.xe))
- indFilterZero = np.where(self.coords[0][:,0,0] > self.xe)
+# indFilterZero = np.where(self.coords[0][:,0,0] > self.xe)
FiltFormula = np.tanh(self.eps * (self.coords[0][indFilter] -
self.xc))
FiltFormula -= np.tanh(self.eps * (self.xe - self.xc))
FiltFormula *= self.coeff
- self._filter[indFilter] = FiltFormula
- self._filter[indFilterZero] = 0.0
+ self._filter[indFilter,:,:] = FiltFormula
+# self._filter[indFilterZero] = 0.0
# Beginning of divergence free vorticity absorption
# for non-periodic inlet BC
for d in xrange(self.vorticity.nbComponents):
- self.vorticity[d][...] *= self._filter[...]
+ self.vorticity.data[d][...] *= self._filter[...]
# Definition of the X-derivative of the filter function
self._filter = npw.zeros_like(self.vorticity.data[0])
@@ -122,11 +122,11 @@ class AbsorptionBC_D(DiscreteOperator):
self._filter[indFilter] = FiltFormula
# End of divergence free vorticity absorption for non-periodic inlet BC
- self.vorticity[YDIR][...] += (- self._filter[...] *
- self.velocity[ZDIR][...]) + \
- (self._filter[...] *
- self.req_flowrate_val[ZDIR])
- self.vorticity[ZDIR][...] += (self._filter[...] *
- self.velocity[YDIR][...]) - \
- (self._filter[...] *
- self.req_flowrate_val[YDIR])
+ self.vorticity.data[YDIR][...] += (- self._filter[...] *
+ self.velocity[ZDIR][...]) + \
+ (self._filter[...] *
+ self.req_flowrate_val[ZDIR])
+ self.vorticity.data[ZDIR][...] += (self._filter[...] *
+ self.velocity[YDIR][...]) - \
+ (self._filter[...] *
+ self.req_flowrate_val[YDIR])