Skip to content
Snippets Groups Projects
Commit d8bd1e28 authored by Chloé Mimeau's avatar Chloé Mimeau
Browse files

update sphere Example

parent 16e1ee27
No related branches found
No related tags found
No related merge requests found
......@@ -31,6 +31,7 @@ from hysop.numerics.finite_differences import FD_C_4, FD_C_2
from hysop.numerics.interpolation import Linear
from hysop.numerics.remeshing import L6_4 as rmsh
import hysop.tools.io_utils as io
import hysop.tools.numpywrappers as npw
from hysop.mpi import main_rank, MPI
from hysop.tools.parameters import Discretization, IO_params
......@@ -50,9 +51,9 @@ dim = 3
Nx = 129
Ny = Nz = 65
g = 2
boxlength = npw.realarray([10.24, 5.12, 5.12])
boxorigin = npw.realarray([-2.0, -2.56, -2.56])
box = Box(dim, length=boxlength, origin=boxorigin)
boxlength = npw.asrealarray([10.24, 5.12, 5.12])
boxorigin = npw.asrealarray([-2.0, -2.56, -2.56])
box = Box(length=boxlength, origin=boxorigin)
# A global discretization with ghost points
d3dg = Discretization([Nx, Ny, Nz], [g, g, g])
......@@ -81,7 +82,6 @@ def computeVort(res, x, y, z, t):
res[2][...] = 0.
return res
# ====== Time-dependant required-flowrate (Variable Parameter) ======
def computeFlowrate(simu):
# === Time-dependant flow rate ===
......@@ -149,6 +149,7 @@ op['diffusion'] = Diffusion(viscosity=VISCOSITY, vorticity=vorti,
rate = VariableParameter(formula=computeFlowrate)
op['poisson'] = Poisson(velo, vorti, discretization=d3d, flowrate=rate)
#op['poisson'] = Poisson(velo, vorti, discretization=d3d)
# ===== Discretization of computational operators ======
for ope in op.values():
......@@ -158,8 +159,10 @@ topofft = op['poisson'].discreteFields[vorti].topology
topoadvec = op['advection'].discreteFields[vorti].topology
# ===== Penalization of the vorticity on a sphere inside the domain =====
from hysop.operator.penalize_vorticity import PenalizeVorticity
op['penalVort'] = PenalizeVorticity(velocity=velo, vorticity=vorti,
#from hysop.operator.penalize_vorticity import PenalizeVorticity
from hysop.operator.penalization_and_curl import PenalizationAndCurl
#op['penalVort'] = PenalizeVorticity(velocity=velo, vorticity=vorti,
op['penalVort'] = PenalizationAndCurl(velocity=velo, vorticity=vorti,
discretization=topo_with_ghosts,
obstacles=[sphere], coeff=1e8,
method={SpaceDiscretisation: FD_C_4})
......@@ -173,7 +176,7 @@ distr['fft2str'] = RedistributeIntra(source=op['poisson'],
target=op['stretching'],
variables=[velo, vorti])
distr['str2fft'] = RedistributeIntra(source=op['stretching'],
target=op['diffusion'],
target=op['poisson'],
variables=[velo, vorti])
distr['fft2advec'] = RedistributeIntra(source=op['poisson'],
target=op['advection'],
......@@ -182,15 +185,17 @@ distr['advec2fft'] = RedistributeIntra(source=op['advection'],
target=op['poisson'],
variables=[velo, vorti])
# ========= Monitoring operators =========
monitors = {}
iop = IO_params('fields', frequency=100)
monitors['writer'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
io_params=iop)
io_ener = IO_params('energy_enstrophy')
monitors['energy'] = EnergyEnstrophy(velo, vorti, discretization=topofft,
io_params=io_ener)
io_params=io_ener, is_normalized=False)
from hysop.domain.subsets.control_box import ControlBox
from hysop.operator.drag_and_lift import MomentumForces, NocaForces
ref_step = topo_with_ghosts.mesh.space_step
cbpos = npw.zeros(dim)
cblength = npw.zeros(dim)
......@@ -198,40 +203,51 @@ cbpos[...] = boxorigin[...]
cbpos += 15 * ref_step
cblength[...] = boxlength[...]
cblength -= 30 * ref_step
cb = ControlBox(box, cbpos, cblength)
cb = ControlBox(parent=box, origin=cbpos, length=cblength)
coeffForce = 1. / (0.5 * uinf ** 2 * pi * RADIUS ** 2)
io_forces=IO_params('drag_and_lift')
monitors['forces'] = DragAndLift(velo, vorti, VISCOSITY, coeffForce,
discretization=topo_with_ghosts, cb,
obstacles=[sphere], io_params=io_forces)
io_forcesPenal=IO_params('drag_and_lift_penal')
monitors['forcesPenal'] = DragAndLiftPenal(velo, vorti, coeffForce,
discretization=topofft,
obstacles=[sphere], factor=[1e8],
io_params=io_forcesPenal)
step_dir = ref_step[0]
io_sliceXY = IO_params('sliceXY', frequency=200)
thickSliceXY = ControlBox(box, origin=[-2.0, -2.56, -2.0 * step_dir],
lengths=[10.24, 5.12, 4.0 * step_dir])
monitors['writerSliceXY'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
io_params=io_sliceXY, subset=thickSliceXY,
xmfalways=True)
io_sliceXZ = IO_params('sliceXZ', frequency=2000)
thickSliceXZ = ControlBox(box, origin=[-2.0, -2.0 * step_dir, -2.56],
lengths=[10.24, 4.0 * step_dir, 5.12])
monitors['writerSliceXZ'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
io_params=io_sliceXZ, subset=thickSliceXZ,
xmfalways=True)
io_subBox = IO_params('subBox', frequency=2000)
subBox = ControlBox(box, origin=[-0.7, -2.0, -2.0], lengths=[8.0, 4.0, 4.0])
monitors['writerSubBox'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
io_params=io_subBox, subset=subBox,
xmfalways=True)
io_forces=IO_params('drag_and_lift_Noca')
#monitors['forcesNoca'] = NocaForces(velo, vorti,
# discretization=topo_with_ghosts,
# nu=VISCOSITY,
# volume_of_control=cb,
# normalization=coeffForce,
# obstacles=[sphere],
# io_params=io_forces)
io_forcesPenal=IO_params('drag_and_lift_Mom')
monitors['forcesMom'] = MomentumForces(velocity=velo,
discretization=topo_with_ghosts,
normalization=coeffForce,
obstacles=[sphere],
io_params=io_forcesPenal)
#io_forcesPenal=IO_params('drag_and_lift_penal')
#monitors['forcesPenal'] = DragAndLiftPenal(velo, vorti, coeffForce,
# discretization=topofft,
# obstacles=[sphere], factor=[1e8],
# io_params=io_forcesPenal)
#step_dir = ref_step[0]
#io_sliceXY = IO_params('sliceXY', frequency=200)
#thickSliceXY = ControlBox(box, origin=[-2.0, -2.56, -2.0 * step_dir],
# lengths=[10.24, 5.12, 4.0 * step_dir])
#monitors['writerSliceXY'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
# io_params=io_sliceXY, subset=thickSliceXY,
# xmfalways=True)
#io_sliceXZ = IO_params('sliceXZ', frequency=2000)
#thickSliceXZ = ControlBox(box, origin=[-2.0, -2.0 * step_dir, -2.56],
# lengths=[10.24, 4.0 * step_dir, 5.12])
#monitors['writerSliceXZ'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
# io_params=io_sliceXZ, subset=thickSliceXZ,
# xmfalways=True)
#io_subBox = IO_params('subBox', frequency=2000)
#subBox = ControlBox(box, origin=[-0.7, -2.0, -2.0], lengths=[8.0, 4.0, 4.0])
#monitors['writerSubBox'] = HDF_Writer(variables={velo: topofft, vorti: topofft},
# io_params=io_subBox, subset=subBox,
# xmfalways=True)
## ========= Monitors discretization=========
for monit in monitors.values():
......@@ -244,7 +260,9 @@ for ope in op.values():
for ope in distr.values():
ope.setup()
for monit in monitors.values():
monit.setUp()
monit.discretize()
for monit in monitors.values():
monit.setup()
print '[', main_rank, '] total time for setup:', MPI.Wtime() - time_setup
......@@ -268,14 +286,13 @@ print '[', main_rank, '] total time for init :', MPI.Wtime() - time_init
fullseq = []
def run(sequence):
op['poisson'].apply(simu) # Poisson
op['correction'].apply(simu) # Velo correction
monitors['forcesPenal'].apply(simu) # Forces Heloise
op['poisson'].apply(simu) # Poisson + correction
monitors['forcesMom'].apply(simu) # Forces Heloise
distr['fft2str'].apply(simu)
distr['fft2str'].wait()
op['penalVort'].apply(simu) # Vorticity penalization
op['stretching'].apply(simu) # Stretching
monitors['forces'].apply(simu) # Forces Noca
# monitors['forcesNoca'].apply(simu) # Forces Noca
distr['str2fft'].apply(simu)
distr['str2fft'].wait()
op['diffusion'].apply(simu) # Diffusion
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment