diff --git a/Examples/levelSet2D.py b/Examples/levelSet2D.py
index fad3310b6a6244b54db64084b3a0fe0ac7cfc485..7983226451d7f03f1a965fc5617eb6cf63cb0794 100755
--- a/Examples/levelSet2D.py
+++ b/Examples/levelSet2D.py
@@ -10,7 +10,13 @@ from parmepy.problem.transport import TransportProblem
from parmepy.operator.monitors.printer import Printer
from parmepy.problem.simulation import Simulation
from parmepy.operator.analytic import Analytic
-# import math
+from parmepy.methods_keys import TimeIntegrator, Interpolation, Remesh,\
+ Support, Splitting
+from parmepy.numerics.integrators.runge_kutta2 import RK2
+from parmepy.numerics.integrators.runge_kutta4 import RK4
+from parmepy.numerics.interpolation import Linear
+from parmepy.numerics.remeshing import L2_1, L2_2, L4_2, L4_4, L6_4, L6_6, L8_4
+import math
# def vitesse(x, y):
# vx = -math.sin(x * math.pi) ** 2 * math.sin(y * math.pi * 2)
@@ -18,76 +24,75 @@ from parmepy.operator.analytic import Analytic
# return vx, vy
-# def vitesse(x, y, t=0):
-# vx = -math.sin(x * math.pi) ** 2 * math.sin(y * math.pi * 2) * \
-# math.cos(t * math.pi / 3.)
-# vy = math.sin(y * math.pi) ** 2 * math.sin(x * math.pi * 2) * \
-# math.cos(t * math.pi / 3.)
-# return vx, vy
+def vitesse(x, y, t=0):
+ vx = -math.sin(x * math.pi) ** 2 * math.sin(y * math.pi * 2) * \
+ math.cos(t * math.pi / 3.)
+ vy = math.sin(y * math.pi) ** 2 * math.sin(x * math.pi * 2) * \
+ math.cos(t * math.pi / 3.)
+ return vx, vy
-# def scalaire(x, y):
-# rr = math.sqrt((x - 0.5) ** 2 + (y - 0.75) ** 2)
-# if rr < 0.15:
-# return 1.
-# else:
-# return 0.
+def scalaire(x, y):
+ rr = math.sqrt((x - 0.5) ** 2 + (y - 0.75) ** 2)
+ if rr < 0.1:
+ return 1.
+ else:
+ return 0.
dim = 2
boxLength = [1., 1.]
boxMin = [0., 0.]
-nbElem = [1025, 1025]
+nbElem = [129, 129]
timeStep = 0.025
finalTime = 3.
outputFilePrefix = './res_2D/levelSet_2D_rect_'
-outputModulo = 0
-simu = Simulation(tinit=0.0, tend=finalTime, timeStep=timeStep)
+outputModulo = 1
+simu = Simulation(tinit=0.0, tend=finalTime, timeStep=timeStep, iterMax=10000)
## Domain
box = Box(dim, length=boxLength, origin=boxMin)
## Fields
-scal = Field(domain=box, name='Scalar')
-velo = Field(domain=box, name='Velocity', isVector=True)
+scal = Field(domain=box, name='Scalar', formula=scalaire)
+velo = Field(domain=box, name='Velocity', isVector=True,formula=vitesse)
## Operators
advec = Advection(velo, scal,
resolutions={velo: nbElem,
scal: nbElem},
- method='gpu_2k_m4prime',
- #method='gpu_1k_m6prime',
- #method='gpu_1k_m8prime',
- #method='gpu_2k_m4prime',
- #method='gpu_2k_m6prime',
- #method='gpu_2k_m8prime',
- #method='scales'
- src=['./levelSet2D.cl'],
- precision=PARMES_REAL_GPU,
- #precision=PARMES_DOUBLE_GPU,
- #splittingConfig='o2'
- #splittingConfig='y_only'
- #splittingConfig='o2_FullHalf'
+ method={TimeIntegrator: RK4,
+ Interpolation: Linear,
+ Remesh: L4_4,
+ Support: '',
+ Splitting: 'o2_FullHalf'},
+#src=['./levelSet2D.cl'],
+#precision=PARMES_REAL_GPU,
)
velocity = Analytic(velo,
resolutions={velo: nbElem},
- method='gpu',
- src=['./levelSet2D.cl'],
- precision=PARMES_REAL_GPU,
+#method={Support: 'gpu'},
+#src=['./levelSet2D.cl'],
+#precision=PARMES_REAL_GPU,
)
##Problem
#pb = TransportProblem([advec],simu,
-pb = TransportProblem([velocity, advec], simu,
- monitors=[Printer(fields=[scal],
- frequency=outputModulo,
- prefix=outputFilePrefix)])
+pb = TransportProblem([velocity, advec], simu, dumpFreq=10000)
## Setting solver to Problem
pb.setUp()
+p = Printer(variables=[scal],
+ topo=scal.topoInit,
+ frequency=outputModulo,
+ prefix=outputFilePrefix,
+ ext='.dat')
+pb.addMonitors([p])
+
## Solve problem
+p._printStep()
pb.solve()
pb.finalize()
diff --git a/Examples/levelSet3D.py b/Examples/levelSet3D.py
index bf1829342ea33bee5383129dd7d3ff41a7db0984..43aa6df168d37472339583125f487c510fb460e6 100755
--- a/Examples/levelSet3D.py
+++ b/Examples/levelSet3D.py
@@ -7,72 +7,80 @@ from parmepy.problem.transport import TransportProblem
from parmepy.operator.monitors.printer import Printer
from parmepy.operator.analytic import Analytic
from parmepy.problem.simulation import Simulation
-# from math import sin, cos, pi
+from parmepy.methods_keys import TimeIntegrator, Interpolation, Remesh,\
+ Support, Splitting
+from parmepy.numerics.integrators.runge_kutta2 import RK2
+from parmepy.numerics.integrators.runge_kutta4 import RK4
+from parmepy.numerics.interpolation import Linear
+from parmepy.numerics.remeshing import L2_1, L2_2, L4_2, L4_4, L6_4, L6_6, L8_4
+from math import sin, cos, pi, sqrt
-# def scalaire(x,y,z):
-# return 1.
+def scalaire(x,y,z):
+ r = sqrt((x - 0.35) ** 2 + (y - 0.35) ** 2 + (z - 0.35) ** 2)
+ if r < 0.15:
+ return 1.
+ else:
+ return 0.
-# def vitesse(x, y, z, t=0):
-# vx = 2. * sin(pi*x)**2*sin(2.*pi*y)*sin(2.*pi*z)*cos(t*pi/3.)
-# vy = -sin(2.*pi*x)*sin(pi*y)**2*sin(2.*pi*z)*cos(t*pi/3.)
-# vz = -sin(2.*pi*x)*sin(2.*pi*y)*sin(pi*z)**2*cos(t*pi/3.)
-# return vx, vy, vz
+def vitesse(x, y, z, t=0):
+ vx = 2. * sin(pi*x)**2*sin(2.*pi*y)*sin(2.*pi*z)*cos(t*pi/3.)
+ vy = -sin(2.*pi*x)*sin(pi*y)**2*sin(2.*pi*z)*cos(t*pi/3.)
+ vz = -sin(2.*pi*x)*sin(2.*pi*y)*sin(pi*z)**2*cos(t*pi/3.)
+ return vx, vy, vz
dim = 3
boxLength = [1., 1., 1.]
boxMin = [0., 0., 0.]
-nbElem = 3 * [129]
+nbElem = 3 * [65]
-timeStep = 0.05
+timeStep = 0.025
finalTime = 3.
-outputFilePrefix = './res3D_new/levelSet_3D_'
+outputFilePrefix = './res3D/levelSet_3D'
outputModulo = 1
-simu = Simulation(tinit=0.0, tend=finalTime, timeStep=timeStep)
+simu = Simulation(tinit=0.0, tend=finalTime, timeStep=timeStep, iterMax=10000)
## Domain
box = pp.Box(dim, length=boxLength, origin=boxMin)
## Fields
-scal = Field(domain=box, name='Scalar')
-velo = Field(domain=box, name='Velocity', isVector=True)
+scal = Field(domain=box, name='Scalar', formula=scalaire)
+velo = Field(domain=box, name='Velocity', isVector=True, formula=vitesse)
## Operators
advec = Advection(velo, scal,
resolutions={velo: nbElem,
scal: nbElem},
- #method='scales',
- method='gpu_1k_m4prime',
- #method='gpu_1k_m6prime',
- #method='gpu_1k_m8prime',
- #method='gpu_2k_m4prime',
- #method='gpu_2k_m6prime',
- #method='gpu_2k_m8prime',
- src=['./levelSet3D.cl'],
- precision=PARMES_REAL_GPU,
- #precision=PARMES_DOUBLE_GPU,
- splittingConfig='o2'
- #splittingConfig='o2_FullHalf'
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L6_4,
+ Support: '',
+ Splitting: 'o2_FullHalf'},
+#src=['./levelSet3D.cl'],
+#precision=PARMES_REAL_GPU,
)
velocity = Analytic(velo,
resolutions={velo: nbElem},
- method='gpu',
- src=['./levelSet3D.cl'],
- precision=PARMES_REAL_GPU,
+#method={Support: 'gpu'},
+#src=['./levelSet3D.cl'],
+#precision=PARMES_REAL_GPU,
)
##Problem
#pb = TransportProblem([advec],
-pb = TransportProblem([velocity, advec], simu,
- monitors=[Printer(fields=[scal],
- frequency=outputModulo,
- prefix=outputFilePrefix)])
-
-## Setting solver to Problem
+pb = TransportProblem([velocity, advec], simu, dumpFreq=10000)
pb.setUp()
+p=Printer(variables=[scal],
+ topo=scal.topoInit,
+ frequency=outputModulo,
+ prefix=outputFilePrefix,
+ ext=".h5")
+pb.addMonitors([p])
+
## Solve problem
+p._printStep()
pb.solve()
pb.finalize()
diff --git a/HySoP/hysop/gpu/gpu_particle_advection_2k.py b/HySoP/hysop/gpu/gpu_particle_advection_2k.py
index f83a804c240a7122ed0cb542514a3615980496cf..8aec4ee4cce0c44a346decaf7a95414c43cc51b5 100644
--- a/HySoP/hysop/gpu/gpu_particle_advection_2k.py
+++ b/HySoP/hysop/gpu/gpu_particle_advection_2k.py
@@ -84,7 +84,7 @@ class GPUParticleAdvection2k(GPUParticleAdvection):
name="Particle_Position",
isVector=False
).discretize(self.advectedFields[0].topology)
- GPUDiscreteField.fromField(self.cl_env,self.part_position,
+ GPUDiscreteField.fromField(self.cl_env, self.part_position,
self.gpu_precision)
## Result scalar
if self.part_advectedFields is None:
@@ -93,7 +93,7 @@ class GPUParticleAdvection2k(GPUParticleAdvection):
name="Particle_AdvectedFields",
isVector=self.advectedFields[0].isVector
).discretize(self.advectedFields[0].topology)]
- GPUDiscreteField.fromField(self.cl_env,self.part_advectedFields[0],
+ GPUDiscreteField.fromField(self.cl_env, self.part_advectedFields[0],
self.gpu_precision, layout=False)
self.variables = [self.advectedFields[0], self.velocity,
diff --git a/HySoP/hysop/gpu/tests/test_advection_randomVelocity.py b/HySoP/hysop/gpu/tests/test_advection_randomVelocity.py
index 39e5e9f918a7fad8a423fd70e4d3f5d904530233..6c146b8ede63206b95cea634b82799fb8e702f61 100644
--- a/HySoP/hysop/gpu/tests/test_advection_randomVelocity.py
+++ b/HySoP/hysop/gpu/tests/test_advection_randomVelocity.py
@@ -2,6 +2,8 @@
@file parmepy.gpu.tests.test_advection_randomVelocity
Testing advection kernels with a random velocity field.
"""
+import parmepy
+parmepy.constants.PARMES_REAL = parmepy.constants.np.float32
from parmepy.domain.box import Box
from parmepy.fields.continuous import Field
from parmepy.operator.advection import Advection
@@ -93,7 +95,7 @@ def assertion_3D_withPython(scal, velo, advec, advec_py):
# M6 testing
-def fixMe_test_2D_m6_1k():
+def test_2D_m6_1k():
"""
Testing M6 remeshing formula in 2D, 1 kernel, simple precision,
o2 splitting.
@@ -122,7 +124,7 @@ def fixMe_test_2D_m6_1k():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m6_2k():
+def test_2D_m6_2k():
"""
Testing M6 remeshing formula in 2D, 2 kernel, simple precision,
o2 splitting.
@@ -151,7 +153,7 @@ def fixMe_test_2D_m6_2k():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m6_1k_sFH():
+def test_2D_m6_1k_sFH():
"""
Testing M6 remeshing formula in 2D, 1 kernel, simple precision,
o2_FullHalf splitting.
@@ -180,7 +182,7 @@ def fixMe_test_2D_m6_1k_sFH():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m6_2k_sFH():
+def test_2D_m6_2k_sFH():
"""
Testing M6 remeshing formula in 2D, 2 kernel, simple precision,
o2_FullHalf splitting.
@@ -190,15 +192,21 @@ def fixMe_test_2D_m6_2k_sFH():
advec = Advection(velo, scal,
resolutions={velo: [33, 33],
scal: [33, 33]},
- method='gpu_2k_l4_2',
- precision=PARMES_REAL_GPU,
- splittingConfig="o2_FullHalf"
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L4_2,
+ Support: 'gpu_2k',
+ Splitting: 'o2_FullHalf'},
+ precision=PARMES_REAL_GPU
)
advec_py = Advection(velo, scal,
resolutions={velo: [33, 33],
scal: [33, 33]},
- method='l4_2',
- splittingConfig="o2_FullHalf"
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L4_2,
+ Support: '',
+ Splitting: 'o2_FullHalf'}
)
assertion_2D_withPython(scal, velo, advec, advec_py)
@@ -320,7 +328,7 @@ def test_3D_m6_2k_sFH():
# M4 testing
-def fixMe_test_2D_m4_1k():
+def test_2D_m4_1k():
"""
Testing M4 remeshing formula in 2D, 1 kernel, simple precision,
o2 splitting.
@@ -349,7 +357,7 @@ def fixMe_test_2D_m4_1k():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m4_2k():
+def test_2D_m4_2k():
"""
Testing M4 remeshing formula in 2D, 2 kernel, simple precision,
o2 splitting.
@@ -378,7 +386,7 @@ def fixMe_test_2D_m4_2k():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m4_1k_sFH():
+def test_2D_m4_1k_sFH():
"""
Testing M4 remeshing formula in 2D, 1 kernel, simple precision,
o2_FullHalf splitting.
@@ -407,7 +415,7 @@ def fixMe_test_2D_m4_1k_sFH():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m4_2k_sFH():
+def test_2D_m4_2k_sFH():
"""
Testing M4 remeshing formula in 2D, 2 kernel, simple precision,
o2_FullHalf splitting.
@@ -553,7 +561,7 @@ def test_3D_m4_2k_sFH():
# M8 testing
-def fixMe_test_2D_m8_1k():
+def test_2D_m8_1k():
"""
Testing M8 remeshing formula in 2D, 1 kernel, simple precision,
o2 splitting.
@@ -582,7 +590,7 @@ def fixMe_test_2D_m8_1k():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m8_2k():
+def test_2D_m8_2k():
"""
Testing M8 remeshing formula in 2D, 2 kernel, simple precision,
o2 splitting.
@@ -611,7 +619,7 @@ def fixMe_test_2D_m8_2k():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m8_1k_sFH():
+def test_2D_m8_1k_sFH():
"""
Testing M8 remeshing formula in 2D, 1 kernel, simple precision,
o2_FullHalf splitting.
@@ -640,7 +648,7 @@ def fixMe_test_2D_m8_1k_sFH():
assertion_2D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_2D_m8_2k_sFH():
+def test_2D_m8_2k_sFH():
"""
Testing M8 remeshing formula in 2D, 2 kernel, simple precision,
o2_FullHalf splitting.
@@ -785,7 +793,7 @@ def test_3D_m8_2k_sFH():
assertion_3D_withPython(scal, velo, advec, advec_py)
-def fixMe_test_rectangular_domain2D():
+def test_rectangular_domain2D():
box = Box(2, length=[1., 1.], origin=[0., 0.])
scal = Field(domain=box, name='Scalar')
velo = Field(domain=box, name='Velocity', isVector=True)
@@ -795,9 +803,9 @@ def fixMe_test_rectangular_domain2D():
scal: [65, 33]},
method={TimeIntegrator: RK2,
Interpolation: Linear,
- Remesh: L4_2,
- Support: 'gpu_1k',
- Splitting: 'o2'},
+ Remesh: L2_1,
+ Support: 'gpu_2k',
+ Splitting: 'x_only'},
precision=PARMES_REAL_GPU,
)
advec_py = Advection(velo, scal,
@@ -805,9 +813,9 @@ def fixMe_test_rectangular_domain2D():
scal: [65, 33]},
method={TimeIntegrator: RK2,
Interpolation: Linear,
- Remesh: L4_2,
+ Remesh: L2_1,
Support: '',
- Splitting: 'o2'},
+ Splitting: 'x_only'},
)
advec.discretize()
advec_py.discretize()
@@ -824,20 +832,21 @@ def fixMe_test_rectangular_domain2D():
velo_d.data[1][...] = np.asarray(
np.random.random(velo_d.data[1].shape),
dtype=PARMES_REAL_GPU, order=ORDER) / (2. * scal_d.resolution[1])
+
scal_d.toDevice()
velo_d.toDevice()
- advec_py.apply(Simulation(0., 0.01, 1))
- advec.apply(Simulation(0., 0.01, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+ advec.apply(Simulation(0., 0.1, 1))
py_res = scal_d.data[0].copy()
scal_d.toHost()
- assert np.allclose(py_res, scal_d.data[0], rtol=5e-02, atol=5e-05)
+ assert np.allclose(py_res, scal_d.data[0], rtol=5e-05, atol=5e-05)
advec.finalize()
-def fixMe_test_rectangular_domain3D():
+def test_rectangular_domain3D():
box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
scal = Field(domain=box, name='Scalar')
velo = Field(domain=box, name='Velocity', isVector=True)
@@ -892,7 +901,7 @@ def fixMe_test_rectangular_domain3D():
advec.finalize()
-def fixMe_test_vector_2D():
+def test_vector_2D():
box = Box(2, length=[1., 1.], origin=[0., 0.])
scal = Field(domain=box, name='Scalar', isVector=True)
velo = Field(domain=box, name='Velocity', isVector=True)
@@ -950,7 +959,7 @@ def fixMe_test_vector_2D():
advec.finalize()
-def fixMe_test_vector_3D():
+def test_vector_3D():
box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
scal = Field(domain=box, name='Scalar', isVector=True)
velo = Field(domain=box, name='Velocity', isVector=True)
diff --git a/HySoP/hysop/operator/advection.py b/HySoP/hysop/operator/advection.py
index d32f39c5d8d4fb64305033af02f5e4cb5d5d34a6..9c1f85f8546f42b5666ea4e16dd5f5d021352944 100644
--- a/HySoP/hysop/operator/advection.py
+++ b/HySoP/hysop/operator/advection.py
@@ -71,18 +71,18 @@ class Advection(Operator):
"""
v = [velocity]
if isinstance(advectedFields, list):
- [v.append(f) for f in advectedFields]
+ self.advectedFields = advectedFields
else:
- v.append(advectedFields)
+ self.advectedFields = [advectedFields]
+ for adF in self.advectedFields:
+ v.append(adF)
Operator.__init__(self, v, method, topo=topo,
ghosts=ghosts)
## Transport velocity
self.velocity = velocity
## Transported fields
- self.advectedFields = [v for v in self.variables
- if not v is self.velocity]
self.output = self.advectedFields
- self.input = [v for v in self.variables]
+ self.input = [var for var in self.variables]
## Resolutions of the fields. Dictionnary with
## fields as keys.
## Example :
diff --git a/HySoP/hysop/operator/discrete/particle_advection.py b/HySoP/hysop/operator/discrete/particle_advection.py
index 199850f294c3d014f451c0615dc92d77e3692170..49d4ada6ec88ba388d9a280398ae128306527eee 100644
--- a/HySoP/hysop/operator/discrete/particle_advection.py
+++ b/HySoP/hysop/operator/discrete/particle_advection.py
@@ -123,25 +123,25 @@ class ParticleAdvection(DiscreteOperator):
assert wk.shape == tuple(memshape)
# Distribute work arrays
- interpol_w = self._workspace[:w_interp]
- interpol_iw = self._iworkspace[:iw_interp]
- integrator_w = self._workspace[w_interp:w_interp + w_rk]
- remesh_w = self._workspace[:w_remesh]
- remesh_iw = self._iworkspace[:iw_remesh]
+ self._interpol_w = self._workspace[:w_interp]
+ self._interpol_iw = self._iworkspace[:iw_interp]
+ self._integrator_w = self._workspace[w_interp:w_interp + w_rk]
+ self._remesh_w = self._workspace[:w_remesh]
+ self._remesh_iw = self._iworkspace[:iw_remesh]
self.num_interpolate = \
Interpol(self.velocity, self.dir,
self.advectedFields[0].topology,
- work=interpol_w,
- iwork=interpol_iw)
- self.num_advec = RKn(1, work=integrator_w,
+ work=self._interpol_w,
+ iwork=self._interpol_iw)
+ self.num_advec = RKn(1, work=self._integrator_w,
f=self.num_interpolate,
topo=self.velocity.topology, optim=WITH_GUESS)
self.num_remesh = Rmsh(
self.velocity.dimension,
self.advectedFields[0].topology,
- self.dir, work=remesh_w,
- iwork=remesh_iw)
+ self.dir, work=self._remesh_w,
+ iwork=self._remesh_iw)
self._isUpToDate = True
@@ -178,7 +178,7 @@ class ParticleAdvection(DiscreteOperator):
# Advect particles
# RK use the first 2 (or 3) works and leave others to interpolation
# First work contains fist evaluation of ode right hand side.
- self._workspace[0][...] = self.velocity.data[self.dir]
+ self._integrator_w[0][...] = self.velocity.data[self.dir][...]
self.part_position.data = self.num_advec(
t, self.part_position.data, dt, result=self.part_position.data)
diff --git a/HySoP/hysop/operator/monitors/printer.py b/HySoP/hysop/operator/monitors/printer.py
index 1d83ea624c43fe329f81626182707ce83027cd5d..07d598682d314d9301c91a84985a72eee818f954 100644
--- a/HySoP/hysop/operator/monitors/printer.py
+++ b/HySoP/hysop/operator/monitors/printer.py
@@ -159,8 +159,6 @@ class Printer(Monitoring):
evtk.imageToVTK(filename, origin=orig, spacing=spacing,
pointData=self._build_vtk_dict())
elif self.ext == '.h5':
- if __VERBOSE__:
- print 'printing output xmf+hdf5 files ...'
filename = self.prefix + "_iter_{0:03d}".format(ite) + self.ext
# Write the h5 file
# (force np.float64, ParaView seems to not be able to read float32)
diff --git a/HySoP/hysop/operator/tests/test_advec_scales.py b/HySoP/hysop/operator/tests/test_advec_scales.py
index fbfedb53bcfc74dd60f3a55fbb15f7b9ad248287..dfd9b92b7802f4d7da7f042425347e451058c283 100755
--- a/HySoP/hysop/operator/tests/test_advec_scales.py
+++ b/HySoP/hysop/operator/tests/test_advec_scales.py
@@ -4,7 +4,7 @@ Testing Scales advection operator.
import numpy as np
from parmepy.constants import PARMES_REAL, ORDER
from parmepy.methods_keys import Scales, TimeIntegrator, Interpolation,\
- Remesh, Support
+ Remesh, Support, Splitting
from parmepy.methods import RK2, L2_1, L4_2, M8Prime, Linear
from parmepy.domain.box import Box
from parmepy.fields.continuous import Field
@@ -160,409 +160,411 @@ def test_nullVelocity_vec_m6():
assert np.allclose(scal_init2, scal_d.data[2])
-# def test_nullVelocity_m8():
-# """Basic test with random velocity. Using M4prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar')
-# scal_ref = Field(domain=box, name='Scalar_ref')
-# velo = Field(domain=box, name='Velocity',
-# formula=lambda x, y, z: (0., 0., 0.), isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M8'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: M8Prime,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-
-# advec.apply(Simulation(0., 0.1, 1))
-# advec_py.apply(Simulation(0., 0.1, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
-
-
-# def test_nullVelocity_vec_m8():
-# """Basic test with random velocity and vector field. Using M4prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar', isVector=True)
-# scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
-# velo = Field(domain=box, name='Velocity',
-# formula=lambda x, y, z: (0., 0., 0.), isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M8'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: M8Prime,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# scal_ref_d.data[1][...] = scal_d.data[1][...]
-# scal_ref_d.data[2][...] = scal_d.data[2][...]
-
-# advec.apply(Simulation(0., 0.075, 1))
-# advec_py.apply(Simulation(0., 0.075, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
-# assert np.allclose(scal_ref_d.data[1], scal_d.data[1])
-# assert np.allclose(scal_ref_d.data[2], scal_d.data[2])
-
-
-# def _randomVelocity_m4():
-# """Basic test with random velocity. Using M4prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar')
-# scal_ref = Field(domain=box, name='Scalar_ref')
-# velo = Field(domain=box, name='Velocity', isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M4'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: L2_1,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-# velo_d = velo.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# velo_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
-# velo_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-# velo_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-
-# advec.apply(Simulation(0., 0.1, 1))
-# advec_py.apply(Simulation(0., 0.1, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
-
-
-# def _randomVelocity_vec_m4():
-# """Basic test with random velocity vector Field. Using M4prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar', isVector=True)
-# scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
-# velo = Field(domain=box, name='Velocity', isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M4'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: L2_1,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-# velo_d = velo.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# scal_ref_d.data[1][...] = scal_d.data[1][...]
-# scal_ref_d.data[2][...] = scal_d.data[2][...]
-# velo_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
-# velo_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-# velo_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-
-# advec.apply(Simulation(0., 0.1, 1))
-# advec_py.apply(Simulation(0., 0.1, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
-# assert np.allclose(scal_ref_d.data[1], scal_d.data[1])
-# assert np.allclose(scal_ref_d.data[2], scal_d.data[2])
-
-
-# def test_randomVelocity_m6():
-# """Basic test with random velocity. Using M6prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar')
-# scal_ref = Field(domain=box, name='Scalar_ref')
-# velo = Field(domain=box, name='Velocity', isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M6'}
-# )
-# methods = {TimeIntegrator: RK2, Interpolation: Linear,
-# Remesh: L4_2, Support: 'CPU'}
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: L4_2,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-# velo_d = velo.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# velo_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
-# velo_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-# velo_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-
-# advec.apply(Simulation(0., 0.1, 1))
-# advec_py.apply(Simulation(0., 0.1, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
-
-
-# def test_randomVelocity_vec_m6():
-# """Basic test with random velocity vector Field. Using M6prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar', isVector=True)
-# scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
-# velo = Field(domain=box, name='Velocity', isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M6'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: L4_2,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-# velo_d = velo.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# scal_ref_d.data[1][...] = scal_d.data[1][...]
-# scal_ref_d.data[2][...] = scal_d.data[2][...]
-# velo_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
-# velo_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-# velo_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-
-# advec.apply(Simulation(0., 0.1, 1))
-# advec_py.apply(Simulation(0., 0.1, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
-# assert np.allclose(scal_ref_d.data[1], scal_d.data[1])
-# assert np.allclose(scal_ref_d.data[2], scal_d.data[2])
-
-
-# def test_randomVelocity_m8():
-# """Basic test with random velocity. Using M8prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar')
-# scal_ref = Field(domain=box, name='Scalar_ref')
-# velo = Field(domain=box, name='Velocity', isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M8'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: M8Prime,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-# velo_d = velo.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# velo_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
-# velo_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-# velo_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-
-# advec.apply(Simulation(0., 0.1, 1))
-# advec_py.apply(Simulation(0., 0.1, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0], atol=1e-07)
-
-
-# def test_randomVelocity_vec_m8():
-# """Basic test with random velocity vector Field. Using M8prime"""
-# box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
-# scal = Field(domain=box, name='Scalar', isVector=True)
-# scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
-# velo = Field(domain=box, name='Velocity', isVector=True)
-# advec = Advection(velo, scal,
-# resolutions={velo: [17, 17, 17],
-# scal: [17, 17, 17]},
-# method={Scales: 'p_M8'}
-# )
-# advec_py = Advection(velo, scal_ref,
-# resolutions={velo: [17, 17, 17],
-# scal_ref: [17, 17, 17]},
-# method={TimeIntegrator: RK2,
-# Interpolation: Linear,
-# Remesh: M8Prime,
-# Support: ''}
-# )
-# advec.discretize()
-# advec_py.discretize()
-# advec.setUp()
-# advec_py.setUp()
-
-# scal_d = scal.discreteFields.values()[0]
-# scal_ref_d = scal_ref.discreteFields.values()[0]
-# velo_d = velo.discreteFields.values()[0]
-
-# scal_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER)
-# scal_ref_d.data[0][...] = scal_d.data[0][...]
-# scal_ref_d.data[1][...] = scal_d.data[1][...]
-# scal_ref_d.data[2][...] = scal_d.data[2][...]
-# velo_d.data[0][...] = np.asarray(
-# np.random.random(scal_d.data[0].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
-# velo_d.data[1][...] = np.asarray(
-# np.random.random(scal_d.data[1].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-# velo_d.data[2][...] = np.asarray(
-# np.random.random(scal_d.data[2].shape),
-# dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
-
-# advec.apply(Simulation(0., 0.01, 1))
-# advec_py.apply(Simulation(0., 0.01, 1))
-
-# assert np.allclose(scal_ref_d.data[0], scal_d.data[0], atol=1e-07)
-# assert np.allclose(scal_ref_d.data[1], scal_d.data[1], atol=1e-07)
-# assert np.allclose(scal_ref_d.data[2], scal_d.data[2], atol=1e-07)
+def test_nullVelocity_m8():
+ """Basic test with random velocity. Using M4prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar')
+ scal_ref = Field(domain=box, name='Scalar_ref')
+ velo = Field(domain=box, name='Velocity',
+ formula=lambda x, y, z: (0., 0., 0.), isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M8'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: M8Prime,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+
+ advec.apply(Simulation(0., 0.1, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
+
+
+def test_nullVelocity_vec_m8():
+ """Basic test with random velocity and vector field. Using M4prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar', isVector=True)
+ scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
+ velo = Field(domain=box, name='Velocity',
+ formula=lambda x, y, z: (0., 0., 0.), isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M8'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: M8Prime,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ scal_ref_d.data[1][...] = scal_d.data[1][...]
+ scal_ref_d.data[2][...] = scal_d.data[2][...]
+
+ advec.apply(Simulation(0., 0.075, 1))
+ advec_py.apply(Simulation(0., 0.075, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
+ assert np.allclose(scal_ref_d.data[1], scal_d.data[1])
+ assert np.allclose(scal_ref_d.data[2], scal_d.data[2])
+
+
+def _randomVelocity_m4():
+ """Basic test with random velocity. Using M4prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar')
+ scal_ref = Field(domain=box, name='Scalar_ref')
+ velo = Field(domain=box, name='Velocity', isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M4'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L2_1,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+ velo_d = velo.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ velo_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
+ velo_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+ velo_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+
+ advec.apply(Simulation(0., 0.1, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
+
+
+def _randomVelocity_vec_m4():
+ """Basic test with random velocity vector Field. Using M4prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar', isVector=True)
+ scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
+ velo = Field(domain=box, name='Velocity', isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M4'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L2_1,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+ velo_d = velo.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ scal_ref_d.data[1][...] = scal_d.data[1][...]
+ scal_ref_d.data[2][...] = scal_d.data[2][...]
+ velo_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
+ velo_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+ velo_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+
+ advec.apply(Simulation(0., 0.1, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
+ assert np.allclose(scal_ref_d.data[1], scal_d.data[1])
+ assert np.allclose(scal_ref_d.data[2], scal_d.data[2])
+
+
+def test_randomVelocity_m6():
+ """Basic test with random velocity. Using M6prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar')
+ scal_ref = Field(domain=box, name='Scalar_ref')
+ velo = Field(domain=box, name='Velocity', isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M6'}
+ )
+ methods = {TimeIntegrator: RK2, Interpolation: Linear,
+ Remesh: L4_2, Support: 'CPU'}
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L4_2,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+ velo_d = velo.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ velo_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
+ velo_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+ velo_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+
+ advec.apply(Simulation(0., 0.1, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
+
+
+def test_randomVelocity_vec_m6():
+ """Basic test with random velocity vector Field. Using M6prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar', isVector=True)
+ scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
+ velo = Field(domain=box, name='Velocity', isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M6'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: L4_2,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+ velo_d = velo.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ scal_ref_d.data[1][...] = scal_d.data[1][...]
+ scal_ref_d.data[2][...] = scal_d.data[2][...]
+ velo_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
+ velo_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+ velo_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+
+ advec.apply(Simulation(0., 0.1, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0])
+ assert np.allclose(scal_ref_d.data[1], scal_d.data[1])
+ assert np.allclose(scal_ref_d.data[2], scal_d.data[2])
+
+
+def test_randomVelocity_m8():
+ """Basic test with random velocity. Using M8prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar')
+ scal_ref = Field(domain=box, name='Scalar_ref')
+ velo = Field(domain=box, name='Velocity', isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M8'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: M8Prime,
+ Support: ''}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+ velo_d = velo.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ velo_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
+ velo_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+ velo_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+
+ advec.apply(Simulation(0., 0.1, 1))
+ advec_py.apply(Simulation(0., 0.1, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0], atol=1e-07)
+
+
+def test_randomVelocity_vec_m8():
+ """Basic test with random velocity vector Field. Using M8prime"""
+ box = Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+ scal = Field(domain=box, name='Scalar', isVector=True)
+ scal_ref = Field(domain=box, name='Scalar_ref', isVector=True)
+ velo = Field(domain=box, name='Velocity', isVector=True)
+ advec = Advection(velo, scal,
+ resolutions={velo: [17, 17, 17],
+ scal: [17, 17, 17]},
+ method={Scales: 'p_M8'}
+ )
+ advec_py = Advection(velo, scal_ref,
+ resolutions={velo: [17, 17, 17],
+ scal_ref: [17, 17, 17]},
+ method={TimeIntegrator: RK2,
+ Interpolation: Linear,
+ Remesh: M8Prime,
+ Support: '',
+ Splitting: 'o2_FullHalf'}
+ )
+ advec.discretize()
+ advec_py.discretize()
+ advec.setUp()
+ advec_py.setUp()
+
+ scal_d = scal.discreteFields.values()[0]
+ scal_ref_d = scal_ref.discreteFields.values()[0]
+ velo_d = velo.discreteFields.values()[0]
+
+ scal_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER)
+ scal_ref_d.data[0][...] = scal_d.data[0][...]
+ scal_ref_d.data[1][...] = scal_d.data[1][...]
+ scal_ref_d.data[2][...] = scal_d.data[2][...]
+ velo_d.data[0][...] = np.asarray(
+ np.random.random(scal_d.data[0].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[0])
+ velo_d.data[1][...] = np.asarray(
+ np.random.random(scal_d.data[1].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+ velo_d.data[2][...] = np.asarray(
+ np.random.random(scal_d.data[2].shape),
+ dtype=PARMES_REAL, order=ORDER) / (2. * scal_d.resolution[1])
+
+ advec.apply(Simulation(0., 0.01, 1))
+ advec_py.apply(Simulation(0., 0.01, 1))
+
+ assert np.allclose(scal_ref_d.data[0], scal_d.data[0], atol=1e-07)
+ assert np.allclose(scal_ref_d.data[1], scal_d.data[1], atol=1e-07)
+ assert np.allclose(scal_ref_d.data[2], scal_d.data[2], atol=1e-07)
+