diff --git a/ci/scripts/build_and_debug.sh b/ci/scripts/build_and_debug.sh
index de66fc1046209834cc08435120916087392c244e..d99b6cb960aeed0db1dc82acd5a0cf3ee8acb6ba 100755
--- a/ci/scripts/build_and_debug.sh
+++ b/ci/scripts/build_and_debug.sh
@@ -27,6 +27,6 @@ cd -
rm -rf build
apt-get update
-apt-get install -y gdb python-dbg
+apt-get install -y gdb python2-dbg
bash
diff --git a/examples/flow_around_sphere/flow_around_sphere.py b/examples/flow_around_sphere/flow_around_sphere.py
index 005fc88435082cbdafd36b1f9d4cc67c7cfc58e7..5a5a90ce71ebb3aa084c872c0e927a72b0593aa3 100644
--- a/examples/flow_around_sphere/flow_around_sphere.py
+++ b/examples/flow_around_sphere/flow_around_sphere.py
@@ -209,9 +209,9 @@ correctFlowrate = FlowRateCorrection(
**extra_op_kwds)
#> outputs
-#penal.dump_inputs(fields=(sphere, ), frequency=outfrea)
-correctFlowrate.dump_inputs(fields=(vorti, ), frequency=outfreq)
-correctFlowrate.dump_outputs(fields=(velo,), frequency=outfreq)
+penal.dump_inputs(fields=(sphere, ), frequency=outfrea, io_params=args.io_params)
+correctFlowrate.dump_inputs(fields=(vorti, ), frequency=outfreq, io_params=args.io_params)
+correctFlowrate.dump_outputs(fields=(velo,), frequency=outfreq, io_params=args.io_params)
#> Operator to compute the infinite norm of the velocity
min_max_U = MinMaxFieldStatistics(name='min_max_U', field=velo,
diff --git a/hysop/operator/base/min_max.py b/hysop/operator/base/min_max.py
index 2ec39a464b69f1803725ca5a87b23e10c55884f3..090637b48dd18f833edf567e239b76d2f84a4f79 100644
--- a/hysop/operator/base/min_max.py
+++ b/hysop/operator/base/min_max.py
@@ -226,7 +226,7 @@ class MinMaxFieldStatisticsBase(object):
else:
comm = self.mpi_params.comm
Fmin, Fmax = self.Fmin, self.Fmax
- if (self.Fmax is not None):
+ if (self.Fmax is not None):
sendbuff = npw.zeros_like(Fmax.value)
recvbuff = npw.zeros_like(Fmax.value)
def _collect_max(val, sendbuff=sendbuff, recvbuff=recvbuff):
@@ -235,7 +235,7 @@ class MinMaxFieldStatisticsBase(object):
return recvbuff.copy()
else:
_collect_max = None
- if (self.Fmin is not None):
+ if (self.Fmin is not None):
sendbuff = npw.zeros_like(Fmin.value)
recvbuff = npw.zeros_like(Fmin.value)
def _collect_min(val, sendbuff=sendbuff, recvbuff=recvbuff):
diff --git a/hysop/tools/debug_dumper.py b/hysop/tools/debug_dumper.py
index 868179b3f0def67cf98a1d5ba8779f1b69122a10..389270cd245722e4ccf397a624f78fde60126195 100644
--- a/hysop/tools/debug_dumper.py
+++ b/hysop/tools/debug_dumper.py
@@ -13,7 +13,8 @@ from hysop.fields.discrete_field import DiscreteScalarFieldView
class DebugDumper(object):
def __init__(self, name, path, force_overwrite=False,
- enable_on_op_apply=False, dump_precision=10,
+ enable_on_op_apply=False,
+ dump_data=False, dump_precision=10,
comm=MPI.COMM_WORLD, io_leader=0):
assert isinstance(name, str), name
assert isinstance(path, str), path
@@ -27,6 +28,7 @@ class DebugDumper(object):
self.blobs_directory = blobs_directory
self.dump_id = 0
self.enable_on_op_apply = enable_on_op_apply
+ self.dump_data = dump_data
self.dump_precision = dump_precision
self.comm = comm
@@ -124,7 +126,7 @@ class DebugDumper(object):
if self.is_io_leader:
self.runfile.write(entry)
- if (comm_size == 1):
+ if (self.dump_data) and (comm_size == 1):
dst = '{}/{}'.format(self.blobs_directory, self.dump_id)
np.savez_compressed(dst, data=d)
diff --git a/hysop_examples/examples/flow_around_sphere/flow_around_sphere.py b/hysop_examples/examples/flow_around_sphere/flow_around_sphere.py
index 858e4bbf71e186cd381554276b68ccc966240156..49402ec47225ea80be1d1e44c1178b0bd3ba60f3 100644
--- a/hysop_examples/examples/flow_around_sphere/flow_around_sphere.py
+++ b/hysop_examples/examples/flow_around_sphere/flow_around_sphere.py
@@ -52,18 +52,18 @@ def compute(args):
elif (impl is Implementation.OPENCL):
# For the OpenCL implementation we need to setup the compute device
# and configure how the code is generated and compiled at runtime.
-
+
# Create an explicit OpenCL context from user parameters
from hysop.backend.device.opencl.opencl_tools import get_or_create_opencl_env, get_device_number
cl_env = get_or_create_opencl_env(
mpi_params=mpi_params,
- platform_id=args.cl_platform_id,
+ platform_id=args.cl_platform_id,
device_id=box.machine_rank%get_device_number() if args.cl_device_id is None else None)
-
+
# Configure OpenCL kernel generation and tuning (already done by HysopArgParser)
from hysop.methods import OpenClKernelConfig
method = { OpenClKernelConfig: args.opencl_kernel_config }
-
+
# Setup opencl specific extra operator keyword arguments
extra_op_kwds['cl_env'] = cl_env
else:
@@ -122,7 +122,7 @@ def compute(args):
def computeFlowrate(t, flowrate):
- fr = np.zeros(3)
+ fr = np.zeros_like(flowrate.value)
fr[-1] = uinf * box.length[1] * box.length[0]
Tstart=3.0
if t() >= Tstart and t() <= Tstart + 1.0:
@@ -171,7 +171,7 @@ def compute(args):
advec_dir = DirectionalAdvection(
implementation=impl,
name='advec',
- velocity = velo,
+ velocity = velo,
advected_fields = (vorti, ),
velocity_cfl = args.cfl,
variables = {velo: npts, vorti: npts},
@@ -182,7 +182,7 @@ def compute(args):
else:
StretchOp = StaticDirectionalStretching
stretch = StretchOp(
- implementation=Implementation.PYTHON if implIsFortran else impl,
+ implementation=Implementation.PYTHON if implIsFortran else impl,
name='stretch',
formulation=StretchingFormulation.CONSERVATIVE,
velocity=velo,
@@ -314,13 +314,13 @@ def compute(args):
simu.write_parameters(t, dt_cfl, dt_advec, dt, enstrophy, flowrate,
min_max_U.Finf, min_max_W.Finf, adapt_dt.equivalent_CFL,
filename='parameters.txt', precision=8)
-
+
problem.initialize_field(vorti, formula=computeVort)
problem.initialize_field(velo, formula=computeVel)
problem.initialize_field(sphere, formula=computeSphere)
# Finally solve the problem
- problem.solve(simu, dry_run=args.dry_run,
+ problem.solve(simu, dry_run=args.dry_run,
debug_dumper=args.debug_dumper,
checkpoint_handler=args.checkpoint_handler)
diff --git a/hysop_examples/examples/particles_above_salt/particles_above_salt_bc_3d.py b/hysop_examples/examples/particles_above_salt/particles_above_salt_bc_3d.py
index f7c1b6c54328948bc6b02751f7daf75df3a88c84..ac2f58ac627f30085bd3d9b028d1af897030d007 100644
--- a/hysop_examples/examples/particles_above_salt/particles_above_salt_bc_3d.py
+++ b/hysop_examples/examples/particles_above_salt/particles_above_salt_bc_3d.py
@@ -1,6 +1,6 @@
## See Meiburg 2012 & 2014
## Sediment-laden fresh water above salt water.
-
+
import numpy as np
import scipy as sp
import sympy as sm
@@ -25,7 +25,7 @@ def delta(Ys, l0):
for Yi in Ys:
Y0 = Y0*Yi
return 0.1*l0*(np.random.rand(*Y0.shape)-0.5)
-
+
def init_concentration(data, coords, l0, component):
assert (component==0)
X = coords[0]
@@ -59,7 +59,7 @@ def compute(args):
from hysop.methods import SpaceDiscretization, Remesh, TimeIntegrator, \
ComputeGranularity, Interpolation
-
+
from hysop.symbolic import sm, space_symbols, local_indices_symbols
from hysop.symbolic.base import SymbolicTensor
from hysop.symbolic.field import curl
@@ -67,7 +67,7 @@ def compute(args):
from hysop.symbolic.misc import Select
from hysop.symbolic.tmp import TmpScalar
from hysop.tools.string_utils import framed_str
-
+
## IO paths
spectral_path = IO.default_path() + '/spectral'
dump_energy_ioparams=IOParams(filepath=spectral_path, filename='E_{fname}', frequency=args.dump_freq)
@@ -96,7 +96,7 @@ def compute(args):
npts=N[::-1]
Xo=Xo[::-1]
Xn=Xn[::-1]
-
+
lboundaries = (BoxBoundaryCondition.PERIODIC,)*(dim-1)+(BoxBoundaryCondition.SYMMETRIC,)
rboundaries = (BoxBoundaryCondition.PERIODIC,)*(dim-1)+(BoxBoundaryCondition.SYMMETRIC,)
@@ -107,7 +107,7 @@ def compute(args):
box = Box(origin=Xo, length=np.subtract(Xn,Xo),
lboundaries=lboundaries, rboundaries=rboundaries)
-
+
# Get default MPI Parameters from domain (even for serial jobs)
mpi_params = MPIParams(comm=box.task_comm,
task_id=box.current_task())
@@ -122,17 +122,17 @@ def compute(args):
elif (impl is Implementation.OPENCL):
# For the OpenCL implementation we need to setup the compute device
# and configure how the code is generated and compiled at runtime.
-
+
# Create an explicit OpenCL context from user parameters
from hysop.backend.device.opencl.opencl_tools import get_or_create_opencl_env
- cl_env = get_or_create_opencl_env(mpi_params=mpi_params,
- platform_id=args.cl_platform_id,
+ cl_env = get_or_create_opencl_env(mpi_params=mpi_params,
+ platform_id=args.cl_platform_id,
device_id=args.cl_device_id)
-
+
# Configure OpenCL kernel generation and tuning (already done by HysopArgParser)
from hysop.methods import OpenClKernelConfig
method = { OpenClKernelConfig: args.opencl_kernel_config }
-
+
# Setup opencl specific extra operator keyword arguments
extra_op_kwds['cl_env'] = cl_env
else:
@@ -145,7 +145,7 @@ def compute(args):
vorti = VorticityField(velocity=velo)
C = Field(domain=box, name='C', dtype=args.dtype, lboundaries=C_lboundaries, rboundaries=C_rboundaries)
S = Field(domain=box, name='S', dtype=args.dtype, lboundaries=S_lboundaries, rboundaries=S_rboundaries)
-
+
# Symbolic fields
frame = velo.domain.frame
Us = velo.s(*frame.vars)
@@ -155,43 +155,43 @@ def compute(args):
dts = dt.s
### Build the directional operators
- #> Directional advection
+ #> Directional advection
advec = DirectionalAdvection(implementation=impl,
name='advec',
- velocity = velo,
+ velocity = velo,
advected_fields = (vorti,S),
velocity_cfl = args.cfl,
variables = {velo: npts, vorti: npts, S: npts},
dt=dt, **extra_op_kwds)
-
+
V0 = [0]*dim
VP = [0]*dim
VP[0] = Vp
advec_C = DirectionalAdvection(implementation=impl,
name='advec_C',
- velocity = velo,
+ velocity = velo,
advected_fields = (C,),
relative_velocity = VP,
velocity_cfl = args.cfl,
variables = {velo: npts, C: npts},
dt=dt, **extra_op_kwds)
-
+
#> Stretch vorticity
if (dim==3):
stretch = DirectionalStretching(implementation=impl,
name='stretch',
pretty_name='stretch',
formulation = args.stretching_formulation,
- velocity = velo,
+ velocity = velo,
vorticity = vorti,
variables = {velo: npts, vorti: npts},
dt=dt, **extra_op_kwds)
elif (dim==2):
- stretch = None
+ stretch = None
else:
msg='Unsupported dimension {}.'.format(dim)
raise RuntimeError(msg)
-
+
#> Diffusion of vorticity, S and C
diffuse_S = Diffusion(implementation=impl,
enforce_implementation=enforce_implementation,
@@ -200,9 +200,9 @@ def compute(args):
nu = nu_S,
Fin = S,
variables = {S: npts},
- dt=dt,
+ dt=dt,
dump_energy=dump_energy_ioparams,
- plot_inout_energy=IOParams(filepath=spectral_path,
+ plot_inout_energy=IOParams(filepath=spectral_path,
filename='E_S_diffusion_{ite}', frequency=args.dump_freq),
**extra_op_kwds)
diffuse_C = Diffusion(implementation=impl,
@@ -214,7 +214,7 @@ def compute(args):
variables = {C: npts},
dt=dt,
dump_energy=dump_energy_ioparams,
- plot_inout_energy=IOParams(filepath=spectral_path,
+ plot_inout_energy=IOParams(filepath=spectral_path,
filename='E_C_diffusion_{ite}', frequency=args.dump_freq),
**extra_op_kwds)
@@ -225,33 +225,33 @@ def compute(args):
lhs = Ws.diff(frame.time)
rhs = curl(Fext, frame)
exprs = Assignment.assign(lhs, rhs)
- external_force = DirectionalSymbolic(name='Fext',
+ external_force = DirectionalSymbolic(name='Fext',
implementation=impl,
exprs=exprs, dt=dt,
variables={vorti: npts,
S: npts,
C: npts},
**extra_op_kwds)
-
- splitting = StrangSplitting(splitting_dim=dim,
+
+ splitting = StrangSplitting(splitting_dim=dim,
order=args.strang_order)
splitting.push_operators(advec, advec_C, stretch, external_force)
-
+
### Build standard operators
#> Poisson operator to recover the velocity from the vorticity
poisson = PoissonCurl(name='poisson',
- velocity=velo, vorticity=vorti,
+ velocity=velo, vorticity=vorti,
variables={velo:npts, vorti: npts},
diffusion=nu_W, dt=dt,
- implementation=impl,
+ implementation=impl,
enforce_implementation=enforce_implementation,
dump_energy=dump_energy_ioparams,
- plot_velocity_energy=IOParams(filepath=spectral_path,
+ plot_velocity_energy=IOParams(filepath=spectral_path,
filename='E_velocity_{ite}', frequency=args.dump_freq),
- plot_inout_vorticity_energy=IOParams(filepath=spectral_path,
+ plot_inout_vorticity_energy=IOParams(filepath=spectral_path,
filename='E_vorticity_{ite}', frequency=args.dump_freq),
**extra_op_kwds)
-
+
#> Operator to compute the infinite norm of the velocity
min_max_U = MinMaxFieldStatistics(name='min_max_U', field=velo,
Finf=True, implementation=impl, variables={velo:npts},
@@ -260,15 +260,15 @@ def compute(args):
min_max_W = MinMaxFieldStatistics(field=vorti,
Finf=True, implementation=impl, variables={vorti:npts},
**extra_op_kwds)
-
+
#> Operators to dump all fields
io_params = IOParams(filename='fields', frequency=args.dump_freq)
dump_fields = HDF_Writer(name='dump',
io_params=io_params,
force_backend=Backend.OPENCL,
variables={vorti: npts,
- velo: npts,
- C: npts,
+ velo: npts,
+ C: npts,
S: npts},
**extra_op_kwds)
@@ -276,23 +276,23 @@ def compute(args):
adapt_dt = AdaptiveTimeStep(dt, equivalent_CFL=True,
name='merge_dt', pretty_name='dt',
max_dt=1.0)
- dt_cfl = adapt_dt.push_cfl_criteria(cfl=args.cfl,
+ dt_cfl = adapt_dt.push_cfl_criteria(cfl=args.cfl,
Fmin=min_max_U.Fmin,
Fmax=min_max_U.Fmax,
- equivalent_CFL=True,
+ equivalent_CFL=True,
relative_velocities=[V0, VP],
name='dt_cfl', pretty_name='CFL')
dt_advec = adapt_dt.push_advection_criteria(lcfl=args.lcfl, Finf=min_max_W.Finf,
criteria=AdvectionCriteria.W_INF,
name='dt_lcfl', pretty_name='LCFL')
-
- ## Create the problem we want to solve and insert our
+
+ ## Create the problem we want to solve and insert our
# directional splitting subgraph and the standard operators.
# The method dictionnary passed to this graph will be dispatched
# accross all operators contained in the graph.
method.update(
- {
+ {
ComputeGranularity: args.compute_granularity,
SpaceDiscretization: args.fd_order,
TimeIntegrator: args.time_integrator,
@@ -301,41 +301,41 @@ def compute(args):
)
problem = Problem(method=method)
- problem.insert(poisson,
+ problem.insert(poisson,
diffuse_S, diffuse_C,
# dump_fields,
- splitting,
+ splitting,
min_max_U, min_max_W, adapt_dt)
problem.build(args)
-
+
# If a visu_rank was provided, and show_graph was set,
# display the graph on the given process rank.
if args.display_graph:
problem.display(args.visu_rank)
-
+
# Create a simulation
# (do not forget to specify the t and dt parameters here)
- simu = Simulation(start=args.tstart, end=args.tend,
+ simu = Simulation(start=args.tstart, end=args.tend,
nb_iter=args.nb_iter,
max_iter=args.max_iter,
dt0=args.dt, times_of_interest=args.times_of_interest,
t=t, dt=dt)
- simu.write_parameters(t, dt_cfl, dt_advec, dt,
+ simu.write_parameters(t, dt_cfl, dt_advec, dt,
min_max_U.Finf, min_max_W.Finf, adapt_dt.equivalent_CFL,
filename='parameters.txt', precision=8)
-
+
# Initialize vorticity, velocity, S and C on all topologies
problem.initialize_field(field=velo, formula=init_velocity)
problem.initialize_field(field=vorti, formula=init_vorticity)
problem.initialize_field(field=C, formula=init_concentration, l0=l0)
problem.initialize_field(field=S, formula=init_salinity, l0=l0)
-
+
# Finally solve the problem
- problem.solve(simu, dry_run=args.dry_run,
+ problem.solve(simu, dry_run=args.dry_run,
debug_dumper=args.debug_dumper,
checkpoint_handler=args.checkpoint_handler)
-
+
# Finalize
problem.finalize()
@@ -346,23 +346,23 @@ if __name__=='__main__':
class ParticleAboveSaltArgParser(HysopArgParser):
def __init__(self):
prog_name = 'particle_above_salt_bc_3d'
- default_dump_dir = '{}/hysop_examples/{}'.format(HysopArgParser.tmp_dir(),
+ default_dump_dir = '{}/hysop_examples/{}'.format(HysopArgParser.tmp_dir(),
prog_name)
description=colors.color('HySoP Particles Above Salt Example: ', fg='blue',
style='bold')
description+=colors.color('[Meiburg 2014]', fg='yellow', style='bold')
- description+=colors.color('\nSediment-laden fresh water above salt water.',
+ description+=colors.color('\nSediment-laden fresh water above salt water.',
fg='yellow')
description+='\n'
description+='\nThis example focuses on a validation study for the '
description+='hybrid particle-mesh vortex method in the Boussinesq approximation.'
-
+
super(ParticleAboveSaltArgParser, self).__init__(
prog_name=prog_name,
description=description,
default_dump_dir=default_dump_dir)
-
+
def _add_main_args(self):
args = super(ParticleAboveSaltArgParser, self)._add_main_args()
args.add_argument('-Sc', '--schmidt', type=float,
@@ -378,7 +378,7 @@ if __name__=='__main__':
dest='Rs',
help='Density expension factor.')
return args
-
+
def _check_main_args(self, args):
super(ParticleAboveSaltArgParser, self)._check_main_args(args)
self._check_default(args, ('schmidt', 'tau', 'Vp', 'Rs'), float, allow_none=False)
@@ -395,11 +395,11 @@ if __name__=='__main__':
parser = ParticleAboveSaltArgParser()
parser.set_defaults(impl='cl', ndim=3, npts=(64,),
- box_origin=(0.0,), box_length=(1.0,),
- tstart=0.0, tend=201.0,
+ box_origin=(0.0,), box_length=(1.0,),
+ tstart=0.0, tend=201.0,
dt=1e-6, cfl=12.00, lcfl=0.95,
dump_times=(25.0, 50.0, 75.0, 100.0, 125.0, 150.0, 175.0, 200.0),
- dump_freq=0,
+ dump_freq=0,
schmidt=7.0, tau=25.0, Vp=0.04, Rs=2.0)
parser.run(compute)