diff --git a/HySoP/CMake/FindFFTW.cmake b/CMake/FindFFTW.cmake
similarity index 100%
rename from HySoP/CMake/FindFFTW.cmake
rename to CMake/FindFFTW.cmake
diff --git a/HySoP/CMake/FindPythonFull.cmake b/CMake/FindPythonFull.cmake
similarity index 100%
rename from HySoP/CMake/FindPythonFull.cmake
rename to CMake/FindPythonFull.cmake
diff --git a/HySoP/CMake/FindPythonModule.cmake b/CMake/FindPythonModule.cmake
similarity index 100%
rename from HySoP/CMake/FindPythonModule.cmake
rename to CMake/FindPythonModule.cmake
diff --git a/HySoP/CMake/HySoPInstallSetup.cmake b/CMake/HySoPInstallSetup.cmake
similarity index 100%
rename from HySoP/CMake/HySoPInstallSetup.cmake
rename to CMake/HySoPInstallSetup.cmake
diff --git a/HySoP/CMake/HySoPTests.cmake b/CMake/HySoPTests.cmake
similarity index 100%
rename from HySoP/CMake/HySoPTests.cmake
rename to CMake/HySoPTests.cmake
diff --git a/HySoP/CMake/LibFindMacros.cmake b/CMake/LibFindMacros.cmake
similarity index 100%
rename from HySoP/CMake/LibFindMacros.cmake
rename to CMake/LibFindMacros.cmake
diff --git a/HySoP/CMake/MyTools.cmake b/CMake/MyTools.cmake
similarity index 100%
rename from HySoP/CMake/MyTools.cmake
rename to CMake/MyTools.cmake
diff --git a/HySoP/CMake/OutOfSourceBuild.cmake b/CMake/OutOfSourceBuild.cmake
similarity index 100%
rename from HySoP/CMake/OutOfSourceBuild.cmake
rename to CMake/OutOfSourceBuild.cmake
diff --git a/HySoP/CMake/PythonSetup.cmake b/CMake/PythonSetup.cmake
similarity index 100%
rename from HySoP/CMake/PythonSetup.cmake
rename to CMake/PythonSetup.cmake
diff --git a/HySoP/CMake/explore_python_config.py b/CMake/explore_python_config.py
similarity index 100%
rename from HySoP/CMake/explore_python_config.py
rename to CMake/explore_python_config.py
diff --git a/Docs/CMake b/Docs/CMake
new file mode 120000
index 0000000000000000000000000000000000000000..2a271a9f072593cf0f7de2713322f58cdc41f6e5
--- /dev/null
+++ b/Docs/CMake
@@ -0,0 +1 @@
+../CMake
\ No newline at end of file
diff --git a/Docs/sphinx_sources/user_guide/domains.rst b/Docs/sphinx_sources/user_guide/domains.rst
index 165b7e4573711fb91cb3e21c32775f7793599434..72404866b06cd7ad304593d7b4fb3405002a6938 100644
--- a/Docs/sphinx_sources/user_guide/domains.rst
+++ b/Docs/sphinx_sources/user_guide/domains.rst
@@ -10,10 +10,16 @@ Box-shaped domain
:class:`~box.Box`
+Data distribution
+-----------------
+.. toctree::
+ :maxdepth: 2
+
+ topologies
+
How to define a subset
----------------------
-
.. toctree::
:maxdepth: 2
diff --git a/Docs/sphinx_sources/user_guide/penalisation.rst b/Docs/sphinx_sources/user_guide/penalisation.rst
index 89a3f85998a1487966edf5d7056ffce80f6689a8..7bd0ba197edff83269f3271e076040217329eadb 100644
--- a/Docs/sphinx_sources/user_guide/penalisation.rst
+++ b/Docs/sphinx_sources/user_guide/penalisation.rst
@@ -3,12 +3,14 @@
Penalisation
============
+.. currentmodule:: hysop.operator.penalization
+
:math:`u(X,t)` is the current velocity field, :math:`\omega(X,t)` is the current vorticity field, :math:`u_D` the velocity of the obstacle, dt the time step.
Velocity penalisation
---------------------
-:class:`hysop.operator.penalization.Penalization`
+:class:`Penalization`
Solves:
@@ -41,7 +43,7 @@ After that, vorticity is computing as:
Vorticity penalisation
----------------------
-:class:`hysop.operator.penalize_vorticity.PenalizeVorticity`
+:class:`PenalizeVorticity`
Solves:
diff --git a/Docs/sphinx_sources/user_guide/subsets.rst b/Docs/sphinx_sources/user_guide/subsets.rst
index 6e07ab2f3017a84a7ff70015a6b7edd573d2fb51..255b8dd846f80d01a0c8fedbb3c5d1dc9b65c308 100644
--- a/Docs/sphinx_sources/user_guide/subsets.rst
+++ b/Docs/sphinx_sources/user_guide/subsets.rst
@@ -3,9 +3,9 @@
Subsets and obstacles
---------------------
-.. currentmodule:: hysop.domain.subsets
+.. currentmodule:: hysop.domain
-The class :class:`~subset.Subset` and its heirs allow user to define some subsets
+The class :class:`subsets.Subset` and its heirs allow user to define some subsets
of a predefined domain. This may be useful for:
* define specific area(s) where an operator applies (penalisation for example)
* reduce the size of i/o by writing hdf output only for a specific area.
@@ -15,12 +15,12 @@ of a predefined domain. This may be useful for:
Porous obstacles
^^^^^^^^^^^^^^^^
-:class:`~porous.BiPole`
+:class:`porous.BiPole`
.. image:: /figures/PolesExample.png
:width: 400pt
-:class:`~porous.QuadriPole`
+:class:`porous.QuadriPole`
.. image:: /figures/QuadriPoleExample.png
:width: 400pt
diff --git a/Docs/sphinx_sources/user_guide/topologies.rst b/Docs/sphinx_sources/user_guide/topologies.rst
new file mode 100644
index 0000000000000000000000000000000000000000..de8d285fae6254f384e19870cdc61a51b1da783d
--- /dev/null
+++ b/Docs/sphinx_sources/user_guide/topologies.rst
@@ -0,0 +1,15 @@
+.. _topologies:
+
+MPI topologies and space discretisation
+=======================================
+
+.. currentmodule:: hysop.mpi
+
+:class:`topology.Cartesian` objects are used to described both the mpi grid layout
+and the space discretisations (global and local to each mpi process).
+
+TODO :
+
+* describe the different ways to build topologies
+* describe the way we define local and global grids.
+
diff --git a/HySoP/CMake b/HySoP/CMake
new file mode 120000
index 0000000000000000000000000000000000000000..6126c513081cd9c08d0031ea93410931a228d572
--- /dev/null
+++ b/HySoP/CMake
@@ -0,0 +1 @@
+../CMake/
\ No newline at end of file
diff --git a/HySoP/CMakeLists.txt b/HySoP/CMakeLists.txt
index c5a4075d6bfdbbf1f88a9de2b72074275b87ac6f..e8fa1c138c9cbdd54671a5c4fd1788dfe3fd8f82 100644
--- a/HySoP/CMakeLists.txt
+++ b/HySoP/CMakeLists.txt
@@ -85,14 +85,14 @@ include(FindPythonModule)
# - python packages -
find_python_module(numpy REQUIRED)
-find_python_module(scipy REQUIRED)
+find_python_module(scipy)
find_python_module(matplotlib)
if(NOT matplotlib_FOUND)
find_python_module(Gnuplot REQUIRED)
endif()
find_python_module(scitools)
find_python_module(h5py REQUIRED)
-find_python_module(sympy REQUIRED)
+find_python_module(sympy)
# --- OpenCL ---
find_python_module(pyopencl REQUIRED)
# --- MPI ---
diff --git a/HySoP/hysop/domain/__init__.py b/HySoP/hysop/domain/__init__.py
index 79b4773973e2ba2129afceeea2ba982ef6edc34f..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/HySoP/hysop/domain/__init__.py
+++ b/HySoP/hysop/domain/__init__.py
@@ -1,5 +0,0 @@
-## @package hysop.domain
-# Physical domains descriptions and related tools.
-# At the time :
-# - Box for box-shaped domains
-# - obstacle
diff --git a/HySoP/hysop/domain/box.py b/HySoP/hysop/domain/box.py
index b67c20794e30d513b349af7b91476146c539b62a..cbf6317fc5d6b267bb700bc9db88947335aaac6f 100644
--- a/HySoP/hysop/domain/box.py
+++ b/HySoP/hysop/domain/box.py
@@ -1,7 +1,4 @@
-"""
-@file box.py
-Box-shaped domains definition.
-
+"""Box-shaped domains definition.
"""
from hysop.domain.domain import Domain
from hysop.constants import PERIODIC, debug
@@ -9,10 +6,9 @@ import hysop.tools.numpywrappers as npw
class Box(Domain):
- """
- Box-shaped domain description.
- @note FP : BC as parameter may be better?
- @todo Have different BC
+ """ Box-shaped domain description.
+
+ todo implement different boundary conditions typesKHave different BC
"""
@debug
@@ -20,8 +16,12 @@ class Box(Domain):
"""
Create a Periodic Box from a dimension, length and origin.
- @param length : Box length. Default [1.0, ...]
- @param origin : Box minimum position. Default [0., ...]
+ Parameters
+ -----------
+ length : list or numpy array of double
+ box sides lengthes. Default = [1.0, ...]
+ origin : list or numpy array of double
+ position of the lowest point of the box. Default [0., ...]
Example:
@@ -59,10 +59,6 @@ class Box(Domain):
assert list(self.boundaries).count(PERIODIC) == self.dimension, msg
def __str__(self):
- """
- Informations display.
- @return Informations
- """
s = str(self.dimension) + \
"D box (parallelepipedic or rectangular) domain : \n"
@@ -74,5 +70,5 @@ class Box(Domain):
c1 = (self.length == other.length).all()
c2 = (self.origin == other.origin).all()
c3 = (self.boundaries == other.boundaries).all()
- c4 = self.currentTask() == other.currentTask()
+ c4 = self.current_task() == other.current_task()
return c1 and c2 and c3 and c4
diff --git a/HySoP/hysop/domain/control_box.py b/HySoP/hysop/domain/control_box.py
index 2c3bd2214030c074db8de62de079d1017820892b..47d285e5a0fc0d522ceaaf583a6302fa13003f59 100644
--- a/HySoP/hysop/domain/control_box.py
+++ b/HySoP/hysop/domain/control_box.py
@@ -12,6 +12,9 @@ class ControlBox(SubBox):
Useful to define control volume to perform integration.
"""
def __init__(self, **kwds):
+ """
+ Same parameters as for :class:`hysop.domain.subsets.SubBox`.
+ """
super(ControlBox, self).__init__(**kwds)
# We must have a real box, not a plane ...
@@ -20,6 +23,9 @@ class ControlBox(SubBox):
self.surf = [None] * len(self.length) * 2
def discretize(self, topo):
+ """Create a sub meshes for topo inside the control box and
+ on its faces.
+ """
# Create the mesh for the whole box
super(ControlBox, self).discretize(topo)
# Create a mesh for each side
@@ -54,17 +60,26 @@ class ControlBox(SubBox):
def integrate_on_faces(self, field, topo, list_dir,
component=0, root=None):
- """
- @param[in] field : a hysop continuous field
- @param[in] topo : a hysop.mpi.topology.Cartesian topology
- @param[in] root : root process used for mpi reduction. If None
- reduction is done on all processes from topo.
- @param[in] list_dir : list of surfaces on which we must integrate.
- 0 : normal to xdir, lower surf,
- 1 : normal to xdir, upper surf (in x dir)
- 2 : normal to ydir, lower surf, and so on.
- the field must be integrated.
- @return a numpy array, with res = sum
+ """Integration of a field on one or more faces of the box
+
+ Parameters
+ ----------
+ field : :class:`~hysop.fields.continuous.Field`
+ a field to be integrated on the box
+ topo : :class:`~hysop.mpi.topology.Cartesian`
+ set mesh/topology used for integration
+ list_dir : list of int
+ indices of faces on which integration is required
+ 0 : normal to xdir, lower surf,
+ 1 : normal to xdir, upper surf (in x dir)
+ 2 : normal to ydir, lower surf, and so on.
+ component : int, optional
+ number of the field component to be integrated
+ root : int, optional
+ rank of the leading mpi process (to collect data)
+ If None reduction is done on all processes from topo.
+
+ Returns a numpy array, with res = sum
of the integrals of a component of field on all surf in list_dir
"""
res = 0.
@@ -79,17 +94,24 @@ class ControlBox(SubBox):
return topo.comm.reduce(res, root=root)
def integrate_on_faces_allc(self, field, topo, list_dir, root=None):
- """
- @param[in] field : a hysop continuous field
- @param[in] topo : a hysop.mpi.topology.Cartesian topology
- @param[in] root : root process used for mpi reduction. If None
- reduction is done on all processes from topo.
- @param[in] list_dir : list of surfaces on which we must integrate.
- 0 : normal to xdir, lower surf,
- 1 : normal to xdir, upper surf (in x dir)
- 2 : normal to ydir, lower surf, and so on.
- the field must be integrated.
- @return a numpy array, with res[i] = sum
+ """Integration of a field on one or more faces of the box
+
+ Parameters
+ ----------
+ field : :class:`~hysop.fields.continuous.Field`
+ a field to be integrated on the box
+ topo : :class:`~hysop.mpi.topology.Cartesian`
+ set mesh/topology used for integration
+ list_dir : list of int
+ indices of faces on which integration is required
+ 0 : normal to xdir, lower surf,
+ 1 : normal to xdir, upper surf (in x dir)
+ 2 : normal to ydir, lower surf, and so on.
+ root : int, optional
+ rank of the leading mpi process (to collect data)
+ If None reduction is done on all processes from topo.
+
+ Returns a numpy array, with res[i] = sum
of the integrals of component i of field on all surf in list_dir
"""
nbc = field.nb_components
diff --git a/HySoP/hysop/domain/domain.py b/HySoP/hysop/domain/domain.py
index cecad94d3803abfac75443f1e9110db504e36697..6868c724dc1f8135af933fe003296014985120cb 100644
--- a/HySoP/hysop/domain/domain.py
+++ b/HySoP/hysop/domain/domain.py
@@ -1,7 +1,4 @@
-"""
-@file domain.py
-
-Abstract interface for physical domains description.
+"""Abstract interface for physical domains description.
"""
from abc import ABCMeta, abstractmethod
@@ -24,21 +21,51 @@ class Domain(object):
@debug
@abstractmethod
def __init__(self, dimension=3, proc_tasks=None, bc=None):
- """ Constructor
- @param dimension integer : domain dimension.
- @param proc_tasks : task - mpi processes connection.
+ """
+
+ Parameters
+ ----------
+ dimension : integer, optional
+ proc_tasks : list or array of int, optional
+ connectivity between mpi process id and task number.
+ See notes below. Default : all procs on task DEFAULT_TASK_ID.
+ bc : list or array of int, optional
+ type of boundary condtions. See notes below.
+
+
+ Attributes
+ ----------
+
+ dimension : int
+ topologies : dictionnary of :class:`~hysop.mpi.topology.Cartesian`
+ all topologies already built for this domain.
+ comm_task : mpi communicator
+ corresponds to the task that owns the current process.
+ boundaries : list or array of int
+ boundary conditions type.
+
+ Notes
+ -----
+ About boundary conditions
+ ==========================
+ At the time, only periodic boundary conditions are implemented.
+ Do not use this parameter (i.e. let default values)
+
+ About MPI Tasks
+ ================
proc_tasks[n] = 12 means that task 12 owns proc n.
+
"""
- ## Domain dimension.
+ # Domain dimension.
self.dimension = dimension
- ## A list of all the topologies defined on this domain.
- ## Each topology is unique in the sense defined by
- ## the comparison operator in the class hysop.mpi.topology.Cartesian.
+ # A list of all the topologies defined on this domain.
+ # Each topology is unique in the sense defined by
+ # the comparison operator in the class hysop.mpi.topology.Cartesian.
self.topologies = {}
- ## Connectivity between mpi tasks and proc numbers :
- ## self.task_of_proc[i] returns the task which is bound to proc of
- ## rank i in main_comm.
+ # Connectivity between mpi tasks and proc numbers :
+ # self.task_of_proc[i] returns the task which is bound to proc of
+ # rank i in main_comm.
# Warning FP: rank in main_comm! We consider that each proc has
# one and only one task.
# Maybe we should change this and allow proc_tasks in subcomm
@@ -53,19 +80,29 @@ class Domain(object):
comm_s = main_comm.Split(color=proc_tasks[main_rank],
key=main_rank)
- ## The sub-communicator corresponding to the task that owns
- ## the current process.
+ # The sub-communicator corresponding to the task that owns
+ # the current process.
self.comm_task = comm_s
- ## Boundary conditions type.
- ## At the moment, only periodic conditions are allowed.
+ # Boundary conditions type.
+ # At the moment, only periodic conditions are allowed.
self.boundaries = bc
- def isOnTask(self, params):
- """
- @param params : a hysop.mpi.MPIParams object
+ def is_on_task(self, params):
+ """Test if the current process corresponds to param task.
+
+ :param params: :class:`~hysop.mpi.MPIParams` or int
+ description of the mpi setup or task number.
or an int (task number)
- @return true if params.task_id (or params if int) corresponds with
- task_id of the current proc.
+ :return bool:
+
+ Example :
+ :code:
+ dom.is_on_task(4)
+ mpi_params = MPIParams(task_id=4)
+ dom.is_on_task(mpi_params)
+ :code:
+
+ returns True if the current process runs on task 4.
"""
if params.__class__ is MPIParams:
return params.task_id == self._tasks_on_proc[main_rank]
@@ -73,147 +110,151 @@ class Domain(object):
return params == self._tasks_on_proc[main_rank]
def tasks_on_proc(self, index):
- """
- @param index: proc number
- @return task number of mpi proc 'index'
+ """Get task number for a given mpi process
+
+ :param index: int
+ proc number
+ :return int:
"""
return self._tasks_on_proc[index]
def tasks_list(self):
- """
- @return the connectivity between proc. numbers and their associated
- task (numbers in main_comm).
+ """Get task/process number connectivity
+
+ :return list or array of int:
+ such that res[proc_number] = task_id
"""
return self._tasks_on_proc
- def currentTask(self):
- """
- @return task number of the current proc.
+ def current_task(self):
+ """Get task number of the current mpi process.
+
+ :return int:
"""
return self._tasks_on_proc[main_rank]
def create_topology(self, discretization, dim=None, mpi_params=None,
shape=None, cutdir=None):
- """
- Create or return an existing hysop.mpi.topology.
-
- Either it gets the topology corresponding to the input arguments
- if it exists (in the sense of the comparison operator defined in
- hysop.mpi.topology.Cartesian)
- or it creates a new topology and register it i n the topology list.
- @param domain : the geometry; it must be a box.
- @param discretization : a hysop.tools.parameters.Discretization
- with:
- - resolution = Number of points in the domain
- in each direction. We assume that first point corresponds
- to origin, and last point to boundary point,
- whatever the boundary type is.
- That is x[0] = domain.origin and
- x[Discretization.resolution-1] = domain.Lengths_x.
- - ghosts = number of points in the ghost layer
- @param dim : dimension of the topology
- @param mpi_params : a hysop.tools.parameters.MPIParams, with:
- - comm : MPI communicator used to create this topology
- (default = main_comm)
- - task_id : id of the task that owns this topology.
- @param isperiodic : periodicity of the topology (in each direction)
- @param cutdir : array of bool, set cutdir[dir] = True if you want
- to distribute data through direction dir.
- @param shape : topology resolution
- (i.e process layout in each direction).
- @return the required topology.
+ """Create or return an existing :class:`~hysop.mpi.topology.Cartesian`.
+
+ Parameters
+ -----------
+ discretization : :class:`~hysop.tools.parameters.Discretization`
+ description of the global space discretization
+ (resolution and ghosts).
+ dim : int, optional
+ mpi topology dimension.
+ mpi_params : :class:`~hysop.tools.parameters.MPIParams`, optional
+ mpi setup (comm, task ...).
+ If None, comm = main_comm, task = DEFAULT_TASK_ID.
+ shape : list or array of int
+ mpi grid layout
+ cutdir : list or array of bool
+ set which directions must (may) be distributed,
+ cutdir[dir] = True to distribute data along dir.
+
+ Returns
+ -------
+ :class:`~hysop.mpi.topology.Cartesian`
+ Either it gets the topology corresponding to the input arguments
+ if it exists (in the sense of the comparison operator defined in
+ :class:`~hysop.mpi.topology.Cartesian`)
+ or it creates a new topology and register it in the topology dict.
+
"""
# set task number
- tid = self.currentTask()
+ tid = self.current_task()
if mpi_params is None:
mpi_params = MPIParams(comm=self.comm_task, task_id=tid)
else:
msg = 'Try to create a topology on a process that does not'
msg += ' belong to the current task.'
assert mpi_params.task_id == tid, msg
- newTopo = Cartesian(self, discretization, dim=dim,
- mpi_params=mpi_params, shape=shape,
- cutdir=cutdir)
- newid = newTopo.get_id()
+ new_topo = Cartesian(self, discretization, dim=dim,
+ mpi_params=mpi_params, shape=shape,
+ cutdir=cutdir)
+ newid = new_topo.get_id()
return self.topologies[newid]
def create_plane_topology_from_mesh(self, localres, global_start,
cdir=None, **kwds):
- """
- Create or return an existing hysop.mpi.topology.
+ """Create 1D topology from predifined :class:`~hysop.domain.mesh.Mesh`.
Define a 'plane' (1D) topology for a given mesh resolution.
This function is to be used when topo/discretization features
come from an external routine (e.g. scales or fftw)
- @param localres : local mesh resolution
- @param global_start : global indices of the lowest point
- of the local mesh.
- @param cdir : direction of cutting (i.e. normal to mpi planes)
- default = last if fortran order, first if C order.
+
+ Parameters
+ ----------
+ localres : list or array of int
+ local mesh resolution
+ global_start : list or array of int
+ indices in the global mesh of the lowest point of the local mesh.
+ cdir : int, optional
+ direction where data must be distributed in priority.
+ Default = last if fortran order, first if C order.
+ **kwds : other args
+ used in :class:`~hysop.mpi.topology.Cartesian` build.
+
"""
- tid = self.currentTask()
+ tid = self.current_task()
if 'mpi_params' not in kwds:
kwds['mpi_params'] = MPIParams(comm=self.comm_task, task_id=tid)
else:
msg = 'Try to create a topology on a process that does not'
msg += ' belong to the current task.'
assert kwds['mpi_params'].task_id == tid, msg
- newTopo = Cartesian.plane_precomputed(localres, global_start, cdir,
- domain=self, **kwds)
+ new_topo = Cartesian.plane_precomputed(localres, global_start, cdir,
+ domain=self, **kwds)
- newid = newTopo.get_id()
+ newid = new_topo.get_id()
return self.topologies[newid]
- def _checkTopo(self, newTopo):
- """
- Return the id of the input topology
+ def _check_topo(self, new_topo):
+ """Returns the id of the input topology
if it exists in the domain list, else return -1.
- @param newTopo : the topology to check
- @return id of the topology if present else -1.
"""
otherid = -1
for top in self.topologies.values():
- if newTopo == top:
+ if new_topo == top:
otherid = top.get_id()
break
return otherid
- def register(self, newTopo):
- """
- Add a new topology in the list of this domain.
+ def register(self, new_topo):
+ """Add a new topology in the list of this domain.
Do nothing if an equivalent topology is already
in the list.
- @param newTopo : the topology to be registered.
- @return the id of the new registered topology
+
+ :Returns the id of the new registered topology
or of the corresponding 'old one' if it already exists.
"""
- otherid = self._checkTopo(newTopo)
+ otherid = self._check_topo(new_topo)
if otherid < 0:
- self.topologies[newTopo.get_id()] = newTopo
- newTopo.isNew = True
- newid = newTopo.get_id()
+ self.topologies[new_topo.get_id()] = new_topo
+ new_topo.isNew = True
+ newid = new_topo.get_id()
else:
# No registration
- newTopo.isNew = False
+ new_topo.isNew = False
newid = otherid
return newid
def remove(self, topo):
- """
- Remove a topology from the list of this domain.
+ """Remove a topology from the list of this domain.
Do nothing if the topology does not exist in the list.
Warning : the object topo is not destroyed, only
removed from self.topologies.
- @param topo : the topology to be removed.
- @return either the od of the removed topology
+
+ Returns either the id of the removed topology
or -1 if nothing is done.
"""
- otherid = self._checkTopo(topo)
+ otherid = self._check_topo(topo)
if otherid >= 0:
self.topologies.pop(otherid)
return otherid
- def printTopologies(self):
+ def print_topologies(self):
"""
Print all topologies of the domain.
"""
@@ -238,7 +279,7 @@ class Domain(object):
def i_periodic_boundaries(self):
"""
- @return list of directions where
+ Returns the list of directions where
boundaries are periodic
"""
return np.where(self.boundaries == PERIODIC)[0]
diff --git a/HySoP/hysop/domain/mesh.py b/HySoP/hysop/domain/mesh.py
index 616be81b42508c5530b63aa8a313ea1fcfce2f25..b3f26498ae42e286c2f19f221d5d7b0c17e75508 100644
--- a/HySoP/hysop/domain/mesh.py
+++ b/HySoP/hysop/domain/mesh.py
@@ -1,5 +1,4 @@
-"""
-@file mesh.py local cartesian grid.
+"""local (mpi process) cartesian grid.
"""
from hysop.constants import debug
import hysop.tools.numpywrappers as npw
@@ -63,12 +62,16 @@ class Mesh(object):
@debug
def __init__(self, parent, discretization, resolution, global_start):
"""
- @param parent : the geometry on which this mesh is defined
- (it must be a rectangle or a rectangular cuboid)
- @param discretization : hysop.tools.parameters.Discretization which
- describes the global discretization of the domain
- (resolution and ghost layer)
- @param global_start
+ Parameters
+ -----------
+ parent : :class:`hysop.domain.domain.Domain`
+ the geometry on which this mesh is defined
+ (it must be a rectangle or a rectangular cuboid)
+ discretization : :class:`hysop.tools.parameters.Discretization`
+ describes the global discretization of the domain
+ (resolution and ghost layer)
+ global_start : list or array of int
+ indices in the global mesh of the lowest point of the local mesh.
"""
# Note : all variables with 'global' prefix are related
@@ -155,9 +158,11 @@ class Mesh(object):
/ self.space_step))
def is_inside(self, point):
- """
- @param point: coordinates of a point
- @return : true if the point is inside the local (mpi proc) mesh
+ """True if the point belongs to volume or surface
+ described by the mesh
+
+ :param point: list or numpy array,
+ coordinates of a point
"""
point = npw.asrealarray(point)
return ((point - self.origin) >= 0.).all() and ((self.end
@@ -175,13 +180,12 @@ class Mesh(object):
return s
def convert2local(self, sl):
- """
- Return the local mesh indices from
- their global indices.
- Ilocal = Iglobal - _global_start + ghost
- @param sl: the list of slices (one per dir) for global indices
+ """convert indices from global mesh to local one
+
+ :param sl: list of slices (one slice per dir) of global indices
something like [(start_dir1, end_dir1), ...]
- @return the same kind of list but in local coordinates
+
+ Returns the same kind of list but in local coordinates
"""
tmp = utils.intersl(sl, self.position)
if tmp is not None:
@@ -192,26 +196,22 @@ class Mesh(object):
return None
def convert2global(self, sl):
- """
- Return the global mesh indices from
- their local indices.
- Iglobal = Ilocal + _global_start - ghost
- @param loc_index an array of size domain.dim*2
- with loc_index = [start_dir1, end_dir1, start_dir2, end_dir2, ...]
- return an array of indices, same 'setup' as loc_index.
+ """convert indices from local mesh to global one
+
+ :param sl: list of slices (one slice per dir) of local indices
+ something like [(start_dir1, end_dir1), ...]
+
+ Returns the same kind of list but in global coordinates
"""
return [slice(sl[i].start - self._shift[i],
sl[i].stop - self._shift[i])
for i in xrange(self._dim)]
def __eq__(self, other):
- """
- = operator for mesh.
- We consider that 2 submeshes are equal if all the conditions
- below are fulfilled:
- - their global resolution is the same
- - the ghost layer is the same
- - their local resolution is the same
+ """True if meshes are equal.
+ Two meshes are equal if they have the same
+ global resolution, the same ghost layer and the same
+ local resolution is the same
"""
if self.__class__ != other.__class__:
return False
@@ -219,8 +219,7 @@ class Mesh(object):
npw.equal(self.resolution, other.resolution).all()
def start(self):
- """
- @return indices of the lowest point of the local mesh
+ """Returns indices of the lowest point of the local mesh
in the global grid
"""
return npw.asdimarray([p for p in
@@ -228,9 +227,9 @@ class Mesh(object):
for d in xrange(self._dim))])
def stop(self):
- """
- @return indices + 1 of the 'highest' point of the local mesh
+ """Returns indices + 1 of the 'highest' point of the local mesh
in the global grid.
+
Warning: made to be used in slices :
slice(mesh.start()[d], mesh.stop()[d]) represents the
points of the current process mesh, indices given
@@ -241,8 +240,7 @@ class Mesh(object):
for d in xrange(self._dim))])
def global_resolution(self):
- """
- @return the resolution of the global grid (on the whole
+ """Returns the resolution of the global grid (on the whole
domain, whatever the mpi distribution is).
Warning : this depends on the type of boundary conditions.
For periodic bc, this corresponds to the user-defined discretization
diff --git a/HySoP/hysop/domain/porous.py b/HySoP/hysop/domain/porous.py
index f57fe564374ede00bb0e3e43653a92c60131262f..5d7a5643b7da85fc99e6bc288aa037685d67dae7 100644
--- a/HySoP/hysop/domain/porous.py
+++ b/HySoP/hysop/domain/porous.py
@@ -1,4 +1,17 @@
-"""Pre-defined subsets as combination of basic geometries.
+"""Porous (multi-layers) subsets:
+
+.. currentmodule hysop.domain.porous
+* :class:`~Porous`
+* :class:`~BiPole`
+* :class:`~QuadriPole`
+* :class:`~RingPole`
+* :class:`~Ring`
+
+See also
+--------
+
+* :class:`~hysop.domain.subsets.Subset` for standard subsets.
+* :ref:`subsets` in HySoP user guide.
"""
from hysop.domain.subsets import Subset, Sphere, HemiSphere
@@ -19,7 +32,7 @@ class Porous(Subset):
"""
Parameters
-----------
- source : :class:`~hysop.domain.subsets.subset.Subset`
+ source : :class:`~hysop.domain.subsets.Subset`
geometry of the object. The only authorized sources are :
Cylinder, Sphere, HemiSphere, HemiCylinder
origin : list or array of coordinates
@@ -27,7 +40,17 @@ class Porous(Subset):
layers : list of real
width of each layer (starting from outside layer)
"""
- super(Porous, self).__init__(**kwds)
+ # the external layer of the porous subset.
+ # Note FP : we need this envelop to set a proper
+ # indicator function for the porous subset, since this
+ # function is used in all the union, intersection ...
+ # methods.
+ # is_inside will return true for any point inside this
+ # envelop.
+ envelop = source(parent=kwds['parent'], origin=origin,
+ radius=sum(layers))
+ assert 'chi' not in kwds
+ super(Porous, self).__init__(chi=envelop.get_chi(), **kwds)
# Center position
self.origin = npw.asrealarray(origin).copy()
self.is_porous = True
@@ -36,7 +59,6 @@ class Porous(Subset):
assert isinstance(layers, list)
# width of porous layers on the source
self._layers = layers
-
self._box = None
def discretize(self, topo):
@@ -117,24 +139,22 @@ class BiPole(Porous):
----------
poles_thickness : real or list or array of real numbers
thickness of layer(s) in direction given by poles_dir.
- Use :
- * poles_thickness = value or [value] if all poles have the
- same width.
- * poles_thickness = [v1, v2] to set bottom pole to v1 and top
- pole to v2
+ Use poles_thickness = value or [value]
+ if all poles have the same width.
+ Or, to set bottom pole to v1 and top pole to v2
+ poles_thickness = [v1, v2].
poles_dir : int
direction for which poles are 'computed'. Default = last dir.
kwds : args for base class
- Example:
- --------
+ Notes
+ -----
- >>> obst = BiPole(..., source=HemiSphere,
- ... poles_thickness=[0.1, 0.2], poles_dir=[2]
+ An example to create poles in direction 2, i.e. z,
+ of thickness 0.2 on top and 0.1 at the bottom of an HemiSphere::
- will create poles in direction 2, i.e. z,
- of thickness 0.2 on top and 0.1 at the
- bottom of an HemiSphere.
+ obst = BiPole(..., source=HemiSphere, poles_thickness=[0.1, 0.2],
+ poles_dir=[2])
"""
super(BiPole, self).__init__(**kwds)
@@ -189,24 +209,22 @@ class QuadriPole(Porous):
----------
poles_thickness : real or list or array of real numbers
thickness of layer(s) in direction given by poles_dir.
- Use :
- * poles_thickness = value or [v1, v2] if all poles in one direction
- have the
- same width.
- * poles_thickness = [v1, v2, v3, v4] to set
- bottom/top poles to v1/v2 in y dir
- and bottom/top poles to v3/v4 in z dir.
+ Use poles_thickness = value or [v1, v2]
+ if all poles in one direction have the same width.
+ Or, to set bottom/top poles to v1/v2 in y dir and
+ bottom/top poles to v3/v4 in z dir,
+ poles_thickness = [v1, v2, v3, v4].
+
kwds : args for base class
- Example:
- --------
+ Notes
+ -----
- >>> obst = QuadriPole(source=HemiSphere,
- ... poles_thickness=[0.1, 0.2])
+ An example to create poles in direction 2, i.e. z,
+ of thickness 0.2 on top and 0.1 at the bottom, on an HemiSphere::
- will create poles in direction 2, i.e. z,
- of thickness 0.2 on top and 0.1 at the
- bottom, on an HemiSphere.
+ obst = QuadriPole(..., source=HemiSphere,
+ poles_thickness=[0.1, 0.2])
"""
super(QuadriPole, self).__init__(**kwds)
diff --git a/HySoP/hysop/domain/subsets.py b/HySoP/hysop/domain/subsets.py
index 572620b8f2ae60aad420419403bfc66ce230a7ff..65baf2d022d0296f1a3c8250335b1a7c059b7739 100644
--- a/HySoP/hysop/domain/subsets.py
+++ b/HySoP/hysop/domain/subsets.py
@@ -1,4 +1,18 @@
-""" subsets of a given domain (sphere ...)
+"""Subsets of a given domain:
+
+* :class:`~hysop.domain.subsets.Sphere`,
+* :class:`~hysop.domain.subsets.HemiSphere`,
+* :class:`~hysop.domain.subsets.Cylinder`,
+* :class:`~hysop.domain.subsets.HemiCylinder`,
+* :class:`~hysop.domain.subsets.SubBox`,
+* :class:`~hysop.domain.subsets.Subset` (abstract base class).
+
+See also
+--------
+
+* :class:`~hysop.domain.porous.Porous` for porous (multi-layers) subsets.
+* :ref:`subsets` in HySoP user guide.
+
"""
from hysop.domain.domain import Domain
import numpy as np
@@ -13,24 +27,18 @@ class Subset(object):
"""
A subset is a geometry defined inside a domain.
Given a topology on the parent domain, the subset must
- provide 'slices' or 'tabs' to allow the computation of a discrete
- field on the subset, with:
+ provide some lists of indices to allow the computation of a discrete
+ field on the subset, with something like
- data[subset.slices] = ...
+ data[subset.ind[topo] = ...
or
data[subset.tab] = ...
- slices : a list of python slices defining the grid points inside the subset
- ind : a tuple of arrays of indices
-
There are two types of subsets:
- 'regular' ones, those on which a regular mesh can be defined
- others, like spheres, cylinders ...
- For regular subsets, one must use 'slices' to compute field inside
- the subset whereas 'tab' is required for others.
-
- Regular subsets have an attribute 'ind' which is a dictionnary of tuples
+ Subsets have an attribute 'ind' which is a dictionnary of tuples
representing the points inside the subset, (keys = topology) such that:
ind[topo] = (i_x, i_y, i_z),
with i_x[i], i_y[i], i_z[i] for each index i being the indices in the
@@ -61,10 +69,10 @@ class Subset(object):
-----
func argument list must start with space coordinates,
- e.g. for a 3D domain something like:
+ e.g. for a 3D domain something like::
- >>> def myfunc(x, y, z, ...):
- ...
+ def myfunc(x, y, z, ...):
+ ...
"""
assert isinstance(parent, Domain)
@@ -98,6 +106,22 @@ class Subset(object):
"""
return self._is_inside(*args) <= 0
+ def _set_chi(self, chi):
+ """Reset indicator function
+
+ Mostly for internal setup (during porous obst. init), so do not
+ use it or use it with care ...
+ """
+ self._is_inside = chi
+
+ def get_chi(self):
+ """Get indicator function
+
+ Mostly for internal setup (during porous obst. init), so do not
+ use it or use it with care ...
+ """
+ return self._is_inside
+
def discretize(self, topo):
"""
Create the list of indices for points inside the subset
@@ -144,7 +168,7 @@ class Subset(object):
Parameters
----------
- slist : a list of :class:`~hysop.subset.subsets.Subset`
+ slist : a list of :class:`~hysop.domain.subsets.Subset`
topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
@@ -162,7 +186,7 @@ class Subset(object):
Parameters
----------
- slist : a list of :class:`~hysop.subset.subsets.Subset`
+ slist : a list of :class:`~hysop.domain.subsets.Subset`
topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
@@ -193,7 +217,7 @@ class Subset(object):
Parameters
----------
- s1, s2 : lists of :class:`~hysop.subset.subsets.Subset`
+ s1, s2 : lists of :class:`~hysop.domain.subsets.Subset`
topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
@@ -211,8 +235,8 @@ class Subset(object):
Parameters
----------
- slist : list of Subset
- topo : `~hysop.mpi.topology.Cartesian`
+ slist : list of :class:`~hysop.domain.subsets.Subset`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
--------
@@ -228,8 +252,8 @@ class Subset(object):
Parameters
----------
- slist : list of Subset
- topo : `~hysop.mpi.topology.Cartesian`
+ slist : list of :class:`~hysop.domain.subsets.Subset`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
--------
@@ -260,8 +284,8 @@ class Subset(object):
Parameters
----------
- s1, s2 : Subset
- topo : `~hysop.mpi.topology.Cartesian`
+ s1, s2 : :class:`~hysop.domain.subsets.Subset`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
--------
@@ -275,8 +299,8 @@ class Subset(object):
Parameters
----------
- s1, s2 : Subset
- topo : `~hysop.mpi.topology.Cartesian`
+ s1, s2 : :class:`~hysop.domain.subsets.Subset`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
--------
@@ -298,8 +322,8 @@ class Subset(object):
Parameters
----------
- s1, s2 : list of Subset
- topo : `~hysop.mpi.topology.Cartesian`
+ s1, s2 : list of :class:`~hysop.domain.subsets.Subset`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
Returns
--------
@@ -314,9 +338,9 @@ class Subset(object):
Parameters
----------
- field : `~hysop.fields.continuous.Field`
+ field : :class:`~hysop.fields.continuous.Field`
a field to be integrated on the box
- topo : `~hysop.mpi.topology.Cartesian`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
set mesh/topology used for integration
root : int, optional
rank of the leading mpi process (to collect data)
@@ -343,9 +367,9 @@ class Subset(object):
Parameters
----------
- field : `~hysop.fields.continuous.Field`
+ field : :class:`~hysop.fields.continuous.Field`
a field to be integrated on the box
- topo : `~hysop.mpi.topology.Cartesian`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
set mesh/topology used for integration
component : int, optional
number of the field component to be integrated
@@ -369,9 +393,9 @@ class Subset(object):
Parameters
----------
- field : `~hysop.fields.continuous.Field`
+ field : :class:`~hysop.fields.continuous.Field`
a field to be integrated on the box
- topo : `~hysop.mpi.topology.Cartesian`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
set mesh/topology used for integration
component : int, optional
number of the field component to be integrated
@@ -395,7 +419,7 @@ class Subset(object):
Parameters
----------
- field : `~hysop.fields.discrete.DiscreteField`
+ field : :class:`~hysop.fields.discrete.DiscreteField`
a field to be integrated on the box
root : int, optional
rank of the leading mpi process (to collect data)
@@ -422,7 +446,7 @@ class Subset(object):
Parameters
----------
- field : `~hysop.fields.discrete.DiscreteField`
+ field : :class:`~hysop.fields.discrete.DiscreteField`
a field to be integrated on the box
component : int, optional
number of the field component to be integrated
@@ -446,7 +470,7 @@ class Subset(object):
Parameters
----------
- field : `~hysop.fields.discrete.DiscreteField`
+ field : :class:`~hysop.fields.discrete.DiscreteField`
a field to be integrated on the box
component : int, optional
number of the field component to be integrated
@@ -474,6 +498,8 @@ class Sphere(Subset):
"""
Parameters
----------
+ origin : :class:`~hysop.domain.subsets.Subset`
+
origin : list or array
position of the center
radius : double, optional
@@ -534,6 +560,13 @@ class HemiSphere(Sphere):
self.LeftBox = left_box
def chi(self, *args):
+ """Indicator function of the subset
+
+ Returns
+ -------
+ array of bool
+
+ """
return np.logical_and(self._is_inside(*args) <= 0,
self.LeftBox(args[self.cutdir]) <= 0)
@@ -583,6 +616,13 @@ class Cylinder(Subset):
length=length, origin=origin)
def chi(self, *args):
+ """Indicator function of the subset
+
+ Returns
+ -------
+ array of bool
+
+ """
return np.logical_and(self._is_inside(*args) <= 0,
args[self.axis] ==
args[self.axis])
@@ -620,6 +660,13 @@ class HemiCylinder(Cylinder):
self.LeftBox = left_box
def chi(self, *args):
+ """Indicator function of the subset
+
+ Returns
+ -------
+ array of bool
+
+ """
return (np.logical_and(
np.logical_and(self._is_inside(*args) <= 0,
self.LeftBox(args[self.cutdir]) <= 0),
@@ -682,7 +729,7 @@ class SubBox(Subset):
"""
Compute a local submesh on the subset, for a given topology
- :param topo: `~hysop.mpi.topology.Cartesian`
+ :param topo: :class:`~hysop.mpi.topology.Cartesian`
"""
assert isinstance(topo, Cartesian)
# Find the indices of starting/ending points in the global_end
@@ -703,9 +750,9 @@ class SubBox(Subset):
Parameters
----------
- field : `~hysop.field.continuous.Field`
+ field : :class:`~hysop.field.continuous.Field`
a field to be integrated on the box
- topo : `~hysop.mpi.topology.Cartesian`
+ topo : :class:`~hysop.mpi.topology.Cartesian`
set mesh/topology used for integration
component : int, optional
number of the field component to be integrated
@@ -732,7 +779,7 @@ class SubBox(Subset):
Parameters
----------
func: python function of space coordinates
- topo: `~hysop.mpi.topology.Cartesian`
+ topo: :class:`~hysop.mpi.topology.Cartesian`
set mesh/topology used for integration
nbc : int
number of components of the return value from func
@@ -759,7 +806,7 @@ class SubBox(Subset):
Parameters
----------
func: python function of space coordinates
- topo: `~hysop.mpi.topology.Cartesian`
+ topo: :class:`~hysop.mpi.topology.Cartesian`
set mesh/topology used for integration
Returns
@@ -783,7 +830,7 @@ class SubBox(Subset):
Parameters
----------
- field : `~hysop.field.discrete.DiscreteField`
+ field : :class:`~hysop.field.discrete.DiscreteField`
a field to be integrated on the box
component : int, optional
number of the field component to be integrated
diff --git a/HySoP/hysop/domain/tests/test_box.py b/HySoP/hysop/domain/tests/test_box.py
index 2fd63b584a4d251bf4fb25011181e08725a85c27..ca231e8ceec97c9634942eda24e8d32ad2fed1e3 100644
--- a/HySoP/hysop/domain/tests/test_box.py
+++ b/HySoP/hysop/domain/tests/test_box.py
@@ -20,11 +20,11 @@ def test_create_box1():
def test_create_box2():
- L = npw.asrealarray([1., 2.5])
+ length = npw.asrealarray([1., 2.5])
ori = npw.asrealarray([1, 2.1])
- dom = Box(length=L, origin=ori)
+ dom = Box(length=length, origin=ori)
assert dom.dimension == 2
- assert allclose(dom.length, L)
+ assert allclose(dom.length, length)
assert allclose(dom.origin, ori)
assert [b == PERIODIC for b in dom.boundaries]
cond = [dom.tasks_on_proc(i) == DEFAULT_TASK_ID for i in range(main_size)]
@@ -33,11 +33,11 @@ def test_create_box2():
def test_create_box3():
- L = [1, 2, 4.]
- dom = Box(length=L)
+ length = [1, 2, 4.]
+ dom = Box(length=length)
assert dom.dimension == 3
- assert allclose(dom.length, npw.asrealarray(L))
- assert allclose(dom.origin, zeros_like(L))
+ assert allclose(dom.length, npw.asrealarray(length))
+ assert allclose(dom.origin, zeros_like(length))
assert [b == PERIODIC for b in dom.boundaries]
cond = [dom.tasks_on_proc(i) == DEFAULT_TASK_ID for i in range(main_size)]
cond = npw.asboolarray(cond)
@@ -45,20 +45,20 @@ def test_create_box3():
def test_create_box4():
- L = [1, 2, 4.]
+ length = [1, 2, 4.]
tasks = [CPU] * main_size
if main_size > 1:
tasks[-1] = GPU
- dom = Box(length=L, proc_tasks=tasks)
+ dom = Box(length=length, proc_tasks=tasks)
last = main_size - 1
if main_size > 1:
if main_rank != last:
- assert dom.currentTask() == CPU
+ assert dom.current_task() == CPU
else:
- assert dom.currentTask() == GPU
+ assert dom.current_task() == GPU
else:
- assert dom.currentTask() == CPU
+ assert dom.current_task() == CPU
# Test topology creation ...
@@ -96,17 +96,17 @@ def test_topo_standard():
def test_topo_multi_tasks():
dom = Box(proc_tasks=proc_tasks)
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
topo = dom.create_topology(discretization=r3D)
- elif dom.isOnTask(GPU):
+ elif dom.is_on_task(GPU):
topo = dom.create_topology(discretization=r3DGh, dim=2)
assert len(dom.topologies) == 1
assert isinstance(topo, Cartesian)
assert topo is dom.topologies.values()[0]
- if dom.isOnTask(CPU):
- assert not topo.hasGhosts()
- elif dom.isOnTask(GPU):
- assert topo.hasGhosts()
+ if dom.is_on_task(CPU):
+ assert not topo.has_ghosts()
+ elif dom.is_on_task(GPU):
+ assert topo.has_ghosts()
def test_topo_plane():
@@ -125,21 +125,21 @@ def test_topo_from_mesh():
# e.g. for fftw
dom = Box(proc_tasks=proc_tasks)
from hysop.f2py import fftw2py
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
localres, global_start = fftw2py.init_fftw_solver(
r3D.resolution, dom.length, comm=comm_s.py2f())
print localres, global_start
topo = dom.create_plane_topology_from_mesh(localres=localres,
global_start=global_start,
discretization=r3D)
- elif dom.isOnTask(GPU):
+ elif dom.is_on_task(GPU):
topo = dom.create_topology(discretization=r3DGh, dim=2)
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
assert (topo.mesh.resolution == localres).all()
assert (topo.mesh.start() == global_start).all()
assert topo.dimension == 1
assert (topo.shape == [1, 1, comm_s.Get_size()]).all()
- elif dom.isOnTask(GPU):
+ elif dom.is_on_task(GPU):
assert topo.size == 1
diff --git a/HySoP/hysop/domain/tests/test_control_box.py b/HySoP/hysop/domain/tests/test_control_box.py
index c4a4ada9ee7fcf716093ad84baa1e96c8ed3de75..19eac3b1d8b2b46418eba1947c06a46bd3364b03 100644
--- a/HySoP/hysop/domain/tests/test_control_box.py
+++ b/HySoP/hysop/domain/tests/test_control_box.py
@@ -1,5 +1,5 @@
-from hysop.domain.subsets.control_box import ControlBox
-from hysop.domain.subsets.boxes import SubBox
+from hysop.domain.control_box import ControlBox
+from hysop.domain.subsets import SubBox
from hysop.tools.parameters import Discretization
from hysop import Field, Box
import numpy as np
@@ -64,7 +64,7 @@ def check_control_box(discr, xr, lr):
return topo, cb
-def test_cb_3D():
+def test_cb_3d():
topo, cb = check_control_box(discr3D, xdef, ldef)
velo = Field(domain=topo.domain, is_vector=True, formula=v3d, name='velo')
velo.discretize(topo)
@@ -84,7 +84,7 @@ def test_cb_3D():
assert np.abs(isurf - sref) < 1e-6
-def test_cb_3D_2():
+def test_cb_3d_2():
topo, cb = check_control_box(discr3D, xdef, ldef)
velo = Field(domain=topo.domain, is_vector=True, formula=v3d, name='velo')
velo.discretize(topo)
@@ -102,7 +102,7 @@ def test_cb_3D_2():
assert np.abs(isurf - sref) < 1e-6
-def test_cb_3D_3():
+def test_cb_3d_3():
topo, cb = check_control_box(discr3D, xdef, ldef)
velo = Field(domain=topo.domain, is_vector=True, formula=v3d, name='velo')
velo.discretize(topo)
@@ -120,7 +120,7 @@ def test_cb_3D_3():
assert (np.abs(isurf - sref) < 1e-6).all()
-def test_cb_2D():
+def test_cb_2d():
topo, cb = check_control_box(discr2D, xdef[:2], ldef[:2])
velo = Field(domain=topo.domain, is_vector=True, formula=v2d, name='velo')
velo.discretize(topo)
@@ -140,7 +140,7 @@ def test_cb_2D():
assert np.abs(isurf - sref) < 1e-6
-def test_cb_2D_2():
+def test_cb_2d_2():
topo, cb = check_control_box(discr2D, xdef[:2], ldef[:2])
velo = Field(domain=topo.domain, is_vector=True, formula=v2d, name='velo')
velo.discretize(topo)
@@ -158,7 +158,7 @@ def test_cb_2D_2():
assert np.abs(isurf - sref) < 1e-6
-def test_cb_2D_3():
+def test_cb_2d_3():
topo, cb = check_control_box(discr2D, xdef[:2], ldef[:2])
velo = Field(domain=topo.domain, is_vector=True, formula=v2d, name='velo')
velo.discretize(topo)
@@ -176,7 +176,7 @@ def test_cb_2D_3():
assert (np.abs(isurf - sref) < 1e-6).all()
-def test_cb_3D_fulldomain():
+def test_cb_3d_fulldomain():
"""
A cb defined on the whole domain.
Integrals on upper surfaces must fail,
@@ -206,10 +206,10 @@ def test_cb_3D_fulldomain():
# This may be useful to run mpi tests
if __name__ == "__main__":
- test_cb_3D()
- test_cb_3D_2()
- test_cb_3D_3()
- test_cb_2D()
- test_cb_2D_2()
- test_cb_2D_3()
- test_cb_3D_fulldomain()
+ test_cb_3d()
+ test_cb_3d_2()
+ test_cb_3d_3()
+ test_cb_2d()
+ test_cb_2d_2()
+ test_cb_2d_3()
+ test_cb_3d_fulldomain()
diff --git a/HySoP/hysop/domain/tests/test_mesh.py b/HySoP/hysop/domain/tests/test_mesh.py
index 993c6df42db5067019120075fc0012b11654a547..7e3e40d36b4a0dcc762bd9a5b3689dfa5b425819 100644
--- a/HySoP/hysop/domain/tests/test_mesh.py
+++ b/HySoP/hysop/domain/tests/test_mesh.py
@@ -60,22 +60,22 @@ def create_mesh(discr):
assert mesh.local_indices(point) is False
-def test_mesh3D():
+def test_mesh3d():
discr = Discretization([Nx + 1, Ny + 1, Nz + 1])
create_mesh(discr)
-def test_mesh3D_ghost():
+def test_mesh3d_ghost():
discr = Discretization([Nx + 1, Ny + 1, Nz + 1], [2, 2, 2])
create_mesh(discr)
-def test_mesh2D():
+def test_mesh2d():
discr = Discretization([Nx + 1, Nz + 1])
create_mesh(discr)
-def test_mesh2D_ghost():
+def test_mesh2d_ghost():
discr = Discretization([Nx + 1, Nz + 1], [2, 2])
create_mesh(discr)
@@ -133,9 +133,9 @@ def test_convert_global():
assert res == sl
if __name__ == '__main__':
- test_mesh3D()
- test_mesh3D_ghost()
- test_mesh2D()
- test_mesh2D_ghost()
+ test_mesh3d()
+ test_mesh3d_ghost()
+ test_mesh2d()
+ test_mesh2d_ghost()
test_convert_local()
test_convert_global()
diff --git a/HySoP/hysop/domain/tests/test_porous.py b/HySoP/hysop/domain/tests/test_porous.py
index 682c56a7ca7c896bd9f994148fc189a85f3a5eb8..b813c1d7b5a8140b50692eac9faa6e018da5bc37 100644
--- a/HySoP/hysop/domain/tests/test_porous.py
+++ b/HySoP/hysop/domain/tests/test_porous.py
@@ -2,11 +2,11 @@
from hysop import Field, Box, Discretization, IOParams
from hysop.mpi.topology import Cartesian
-from hysop.domain.subsets.sphere import HemiSphere, Sphere
-from hysop.domain.subsets.cylinder import Cylinder, HemiCylinder
+from hysop.domain.subsets import HemiSphere, Sphere
+from hysop.domain.subsets import Cylinder, HemiCylinder
import hysop.tools.numpywrappers as npw
import numpy as np
-from hysop.domain.subsets.porous import BiPole, QuadriPole,\
+from hysop.domain.porous import BiPole, QuadriPole,\
Ring, RingPole, Porous
from hysop.operator.hdf_io import HDF_Reader
diff --git a/HySoP/hysop/domain/tests/test_regular_subset.py b/HySoP/hysop/domain/tests/test_regular_subset.py
index cfb14a9485321f34d146bc0eafab2c9e3949b52b..da1b44652d5184dca3bb496bc68251a094d431b5 100644
--- a/HySoP/hysop/domain/tests/test_regular_subset.py
+++ b/HySoP/hysop/domain/tests/test_regular_subset.py
@@ -1,4 +1,4 @@
-from hysop.domain.subsets.boxes import SubBox
+from hysop.domain.subsets import SubBox
from hysop.tools.parameters import Discretization
from hysop import Field, Box
import hysop.tools.numpywrappers as npw
@@ -182,7 +182,7 @@ def test_integ_fulldom():
assert (np.abs(i1 - vref) < 1e-6).all()
-def test_integ_fulldom_2D():
+def test_integ_fulldom_2d():
topo, cub = check_subset(discr2D, xdom[:2], ldom[:2])
velo = Field(domain=topo.domain, is_vector=True, formula=v2d, name='velo')
vd = velo.discretize(topo)
@@ -228,4 +228,4 @@ if __name__ == "__main__":
test_2d_subbox()
test_2d_line()
test_integ_fulldom()
- test_integ_fulldom_2D()
+ test_integ_fulldom_2d()
diff --git a/HySoP/hysop/domain/tests/test_submesh.py b/HySoP/hysop/domain/tests/test_submesh.py
index f599b6fea722fb29d778c8f52b561e8c4d81c1e8..bec87e49d42d5fc1dd570f74a87c4abe64f326de 100644
--- a/HySoP/hysop/domain/tests/test_submesh.py
+++ b/HySoP/hysop/domain/tests/test_submesh.py
@@ -4,7 +4,7 @@ Testing regular grids.
"""
from hysop.domain.box import Box
from hysop.tools.parameters import Discretization
-from hysop.domain.subsets.submesh import SubMesh
+from hysop.domain.submesh import SubMesh
from hysop.tools.misc import utils
import numpy as np
@@ -75,7 +75,7 @@ def check_coords(m1, m2):
assert np.allclose(xend, req_end)
-def test_submesh_3D():
+def test_submesh_3d():
"""
3D Periodic mesh
"""
@@ -83,7 +83,7 @@ def test_submesh_3D():
check_submesh(discr)
-def test_submesh_2D():
+def test_submesh_2d():
"""
2D Periodic mesh
"""
@@ -91,5 +91,5 @@ def test_submesh_2D():
check_submesh(discr)
if __name__ == '__main__':
- test_submesh_3D()
- test_submesh_2D()
+ test_submesh_3d()
+ test_submesh_2d()
diff --git a/HySoP/hysop/domain/tests/test_subset.py b/HySoP/hysop/domain/tests/test_subset.py
index 73cc589ceefaa54b5ddbb5ded91c0806d277af62..0b11f5288725e28b070e7f4049f5615a4462738e 100644
--- a/HySoP/hysop/domain/tests/test_subset.py
+++ b/HySoP/hysop/domain/tests/test_subset.py
@@ -1,7 +1,6 @@
-from hysop.domain.subsets.sphere import Sphere, HemiSphere
-from hysop.domain.subsets.subset import Subset
-from hysop.domain.subsets.boxes import SubBox
-from hysop.domain.subsets.cylinder import Cylinder, HemiCylinder
+from hysop.domain.subsets import Sphere, HemiSphere
+from hysop.domain.subsets import Subset, SubBox
+from hysop.domain.subsets import Cylinder, HemiCylinder
from hysop.operator.hdf_io import HDF_Reader
from hysop.tools.parameters import Discretization, IOParams
from hysop import Field, Box
@@ -55,8 +54,8 @@ def check_subset(discr, fileref, class_name):
scal = Field(domain=topo.domain, is_vector=False)
sd = scal.discretize(topo)
sd[0][subs.ind[topo][0]] = 2.
- iCompute = topo.mesh.iCompute
- return np.allclose(sd[0][iCompute], vd[0][iCompute])
+ icompute = topo.mesh.iCompute
+ return np.allclose(sd[0][icompute], vd[0][icompute])
def test_sphere_3d():
@@ -117,8 +116,8 @@ def init_obst_list(discr, fileref):
sd = Field(domain=topo.domain, is_vector=False).discretize(topo)
indices = Subset.union(obs_list, topo)
sd[0][indices] = 2.
- iCompute = topo.mesh.iCompute
- return np.allclose(sd[0][iCompute], vd[0][iCompute])
+ icompute = topo.mesh.iCompute
+ return np.allclose(sd[0][icompute], vd[0][icompute])
def test_union_3d():
@@ -142,8 +141,8 @@ def init_subtract(discr, fileref):
sd = Field(domain=topo.domain, is_vector=False).discretize(topo)
indices = Subset.subtract(box, s1, topo)
sd[0][indices] = 2.
- iCompute = topo.mesh.iCompute
- return np.allclose(sd[0][iCompute], vd[0][iCompute])
+ icompute = topo.mesh.iCompute
+ return np.allclose(sd[0][icompute], vd[0][icompute])
def test_subtract_3d():
diff --git a/HySoP/hysop/mpi/__init__.py b/HySoP/hysop/mpi/__init__.py
index f990cd8c790a9effd03543088d92b99edf09e8a5..94aa37bd9ca80eb741e837d584331b2efd9a8a29 100644
--- a/HySoP/hysop/mpi/__init__.py
+++ b/HySoP/hysop/mpi/__init__.py
@@ -1,33 +1,36 @@
-"""
-@package hysop.mpi
-hysop interface to the mpi implementation.
+"""Hysop interface to the mpi implementation.
It contains :
-- mpi basic variables (main communicator, rank, size ...)
-- topology.Cartesian : mpi process distribution + local mesh
-- mesh.SubMesh : local mesh
+
+* mpi basic variables (main communicator, rank, size ...)
+* :class:`hysop.mpi.topology.Cartesian` : mpi process distribution + local mesh
This package is used to hide the underlying mpi interface
in order to make any change of this interface, if required, easiest.
+
At the time we use mpi4py : http://mpi4py.scipy.org
"""
-## Everything concerning the chosen mpi implementation is hidden in main_var
+# Everything concerning the chosen mpi implementation is hidden in main_var
# Why? --> to avoid that things like mpi4py. ... spread everywhere in the
# soft so to ease a change of this implementation (if needed).
from hysop.mpi import main_var
-## A list of mpi variables that can be "seen" by user
-## MPI underlying implementation
+
MPI = main_var.MPI
+"""MPI underlying implementation
+"""
-## Main communicator (copy of comm_world)
main_comm = main_var.main_comm
+"""Main communicator (copy of comm_world)
+"""
-## Rank of the current process in the main communicator
main_rank = main_var.main_rank
+"""Rank of the current process in the main communicator
+"""
-## Size of the main communicator
main_size = main_var.main_size
+"""Size of the main communicator
+"""
diff --git a/HySoP/hysop/mpi/bridge.py b/HySoP/hysop/mpi/bridge.py
index 33957e1fe76f71c64570498e695ef5ae4e39bba8..b9765d94dbaa751aac8fd52a05b51889bdaa85b8 100644
--- a/HySoP/hysop/mpi/bridge.py
+++ b/HySoP/hysop/mpi/bridge.py
@@ -3,7 +3,7 @@
Tools to compute the intersection between
two HySoP topologies.
"""
-from hysop.mpi.topology import Cartesian, topotools
+from hysop.mpi.topology import Cartesian, TopoTools
from hysop.tools.misc import utils
@@ -51,7 +51,7 @@ class Bridge(object):
msg = 'Bridge error, both source/target topologies'
msg += ' must have the same parent communicator.'
- assert topotools.compare_comm(self._source.parent(),
+ assert TopoTools.compare_comm(self._source.parent(),
self._target.parent()), msg
# The assert above ensure that source and target hold the same
# group of process in the same communication context.
@@ -63,10 +63,10 @@ class Bridge(object):
# are on both source and target mesh.
# Get global indices of the mesh on source for all mpi processes.
- indices_source = topotools.gatherGlobalIndices(self._source)
+ indices_source = TopoTools.gather_global_indices(self._source)
# Get global indices of the mesh on target for all mpi processes.
- indices_target = topotools.gatherGlobalIndices(self._target)
+ indices_target = TopoTools.gather_global_indices(self._target)
# From now on, we have indices_source[rk] = global indices (slice)
# of grid points of the source on process number rk in parent.
# And the same thing for indices_target.
@@ -128,8 +128,8 @@ class Bridge(object):
"""
if self._recv_types is None:
data_shape = self._target.mesh.resolution
- self._recv_types = topotools.createSubArray(self._recv_indices,
- data_shape)
+ self._recv_types = TopoTools.create_subarray(self._recv_indices,
+ data_shape)
return self._recv_types
def sendTypes(self):
@@ -141,6 +141,6 @@ class Bridge(object):
"""
if self._send_types is None:
data_shape = self._source.mesh.resolution
- self._send_types = topotools.createSubArray(self._send_indices,
- data_shape)
+ self._send_types = TopoTools.create_subarray(self._send_indices,
+ data_shape)
return self._send_types
diff --git a/HySoP/hysop/mpi/bridge_inter.py b/HySoP/hysop/mpi/bridge_inter.py
index 6516bd7abcf6071264448b544295b13007b1cc40..b800d6988b0d79f5ee0c89692d020146df3d2fb0 100644
--- a/HySoP/hysop/mpi/bridge_inter.py
+++ b/HySoP/hysop/mpi/bridge_inter.py
@@ -3,7 +3,7 @@
Tools to compute the intersection between
two HySoP topologies.
"""
-from hysop.mpi.topology import Cartesian, topotools
+from hysop.mpi.topology import Cartesian, TopoTools
from hysop.tools.misc import utils
from hysop.mpi import MPI
import hysop.tools.numpywrappers as npw
@@ -40,7 +40,7 @@ class BridgeInter(object):
assert isinstance(parent, MPI.Intracomm)
self._topology = current
# current task id
- current_task = self._topology.domain.currentTask()
+ current_task = self._topology.domain.current_task()
# True if current process is in the 'from' group'
task_is_source = current_task == self.source_id
@@ -87,7 +87,7 @@ class BridgeInter(object):
def _swap_indices(self):
# First, we need to collect the global indices, as arrays
# since we need to broadcast them later.
- current_indices = topotools.gatherGlobalIndices(self._topology,
+ current_indices = TopoTools.gather_global_indices(self._topology,
toslice=False)
# To allocate remote_indices array, we need the size of
# the remote communicator.
@@ -96,7 +96,7 @@ class BridgeInter(object):
remote_indices = npw.dim_zeros((dimension * 2, remote_size))
# Then they are broadcasted to the remote communicator
rank = self._topology.rank
- current_task = self._topology.domain.currentTask()
+ current_task = self._topology.domain.current_task()
if current_task is self.source_id:
# Local 0 broadcast current_indices to remote comm
if rank == 0:
@@ -130,6 +130,6 @@ class BridgeInter(object):
"""
if self._transfer_types is None:
data_shape = self._topology.mesh.resolution
- self._transfer_types = topotools.createSubArray(
+ self._transfer_types = TopoTools.create_subarray(
self._tranfer_indices, data_shape)
return self._transfer_types
diff --git a/HySoP/hysop/mpi/bridge_overlap.py b/HySoP/hysop/mpi/bridge_overlap.py
index 64585e68b31fe4e87a31e96952fe6cbe6c9772b6..4dfa0b904ed6c83bdd59d0f250ed00c4b41c22dc 100644
--- a/HySoP/hysop/mpi/bridge_overlap.py
+++ b/HySoP/hysop/mpi/bridge_overlap.py
@@ -4,7 +4,7 @@ Tools to compute the intersection between
two HySoP topologies defined on the same comm but for a
different number of processes.
"""
-from hysop.mpi.topology import Cartesian, topotools
+from hysop.mpi.topology import Cartesian, TopoTools
from hysop.tools.misc import utils
from hysop.mpi.bridge import Bridge
@@ -48,21 +48,21 @@ class BridgeOverlap(Bridge):
msg = 'BridgeOverlap error: mpi group from '
msg += 'source and topo must overlap. If not '
msg += 'BridgeInter will probably suits better.'
- assert topotools.intersection_size(self._source.comm,
+ assert TopoTools.intersection_size(self._source.comm,
self._target.comm) > 0, msg
assert self._source.domain == self._target.domain
if self._source is not None:
assert isinstance(self._source, Cartesian)
s_size = self._source.size
- assert topotools.intersection_size(self._source.comm,
+ assert TopoTools.intersection_size(self._source.comm,
self.comm) == s_size
self.domain = self._source.domain
if self._target is not None:
assert isinstance(self._target, Cartesian)
t_size = self._target.size
- assert topotools.intersection_size(self._target.comm,
+ assert TopoTools.intersection_size(self._target.comm,
self.comm) == t_size
self.domain = self._target.domain
@@ -71,12 +71,12 @@ class BridgeOverlap(Bridge):
def _build_send_recv_dict(self):
# Compute local intersections : i.e. find which grid points
# are on both source and target mesh.
- indices_source = topotools.gatherGlobalIndicesOverlap(self._source,
- self.comm,
- self.domain)
- indices_target = topotools.gatherGlobalIndicesOverlap(self._target,
- self.comm,
- self.domain)
+ indices_source = TopoTools.gather_global_indices_overlap(self._source,
+ self.comm,
+ self.domain)
+ indices_target = TopoTools.gather_global_indices_overlap(self._target,
+ self.comm,
+ self.domain)
# From now on, we have indices_source[rk] = global indices (slice)
# of grid points of the source on process number rk in parent.
diff --git a/HySoP/hysop/mpi/tests/test_bridge.py b/HySoP/hysop/mpi/tests/test_bridge.py
index fa02f3b69ca84d28f565c14b1abbb79821a49d23..0621466c23680c0181c82c2c6a445d8aebbb20b0 100644
--- a/HySoP/hysop/mpi/tests/test_bridge.py
+++ b/HySoP/hysop/mpi/tests/test_bridge.py
@@ -88,9 +88,9 @@ def test_bridgeInter2D():
# The real tests are done in test_redistribute.py
# Bridge from topo2 on GPU to topo1 on CPU:
- if dom.isOnTask(GPU):
+ if dom.is_on_task(GPU):
bridge2 = BridgeInter(topo2, main_comm, source_id=GPU, target_id=CPU)
- elif dom.isOnTask(CPU):
+ elif dom.is_on_task(CPU):
bridge2 = BridgeInter(topo1, main_comm, source_id=GPU, target_id=CPU)
assert bridge2 is not None
@@ -114,9 +114,9 @@ def test_bridgeInter3D():
# so we just create a bridge.
# The real tests are done in test_redistribute.py
# Bridge from topo2 on GPU to topo1 on CPU:
- if dom.isOnTask(GPU):
+ if dom.is_on_task(GPU):
bridge2 = BridgeInter(topo2, main_comm, source_id=GPU, target_id=CPU)
- elif dom.isOnTask(CPU):
+ elif dom.is_on_task(CPU):
bridge2 = BridgeInter(topo1, main_comm, source_id=GPU, target_id=CPU)
assert bridge2 is not None
diff --git a/HySoP/hysop/mpi/tests/test_topology.py b/HySoP/hysop/mpi/tests/test_topology.py
index 926eb3af659c3fe581a6827be56b1d303674660b..23e6fcb0095075e99cbee9ae29a5123900e21afb 100644
--- a/HySoP/hysop/mpi/tests/test_topology.py
+++ b/HySoP/hysop/mpi/tests/test_topology.py
@@ -62,9 +62,9 @@ def test_create_default_topology2_2d():
topo = dom.create_topology(r2D)
assert topo.domain == dom
assert topo.size == dom2D.comm_task.Get_size()
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
assert topo.task_id() == CPU
- if dom.isOnTask(GPU):
+ if dom.is_on_task(GPU):
assert topo.task_id() == GPU
@@ -151,9 +151,9 @@ def test_create_default_topology2():
topo = dom.create_topology(r3D)
assert topo.domain == dom
assert topo.size == dom3D.comm_task.Get_size()
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
assert topo.task_id() == CPU
- if dom.isOnTask(GPU):
+ if dom.is_on_task(GPU):
assert topo.task_id() == GPU
diff --git a/HySoP/hysop/mpi/topology.py b/HySoP/hysop/mpi/topology.py
index b11214c134f513157ab7a504bf77840df3dbaed3..41dd5d78ac41e906a9c94b3cae8cbb968f339d8a 100644
--- a/HySoP/hysop/mpi/topology.py
+++ b/HySoP/hysop/mpi/topology.py
@@ -1,8 +1,9 @@
-"""
-@file topology.py
-Tools and definitions for HySoP topologies
+"""Tools and definitions for HySoP topologies
(MPI Processes layout + space discretization)
+* :class:`~hysop.mpi.topology.Cartesian`
+* :class:`~hysop.mpi.topology.TopoTools`
+
"""
from hysop.constants import debug, ORDER, PERIODIC
@@ -51,34 +52,38 @@ class Cartesian(object):
def __init__(self, domain, discretization, dim=None, mpi_params=None,
isperiodic=None, cutdir=None, shape=None, mesh=None):
"""
- Required parameters : domain, discretization
- Others are optional. You must choose one and only one param
- among dim, cutdir and shape.
+
+ Parameters
+ -----------
+ domain : :class:`~hysop.domain.domain.Domain`
+ the geometry on which the topology is defined.
+ discretization : :class:`~hysop.tools.parameters.Discretization`
+ description of the global space discretization
+ (resolution and ghosts).
+ dim : int, optional
+ mpi topology dimension.
+ mpi_params : :class:`~hysop.tools.parameters.MPIParams`, optional
+ mpi setup (comm, task ...).
+ If None, comm = main_comm, task = DEFAULT_TASK_ID.
+ isperiodic : list or array of bool, optional
+ mpi grid periodicity
+ cutdir : list or array of bool
+ set which directions must (may) be distributed,
+ cutdir[dir] = True to distribute data along dir.
+ shape : list or array of int
+ mpi grid layout
+ mesh : :class:`~hysop.domain.mesh.Mesh`, optional
+ a predifined mesh (it includes local and global grids.)
+
+ Notes:
+ ------
+ Almost all parameters above are optional.
+ Only one must be chosen among dim, cutdir and shape.
+ See :ref:`topologies` in the Hysop User Guide for details.
See hysop.mpi.topology.Cartesian.plane_precomputed
details to build a plane topology from a given local discretization
(e.g. from fftw or scales precomputation).
- @param domain : the geometry; it must be a box.
- @param discretization : a hysop.tools.parameters.Discretization
- with:
- - resolution = Number of points in the domain
- in each direction. We assume that first point corresponds
- to origin, and last point to boundary point,
- whatever the boundary type is.
- That is x[0] = domain.origin and
- x[Discretization.resolution-1] = domain.Lengths_x.
- - ghosts = number of points in the ghost layer
- @param dim : dimension of the topology
- @param mpi_params : a hysop.tools.parameters.MPIParams, with:
- - comm : MPI communicator used to create this topology
- (default = main_comm)
- - task_id : id of the task that owns this topology.
- @param isperiodic : periodicity of the topology (in each direction)
- @param cutdir : array of bool, set cutdir[dir] = True if you want
- to distribute data through direction dir.
- @param shape : topology resolution
- (i.e process layout in each direction).
- @param mesh : a predefined hysop.mpi.mesh.SubMesh
"""
# ===== 1 - Required parameters : domain and mpi (comm, task) ====
# An id for the topology
@@ -89,10 +94,10 @@ class Cartesian(object):
# - MPI topo params :
if mpi_params is None:
mpi_params = MPIParams(comm=domain.comm_task,
- task_id=domain.currentTask())
+ task_id=domain.current_task())
- ## mpi parameters : Communicator used to build the topology
- ## and task_id
+ # mpi parameters : Communicator used to build the topology
+ # and task_id
self._mpis = mpi_params
# Each topo must be associated to one and only one topology
assert self._mpis.task_id is not None
@@ -105,50 +110,50 @@ class Cartesian(object):
# (c) - from dimension of the topology ==> let MPI
# choose the 'best' layout.
- ## Layout of the grid of mpi processes.
+ # Layout of the grid of mpi processes.
self.shape = None
- ## Dim of the cartesian grid for MPI processes.
+ # Dim of the cartesian grid for MPI processes.
self.dimension = None
# MPI grid periodicity
self._isperiodic = None
# directions where data are distributed
self.cutdir = None
- ## mpi communicator (cartesian topo)
+ # mpi communicator (cartesian topo)
self.comm = None
self._build_mpi_topo(shape, cutdir, dim, isperiodic)
# ===== 3 - Get features of the mpi processes grid =====
- ## Size of the topology (i.e. total number of mpi processes)
+ # Size of the topology (i.e. total number of mpi processes)
self.size = self.comm.Get_size()
- ## Rank of the current process in the topology
+ # Rank of the current process in the topology
self.rank = self.comm.Get_rank()
- ## Coordinates of the current process
+ # Coordinates of the current process
reduced_coords = npw.asdimarray(self.comm.Get_coords(self.rank))
- ## Coordinates of the current process
- ## What is different between proc_coords and reduced_coords?
- ## --> proc_coords has values even for directions that
- ## are not distributed. If cutdir = [False, False, True]
- ## then reduced_coords = [ nx ] and proc_coords = [0, 0, nx]
+ # Coordinates of the current process
+ # What is different between proc_coords and reduced_coords?
+ # --> proc_coords has values even for directions that
+ # are not distributed. If cutdir = [False, False, True]
+ # then reduced_coords = [ nx ] and proc_coords = [0, 0, nx]
self.proc_coords = npw.dim_zeros(self.domain.dimension)
self.proc_coords[self.cutdir] = reduced_coords
- ## Neighbours : self.neighbours[0,i] (resp. [1,i])
- ## previous (resp. next) neighbour in direction i
- ## (warning : direction in the grid of process).
+ # Neighbours : self.neighbours[0,i] (resp. [1,i])
+ # previous (resp. next) neighbour in direction i
+ # (warning : direction in the grid of process).
self.neighbours = npw.dim_zeros((2, self.dimension))
for direction in range(self.dimension):
self.neighbours[:, direction] = self.comm.Shift(direction, 1)
# ===== 4 - Computation of the local meshes =====
- ## Local mesh on the current mpi process.
+ # Local mesh on the current mpi process.
# mesh from external function (e.g. fftw, scales ...)
self.mesh = mesh
# If mesh is None, we must compute local resolution and other
# parameters, using discretization.
if mesh is None:
- self._computeMesh(discretization)
+ self._compute_mesh(discretization)
# ===== 5 - Final setup ====
@@ -168,11 +173,11 @@ class Cartesian(object):
self.mesh = topo.mesh
self.comm = topo.comm
-# self.is_registered = False
-
- def _build_mpi_topo(self, shape, cutdir, dim, isperiodic):
- """
- Build mpi topology
+ def _build_mpi_topo(self, shape=None, cutdir=None,
+ dim=None, isperiodic=None):
+ """ default builder
+ See :ref:`topologies` for details.
+
"""
# number of process in parent comm
origin_size = self._mpis.comm.Get_size()
@@ -197,7 +202,7 @@ class Cartesian(object):
self.dimension = self.cutdir[self.cutdir].size
shape = npw.asdimarray(MPI.Compute_dims(origin_size,
self.dimension))
- self.optimizeshape(shape)
+ self._optimizeshape(shape)
self.shape = npw.dim_ones(self.domain.dimension)
self.shape[self.cutdir] = shape
@@ -227,10 +232,10 @@ class Cartesian(object):
# For C-like (row major) arrays, first dir is the
# first to be distributed, and so on.
self.shape[:self.dimension] = shape
- self.optimizeshape(self.shape)
+ self._optimizeshape(self.shape)
self.cutdir = self.shape != 1
- ## MPI processes grid periodicity. Default is true.
+ # MPI processes grid periodicity. Default is true.
if isperiodic is None:
self._isperiodic = npw.ones((self.domain.dimension),
dtype=npw.bool)
@@ -257,9 +262,8 @@ class Cartesian(object):
reorder=True)
@staticmethod
- def optimizeshape(shape):
- """
- Reorder 'shape' according to the chosen
+ def _optimizeshape(shape):
+ """Reorder 'shape' according to the chosen
data layout to optimize data distribution.
"""
shape.sort()
@@ -267,34 +271,36 @@ class Cartesian(object):
shape[:] = shape[::-1]
def parent(self):
- """
- returns the communicator used to build this topology
+ """returns the communicator used to build this topology
"""
return self._mpis.comm
def ghosts(self):
- """
- Get ghost layer width.
+ """returns ghost layer width.
"""
return self.mesh.discretization.ghosts
def task_id(self):
- """
- @return id of the task that owns this topology
+ """returns id of the task that owns this topology
"""
return self._mpis.task_id
@classmethod
def plane_precomputed(cls, localres, global_start, cdir=None, **kwds):
- """
- Define a 'plane' (1D) topology for a given mesh resolution.
+ """Defines a 'plane' (1D) topology for a given mesh resolution.
+
This function is to be used when topo/discretization features
come from an external routine (e.g. scales or fftw)
- @param localres : local mesh resolution
- @param global_start : global indices of the lowest point
- of the local mesh
- @param cdir : direction of cutting (i.e. normal to mpi planes)
- default = last if fortran order, first if C order.
+
+ Parameters
+ ----------
+ localres : list or array of int
+ local mesh resolution
+ global_start : list or array of int
+ indices in the global mesh of the lowest point of the local mesh.
+ cdir : int, optional
+ direction where data must be distributed in priority.
+ Default = last if fortran order, first if C order.
"""
msg = 'parameter is not required for plane_precomputed'
msg += ' topology construction.'
@@ -320,7 +326,7 @@ class Cartesian(object):
return cls(mesh=mesh, cutdir=cutdir, **kwds)
- def _computeMesh(self, discretization):
+ def _compute_mesh(self, discretization):
assert isinstance(discretization, Discretization)
assert discretization.resolution.size == self.domain.dimension
assert self.domain.dimension == discretization.resolution.size, \
@@ -334,28 +340,28 @@ class Cartesian(object):
is_periodic = len(np.where(self.domain.boundaries == PERIODIC)[0])\
== self.domain.dimension
if is_periodic:
- resolCalc = discretization.resolution - 1
+ computational_grid_resolution = discretization.resolution - 1
else:
raise AttributeError('Unknwon boundary conditions.')
- pts_noghost[:] = resolCalc // self.shape
+ pts_noghost[:] = computational_grid_resolution // self.shape
# If any, remaining points are
# added on the mesh of the last process.
- remainingPts = npw.dim_zeros(self.domain.dimension)
- remainingPts[:] = resolCalc % self.shape
+ remaining_points = npw.dim_zeros(self.domain.dimension)
+ remaining_points[:] = computational_grid_resolution % self.shape
# Total number of points (size of arrays to be allocated)
nbpoints = pts_noghost.copy()
for i in range(self.domain.dimension):
if self.proc_coords[i] == self.shape[i] - 1:
- nbpoints[i] += remainingPts[i]
+ nbpoints[i] += remaining_points[i]
- local_resolution = resolCalc.copy()
+ local_resolution = computational_grid_resolution.copy()
local_resolution[:] = nbpoints[:]
local_resolution[:] += 2 * discretization.ghosts[:]
- ## Global indices for the local mesh points
+ # Global indices for the local mesh points
global_start = npw.dim_zeros((self.domain.dimension))
global_start[:] = self.proc_coords[:] * pts_noghost[:]
@@ -363,12 +369,9 @@ class Cartesian(object):
local_resolution, global_start)
def __eq__(self, other):
- """
- Operator for comparison. Based on:
- @li mesh
- @li shape
- @li domain
- @param other : Topology to compare with.
+ """Comparison of two topologies.
+
+ Two topos are equal if they have the same mesh, shape and domain.
"""
if self.__class__ != other.__class__:
return False
@@ -377,11 +380,10 @@ class Cartesian(object):
self.domain == other.domain
def __ne__(self, other):
- """
- Not equal operator.
+ """Not equal operator.
+
Seems to be required in addition to __eq__ to
avoid 'corner-case' behaviors.
- @param other : Topology to compare with.
"""
result = self.__eq__(other)
if result is NotImplemented:
@@ -389,7 +391,6 @@ class Cartesian(object):
return not result
def __str__(self):
- """ Topology info display """
s = '======== Topology id : ' + str(self.__id) + ' ========\n'
s += ' - on task : ' + str(self.task_id()) + '\n'
s += ' - shape :' + str(self.shape) + '\n'
@@ -402,7 +403,7 @@ class Cartesian(object):
s += '\n=================================\n'
return s
- def hasGhosts(self):
+ def has_ghosts(self):
"""
True if ghost layer length is not zero.
"""
@@ -415,20 +416,21 @@ class Cartesian(object):
"""
return self.__id
- def isConsistentWith(self, topo):
- same_parent = self.parent() == topo.parent()
+ def is_consistent_with(self, target):
+ """True if target and current object are equal and
+ have the same parent. Equal means same mesh, same shape and
+ same domain.
+ """
+ same_parent = self.parent() == target.parent()
# Note FP. Is it really required to have the
# same parent? Inclusion of all proc may be enough?
- return npw.equal(self.shape, topo.shape).all() and same_parent
+ return npw.equal(self.shape, target.shape).all() and same_parent
- def canCommunicateWith(self, target):
- """
- return True if current topo is complient with topo.
- (See below for what 'complient' implies)
- @param target : the targeted topology
+ def can_communicate_with(self, target):
+ """True if current topo is complient with target.
- Complient if :
- - all processes in current are in topo
+ Complient means :
+ - all processes in current are in target
- both topologies belong to the same mpi task
"""
if self == target:
@@ -438,44 +440,51 @@ class Cartesian(object):
msg += ' InterBridge.'
assert self.task_id() == target.task_id(), msg
- ## Parent communicator
- ## Todo : define some proper conditions for compatibility
- ## between topo_from, topo_to and parent:
- ## - same size
- ## - same domain
- ## - common processus ...
- ## At the time we check that both topo have
- ## the same comm_origin.
- return self.isConsistentWith(target)
+ # Parent communicator
+ # Todo : define some proper conditions for compatibility
+ # between topo_from, topo_to and parent:
+ # - same size
+ # - same domain
+ # - common processus ...
+ # At the time we check that both topo have
+ # the same comm_origin.
+ return self.is_consistent_with(target)
@staticmethod
def reset_counter():
Cartesian.__ids = count(0)
-class topotools(object):
+class TopoTools(object):
@staticmethod
- def gatherGlobalIndices(topo, toslice=True, root=None, comm=None):
- """
- Collect global indices of local meshes on each process of topo
- @param topo : a hysop.mpi.topology.Cartesian
- @param toslice : true (default) if you want a dict of slice
- as return value, false if you need a numpy array.
- @param root : rank (in topo.parent()) of the gathering process.
- If root is None (default) indices are gathered
- on all processes of topo.
- @param comm : communicator used for global communications.
- Default = topo.parent().
- @return : either :
- - a dictionnary which maps rank number with
+ def gather_global_indices(topo, toslice=True, root=None, comm=None):
+ """Collect global indices of local meshes on each process of topo
+
+ Parameters
+ ----------
+ topo : :class:`~ hysop.mpi.topology.Cartesian`
+ topology on which indices are collected.
+ toslice : boolean, optional
+ true to return the result a dict of slices, else
+ return a numpy array. See notes below.
+ root : int, optional
+ rank of the root mpi process. If None, reduce operation
+ is done on all processes.
+ comm : mpi communicator, optional
+ Communicator used to reduce indices. Default = topo.parent
+
+ Returns
+ -------
+ either :
+ * a dictionnary which maps rank number with
a list of slices such that res[rank][i] = a slice
defining the indices of the points of the local mesh,
in direction i, in global notation.
- - a numpy array, each column corresponds to a rank number,
+ * or a numpy array where each column corresponds to a rank number,
with column = [start_x, end_x, start_y, end_y ...]
+ Ranks number are the processes numbers in comm.
- Ranks number are the processes numbers in topo.parent().
"""
if comm is None:
comm = topo.parent()
@@ -502,11 +511,11 @@ class topotools(object):
return iglob
@staticmethod
- def gatherGlobalIndicesOverlap(topo=None, comm=None, dom=None,
- toslice=True, root=None):
- """
- This functions does the same thing as gatherGlobalIndices but
+ def gather_global_indices_overlap(topo=None, comm=None, dom=None,
+ toslice=True, root=None):
+ """This functions does the same thing as gather_global_indices but
may also work when topo is None.
+
The function is usefull if you need to collect global indices
on a topo define only on a subset of comm,
when for the procs not in this subset, topo will be
@@ -516,9 +525,36 @@ class topotools(object):
between the two groups of processes of the topologies.
In that case, a call to
- gatherGlobalIndices(topo, comm, dom)
+ gather_global_indices(topo, comm, dom)
will work on all processes belonging to comm, topo being None or not.
The values corresponding to ranks not in topo will be empty slices.
+
+ Parameters
+ ----------
+ topo : :class:`~ hysop.mpi.topology.Cartesian`, optional
+ topology on which indices are collected.
+ toslice : boolean, optional
+ true to return the result a dict of slices, else
+ return a numpy array. See notes below.
+ root : int, optional
+ rank of the root mpi process. If None, reduce operation
+ is done on all processes.
+ comm : mpi communicator, optional
+ Communicator used to reduce indices. Default = topo.parent
+ dom : :class:`~hysop.domain.domain.Domain`
+ current domain.
+
+ Returns
+ -------
+ either :
+ * a dictionnary which maps rank number with
+ a list of slices such that res[rank][i] = a slice
+ defining the indices of the points of the local mesh,
+ in direction i, in global notation.
+ * or a numpy array where each column corresponds to a rank number,
+ with column = [start_x, end_x, start_y, end_y ...]
+ Ranks number are the processes numbers in comm.
+
"""
if topo is None:
assert comm is not None and dom is not None
@@ -538,12 +574,12 @@ class topotools(object):
return iglob
else:
- return topotools.gatherGlobalIndices(topo, toslice, root, comm)
+ return TopoTools.gather_global_indices(topo, toslice, root, comm)
@staticmethod
def is_parent(child, parent):
"""
- Return true if all mpi processes of child belongs to parent
+ Return true if all mpi processes of child belong to parent
"""
# Get the list of processes
assert child is not None
@@ -556,6 +592,8 @@ class topotools(object):
@staticmethod
def intersection_size(comm_1, comm_2):
+ """Number of processess common to comm_1 and comm_2
+ """
if comm_1 == MPI.COMM_NULL or comm_2 == MPI.COMM_NULL:
return None
group_1 = comm_1.Get_group()
@@ -565,12 +603,11 @@ class topotools(object):
@staticmethod
def compare_comm(comm_1, comm_2):
- """
- Compare two mpi communicators.
- Return true if the two communicators are handles for the same
+ """Compare two mpi communicators.
+
+ Returns true if the two communicators are handles for the same
group of proc and for the same communication context.
- @param comm_1 : a mpi communicator
- @param comm_2 : a mpi communicator
+
Warning : if comm_1 or comm_2 is invalid, the
function will fail.
"""
@@ -582,12 +619,11 @@ class topotools(object):
@staticmethod
def compare_groups(comm_1, comm_2):
- """
- Compare two mpi communicators.
- Return true if each comm handles the
+ """Compare the groups of two mpi communicators.
+
+ Returns true if each comm handles the
same group of mpi processes.
- @param comm_1 : a mpi communicator
- @param comm_2 : a mpi communicator
+
Warning : if comm_1 or comm_2 is invalid, the
function will fail.
"""
@@ -599,11 +635,13 @@ class topotools(object):
@staticmethod
def convert_ranks(source, target):
- """
- find the values of ranks in target from ranks in source.
- @param source : a mpi communicator
- @param target : a mpi communicator
- @return a list 'ranks' such that ranks[i] = rank in target
+ """Find the values of ranks in target from ranks in source.
+
+ Parameters
+ ----------
+ source, target : mpi communicators
+
+ Returns a list 'ranks' such that ranks[i] = rank in target
of process of rank i in source.
"""
assert source != MPI.COMM_NULL and target != MPI.COMM_NULL
@@ -615,15 +653,18 @@ class topotools(object):
return {r_source[i]: res[i] for i in xrange(size_source)}
@staticmethod
- def createSubArray(sl_dict, data_shape):
- """
- Create a MPI subarray mask to be used in send/recv operations
+ def create_subarray(sl_dict, data_shape):
+ """Create a MPI subarray mask to be used in send/recv operations
between some topologies.
- @param[in] sl_dict : dictionnary which contains mesh indices
- of the subarray for each rank,
- such that sl_dict[rk] = (slice(...), slice(...), ...)
- @param[in] data_shape : shape (numpy-like) of the original array
- @return : dictionnary of MPI derived types.
+
+ Parameters
+ ----------
+ sl_dict : dictionnary
+ indices of the subarray for each rank,
+ such that sl_dict[rk] = (slice(...), slice(...), ...)
+ data_shape : shape (numpy-like) of the original array
+
+ :Returns : dictionnary of MPI derived types.
Keys = ranks in parent communicator.
"""
from hysop.constants import HYSOP_MPI_REAL, ORDERMPI
diff --git a/HySoP/hysop/operator/absorption_BC.py b/HySoP/hysop/operator/absorption_BC.py
index e6ed682ac9c095ae2516b67329c291279675dd9e..303ac78ae6344c1f3b70a3a05d62c58e525f5c59 100755
--- a/HySoP/hysop/operator/absorption_BC.py
+++ b/HySoP/hysop/operator/absorption_BC.py
@@ -1,17 +1,14 @@
# -*- coding: utf-8 -*-
-"""
-@file operator/absorption_BC.py
-
-Operator to kill the vorticity at the outlet boundary
-(i.e. removal of the periodic BC in the flow direction
-by vorticity absorption in order to set the far field
+"""Operator to kill the vorticity at the outlet boundary
+(i.e. removal of the periodic BC in the flow direction
+by vorticity absorption in order to set the far field
velocity to u_inf at the inlet)
"""
from hysop.constants import debug
from hysop.operator.discrete.absorption_BC import AbsorptionBC_D
from hysop.operator.computational import Computational
-from hysop.domain.subsets.control_box import ControlBox
+from hysop.domain.control_box import ControlBox
from hysop.operator.continuous import opsetup
@@ -25,7 +22,7 @@ class AbsorptionBC(Computational):
"""
@debug
- def __init__(self, velocity, vorticity, req_flowrate,
+ def __init__(self, velocity, vorticity, req_flowrate,
x_coords_absorp, **kwds):
"""
@param[in] velocity field
diff --git a/HySoP/hysop/operator/adapt_timestep.py b/HySoP/hysop/operator/adapt_timestep.py
index 9f21d3ff7a005d52262026d676d53b3f3417d1cc..0079aa1eb6c29f4f82846e05cfcef484cf790c23 100755
--- a/HySoP/hysop/operator/adapt_timestep.py
+++ b/HySoP/hysop/operator/adapt_timestep.py
@@ -75,7 +75,7 @@ class AdaptTimeStep(Computational):
Create intercommunicators, if required (i.e. if there are several
tasks defined in the domain).
"""
- task_is_source = self._mpis.task_id == self.domain.currentTask()
+ task_is_source = self._mpis.task_id == self.domain.current_task()
tasks_list = self.domain.tasks_list()
others = (v for v in tasks_list if v != self._mpis.task_id)
if task_is_source:
@@ -127,7 +127,7 @@ class AdaptTimeStep(Computational):
self._is_uptodate = True
def wait(self):
- task_is_source = self._mpis.task_id == self.domain.currentTask()
+ task_is_source = self._mpis.task_id == self.domain.current_task()
rank = self._mpis.rank
dt = self.simulation.timeStep
for rk in self._intercomms:
diff --git a/HySoP/hysop/operator/continuous.py b/HySoP/hysop/operator/continuous.py
index a45d6cf377faa3f6fbe90139a213e9a98eeec8f3..a470261f34cb9ec53dd60b72ae4e45c9f985135c 100644
--- a/HySoP/hysop/operator/continuous.py
+++ b/HySoP/hysop/operator/continuous.py
@@ -9,8 +9,7 @@ import hysop.tools.io_utils as io
class Operator(object):
- """
- Abstract interface to continuous operators.
+ """Abstract interface to continuous operators.
"""
__metaclass__ = ABCMeta
@@ -38,21 +37,23 @@ class Operator(object):
The elements of the list or the keys of the dict
are :class:`hysop.fields.continuous.Fields`.
- The values of the dict can be either :class:`hysop.mpi.topology.Cartesian`
- or :class:`hysop.tools.parameters.Discretization`.
+ The values of the dict can be either
+ :class:`hysop.mpi.topology.Cartesian`
+ or :class:`hysop.tools.parameters.Discretization`::
- op = Analytic(variables = [velo, vorti], ...
+ op = Analytic(variables = [velo, vorti], ...)
- or
- op = Analytic(variables = {velo: topo, vorti: discr3D}, ...
+ or::
- TODO : write a full explanation for field discretization.
+ op = Analytic(variables = {velo: topo, vorti: discr3D}, ...)
Attributes
----------
variables : dict or list
- fields with their discretisation
- domain : the physical domain on which this operator applies
+ :class:`~hysop.fields.continuous.Continuous`
+ with their discretisation
+ domain : :class:`~hysop.domain.domain.Domain`,
+ the geometry on which this operator applies
"""
# 1 ---- Variables setup ----
@@ -125,7 +126,7 @@ class Operator(object):
if isinstance(self.variables, list):
self.domain = self.variables[0].domain
elif isinstance(self.variables, dict):
- self.domain = self.variables.keys()[0].domain
+ self.domain = self.variables.keys()[0].domain
# Check if all variables have the same domain
for v in self.variables:
@@ -134,7 +135,7 @@ class Operator(object):
# Set/check mpi context
if self._mpis is None:
self._mpis = MPIParams(comm=self.domain.comm_task,
- task_id=self.domain.currentTask())
+ task_id=self.domain.current_task())
# Set profiler
self.profiler = Profiler(self, self.domain.comm_task)
@@ -143,9 +144,9 @@ class Operator(object):
raise RuntimeError("This operator is not defined for the current task")
def wait_for(self, op):
- """
+ """MPI synchronisation
- :param op: Continuous
+ :param op: :class:`~hysop.operator.continuous.Continuous`
Add an operator into 'wait' list of the present object.
It means that before any apply of this operator, all
(mpi) operations in op must be fulfilled, which implies
@@ -155,7 +156,8 @@ class Operator(object):
self._wait_list.append(op)
def wait_list(self):
- """
+ """Get MPI running comm. list
+
Returns
-------
python list of all operators that must be fulfilled
@@ -173,10 +175,11 @@ class Operator(object):
pass
def test_requests(self):
- """
+ """Checks for unfinished mpi communications.
+
Returns
-------
- bool : MPI send/recv test for synchronisation: when this function is
+ bool : MPI send/recv test for synchronisation, when this function is
called, the programm checks if this operator handles some uncomplete
mpi requests (if so return true, else false).
This is a non-blocking call.
@@ -190,7 +193,7 @@ class Operator(object):
ready to apply.
In derived classes, called through @opsetup decorator.
"""
- if not self.domain.currentTask() == self._mpis.task_id:
+ if not self.domain.current_task() == self._mpis.task_id:
self.ontask = False
self._error_()
@@ -225,9 +228,8 @@ class Operator(object):
"""
def is_up(self):
- """
- Returns True if ready to be applied (--> setup function has
- been called succesfully)
+ """Returns True if ready to be applied
+ (--> setup function has been called succesfully)
"""
return self._is_uptodate
diff --git a/HySoP/hysop/operator/discrete/discrete.py b/HySoP/hysop/operator/discrete/discrete.py
index e24b24701eebb8d599532e1fdcc76e997053e807..baff4c8f8d71adf8672c2b866a33d4c72c7c7680 100644
--- a/HySoP/hysop/operator/discrete/discrete.py
+++ b/HySoP/hysop/operator/discrete/discrete.py
@@ -95,9 +95,9 @@ class DiscreteOperator(object):
self._writer = None
# Check topologies consistency
if self.variables is not None:
- topoRef = self.variables[0].topology
+ toporef = self.variables[0].topology
for v in self.variables:
- assert v.topology.isConsistentWith(topoRef)
+ assert v.topology.is_consistent_with(toporef)
def get_work_properties(self):
"""
diff --git a/HySoP/hysop/operator/discrete/drag_and_lift.py b/HySoP/hysop/operator/discrete/drag_and_lift.py
index 2eea997d9a7676f6a53610a6dd42ed1eec51a872..a85de1030d3fbf3c4208af159b2ab7857780c4ac 100644
--- a/HySoP/hysop/operator/discrete/drag_and_lift.py
+++ b/HySoP/hysop/operator/discrete/drag_and_lift.py
@@ -7,7 +7,8 @@ import hysop.tools.numpywrappers as npw
from abc import ABCMeta, abstractmethod
from hysop.numerics.utils import Utils
from hysop.constants import XDIR, YDIR, ZDIR
-import numpy as np
+from hysop.domain.control_box import ControlBox
+from hysop.domain.subsets import Subset
class Forces(DiscreteOperator):
@@ -23,7 +24,7 @@ class Forces(DiscreteOperator):
"""
Parameters
----------
- obstacles : list of :class:`~hysop.domain.subsets.subset.Subset`
+ obstacles : list of :class:`~hysop.domain.subsets.Subset`
list of bodies inside the flow
normalization : double, optional
a normalization coefficient applied to the force, default = 1.
@@ -85,7 +86,7 @@ class Forces(DiscreteOperator):
"""
Parameters
-----------
- obstacles : a list of :class:`hysop.domain.subsets.subset.Subset`
+ obstacles : a list of :class:`hysop.domain.subsets.Subset`
Returns
-------
@@ -186,12 +187,6 @@ class MomentumForces(Forces):
self._buff = None
self._init_buffers()
- # A local discrete operator, to penalize temporary var.
- ## from hysop.operator.discrete.penalization import Penalization
- ## self._penalize = Penalization(variables=self._buff,
- ## obstacles=self._obstacles,
- ## coeff=self._coeff)
-
def _init_buffers(self):
"""
Allocate memory for local buffers (used to compute
@@ -216,10 +211,7 @@ class MomentumForces(Forces):
# mpi communicator
#self._subcomm = obstacles[0].subcomm[self._topology]
# Return the list of indices of points inside the obstacle
- if obst.is_porous:
- return obst.ind[self._topology]
- else:
- return [obst.ind[self._topology]]
+ return obst.ind[self._topology]
def _momentum(self, dt):
# -- Integration over the obstacle --
@@ -233,21 +225,9 @@ class MomentumForces(Forces):
for i in xrange(len(self._indices)):
self._buff[...] = 0.0
icurrent = self._indices[i]
-# print 'shape', np.shape(icurrent), 'r =', self._topology.rank, 'i =', i
- current = self._buff[icurrent]
coeff = self._coeff[i] * dt / (1. + self._coeff[i] * dt)
- current[...] = self.velocity.data[d][icurrent] * coeff
- self.force[d] += npw.real_sum(current)
-# print 'd=', d, 'force =', self.force[d]
-
-# i = 0
-# for ind in self._indices:
-# self._buff[...] = 0.0
-# self._buff[ind] = self.velocity.data[d][ind] / (1.0 + self._coeff[i] * dt)
-# self._buff[ind] *= -1.0
-# self._buff[ind] += self.velocity.data[d][ind]
-# i += 1
-# self.force[d] += npw.real_sum(self._buff[ind])
+ self._buff[icurrent] = self.velocity.data[d][icurrent] * coeff
+ self.force[d] += npw.real_sum(self._buff[icurrent])
self.force *= self._dvol
self.force /= dt
@@ -338,7 +318,7 @@ class NocaForces(Forces):
Parameters
----------
- obstacles: list of :class:`~hysop.domain.subsets.subset.Subset`
+ obstacles: list of :class:`~hysop.domain.subsets.Subset`
obstacles in the flow
Returns
@@ -346,12 +326,10 @@ class NocaForces(Forces):
list of indices of points inside those obstacles
"""
- from hysop.domain.subsets.control_box import ControlBox
assert isinstance(self._voc, ControlBox)
self._on_proc = self._voc.on_proc[self._topology]
msg = 'obstacles arg must be a list.'
assert isinstance(obstacles, list), msg
- from hysop.domain.subsets.subset import Subset
# no obstacle in the box, just for test purpose.
if obstacles is None or len(obstacles) == 0:
return self._voc.ind[self._topology]
diff --git a/HySoP/hysop/operator/discrete/penalization.py b/HySoP/hysop/operator/discrete/penalization.py
index 983faf515316aa379025051f2e237fa016400424..b8113bb13928dea8a1506d74a48fcba598b721a4 100644
--- a/HySoP/hysop/operator/discrete/penalization.py
+++ b/HySoP/hysop/operator/discrete/penalization.py
@@ -1,18 +1,18 @@
# -*- coding: utf-8 -*-
-"""
-@file operator/discrete/penalization.py
-Discrete operator for penalization problem.
+"""Discrete operators for penalization problem.
+.. currentmodule:: hysop.operator.discrete.penalization
+* :class:`Penalization` : standard penalisation
+* :class:`PenalizeVorticity` : vorticity formulation
+
"""
from hysop.constants import debug
from hysop.operator.discrete.discrete import DiscreteOperator
from hysop.tools.profiler import profile
-from hysop.domain.subsets.subset import Subset
+from hysop.domain.subsets import Subset
class Penalization(DiscreteOperator):
- """
- Discretized penalisation operator.
- See details in hysop.operator.penalization
+ """Discretized penalisation operator.
"""
@debug
@@ -25,19 +25,6 @@ class Penalization(DiscreteOperator):
coeff : double, optional
penalization factor applied to all geometries.
- Set :
- - obstacles = {obs1: coeff1, obs2: coeff2, ...} and coeff=None
- to apply a different coefficient on each subset.
- - obstacles = [obs1, obs2, ...], coeff=coeff1 to apply the
- same penalization on all subsets.
- obs1, ob2 ... must be some `~hysop.domain.subsets.Subset`
- and coeff1 must be either a real scalar or a function of the
- coordinates like
- def coeff(*args):
- return 3 * args[0]
-
- with args[0,1,...] = x,y,...
- Warning : coeff as a function is not yet implemented!!
"""
super(Penalization, self).__init__(**kwds)
@@ -123,3 +110,72 @@ class Penalization(DiscreteOperator):
"Simulation parameter is required."
dt = simulation.timeStep
self._apply(dt)
+
+
+class PenalizeVorticity(Penalization):
+ """
+ Discretized penalisation operator.
+ See details in hysop.operator.penalization
+ """
+
+ @debug
+ def __init__(self, vorticity, velocity, curl, **kwds):
+ """
+ Parameters
+ ----------
+ velocity, vorticity: :class:`~hysop.fields.continuous.Field`
+ curl : :class:`~hysop..operator.differential`
+ internal operator to compute the curl of the penalised velocity
+ **kwds : extra parameters for parent class.
+
+ Notes
+ -----
+ velocity is not modified by this operator.
+ vorticity is an in-out parameter.
+ input and ouput variables of the curl are some local buffers.
+ """
+ assert 'variables' not in kwds, 'variables parameter is useless.'
+ super(PenalizeVorticity, self).__init__(variables=[vorticity,
+ velocity],
+ **kwds)
+ self.velocity = velocity
+ self.vorticity = vorticity
+ # warning : a buffer is added for invar variable in curl
+ topo = self.velocity.topology
+ msg = 'Multiresolution not implemented for penalization.'
+ assert self.vorticity.topology == topo, msg
+ self._curl = curl
+
+ def _apply_single_coeff(self, dt):
+ # Vorticity penalization
+ # warning : the buff0 array ensures "invar" to be 0
+ # outside the obstacle for the curl evaluation
+ invar = self._curl.invar
+ nbc = invar.nb_components
+ for d in xrange(nbc):
+ invar.data[d][...] = 0.0
+ coeff = -dt * self._coeff / (1.0 + dt * self._coeff)
+ for d in xrange(nbc):
+ invar.data[d][self._cond] = \
+ self.velocity[d][self._cond] * coeff
+ self._curl.apply()
+ for d in xrange(self.vorticity.nb_components):
+ self.vorticity[d][...] += self._curl.outvar[d][...]
+
+ def _apply_multi_coeff(self, dt):
+ invar = self._curl.invar
+ nbc = invar.nb_components
+
+ for d in xrange(nbc):
+ invar.data[d][...] = 0.0
+
+ for i in xrange(len(self._cond)):
+ coeff = -dt * self._coeff[i] / (1.0 + dt * self._coeff[i])
+ cond = self._cond[i]
+ for d in xrange(nbc):
+ invar.data[d][cond] = self.velocity[d][cond] * coeff
+
+ self._curl.apply()
+
+ for d in xrange(self.vorticity.nb_components):
+ self.vorticity[d][...] += self._curl.outvar[d][...]
diff --git a/HySoP/hysop/operator/discrete/penalize_vorticity.py b/HySoP/hysop/operator/discrete/penalize_vorticity.py
index a33cead92f98479a966852f3da39ad7a95d2f405..279bf74b399ec6dc2d9a0aaccfb6e2847dfed021 100644
--- a/HySoP/hysop/operator/discrete/penalize_vorticity.py
+++ b/HySoP/hysop/operator/discrete/penalize_vorticity.py
@@ -1,75 +1,7 @@
# -*- coding: utf-8 -*-
-"""
-@file operator/discrete/penalization.py
-Discrete operator for penalization problem.
+"""Discrete operators for penalization problem.
"""
from hysop.constants import debug
from hysop.operator.discrete.penalization import Penalization
-class PenalizeVorticity(Penalization):
- """
- Discretized penalisation operator.
- See details in hysop.operator.penalization
- """
-
- @debug
- def __init__(self, vorticity, velocity, curl, **kwds):
- """
- @param[in,out] vorticity : vorticity field
- @param[in] curl : curl operator used to compute
- somebuffer = curl(penalized vorticity)
- """
- assert 'variables' not in kwds, 'variables parameter is useless.'
- super(PenalizeVorticity, self).__init__(variables=[vorticity,
- velocity],
- **kwds)
- self.velocity = velocity
- self.vorticity = vorticity
- # \warning : a buffer is added for invar variable in curl
- topo = self.velocity.topology
- msg = 'Multiresolution not implemented for penalization.'
- assert self.vorticity.topology == topo, msg
- self._curl = curl
-
- def _apply_single_coeff(self, dt):
- # Vorticity penalization
- # \warning : the buff0 array ensures "invar" to be 0
- # outside the obstacle for the curl evaluation
- invar = self._curl.invar
- nbc = invar.nb_components
- for d in xrange(nbc):
- invar.data[d][...] = 0.0
- coeff = -dt * self._coeff / (1.0 + dt * self._coeff)
- for d in xrange(nbc):
- invar.data[d][self._cond] = \
- self.velocity[d][self._cond] * coeff
- self._curl.apply()
- for d in xrange(self.vorticity.nb_components):
- self.vorticity[d][...] += self._curl.outvar[d][...]
- # Velocity penalization (nor used for the moment since
- # it implies a non divergence free velocity field)
-# coeff = 1. / (1. + dt * self._coeff)
-# self._penalize(self.velocity, coeff, self._cond)
-
- def _apply_multi_coeff(self, dt):
- invar = self._curl.invar
- nbc = invar.nb_components
-
- for d in xrange(nbc):
- invar.data[d][...] = 0.0
-
- for i in xrange(len(self._cond)):
- coeff = -dt * self._coeff[i] / (1.0 + dt * self._coeff[i])
- cond = self._cond[i]
- for d in xrange(nbc):
- invar.data[d][cond] = self.velocity[d][cond] * coeff
-
- self._curl.apply()
-
- for d in xrange(self.vorticity.nb_components):
- self.vorticity[d][...] += self._curl.outvar[d][...]
-
- def _penalize(self, field, coef, cond):
- for d in xrange(field.nb_components):
- field[d][cond] *= coef
diff --git a/HySoP/hysop/operator/drag_and_lift.py b/HySoP/hysop/operator/drag_and_lift.py
index a9100f7f59b8a20e0a1f1f660085d6555cefabf5..e8ecb741db2cb59b6017d01b05850c3230c8f5bc 100644
--- a/HySoP/hysop/operator/drag_and_lift.py
+++ b/HySoP/hysop/operator/drag_and_lift.py
@@ -7,6 +7,7 @@
from hysop.operator.computational import Computational
from hysop.operator.continuous import opsetup
from abc import ABCMeta, abstractmethod
+from hysop.domain.control_box import ControlBox
class Forces(Computational):
@@ -201,7 +202,6 @@ class NocaForces(Forces):
if volume_of_control is None:
lr = self.domain.length * 0.9
xr = self.domain.origin + 0.04 * self.domain.length
- from hysop.domain.subsets.control_box import ControlBox
volume_of_control = ControlBox(parent=self.domain,
origin=xr, length=lr)
self._voc = volume_of_control
diff --git a/HySoP/hysop/operator/penalization.py b/HySoP/hysop/operator/penalization.py
index 9d164a41253c8756ccde375c56e5b5e1d1432467..81b451ef8960972617c235ec49c7fa85c8f0b55d 100644
--- a/HySoP/hysop/operator/penalization.py
+++ b/HySoP/hysop/operator/penalization.py
@@ -1,13 +1,27 @@
# -*- coding: utf-8 -*-
-"""
-@file operator/penalization.py
+"""Operators for penalization problem.
+
+.. currentmodule:: hysop.operator.penalization
+
+* :class:`Penalization` : standard penalisation
+* :class:`PenalizeVorticity` : vorticity formulation
+
+See details in :ref:`penalisation` section of HySoP user guide.
-Penalisation of a given field
"""
+
from hysop.operator.computational import Computational
-from hysop.operator.discrete.penalization import Penalization as DiscrPenal
+from hysop.operator.discrete.penalization import Penalization as DPenalV
+from hysop.operator.discrete.penalization import PenalizeVorticity as DPenalW
from hysop.constants import debug
from hysop.operator.continuous import opsetup
+from hysop.domain.subsets import Subset
+import hysop.default_methods as default
+from hysop.methods_keys import SpaceDiscretisation
+from hysop.numerics.finite_differences import FD_C_4,\
+ FD_C_2
+from hysop.operator.differential import Curl
+from hysop.fields.continuous import Field
class Penalization(Computational):
@@ -25,30 +39,45 @@ class Penalization(Computational):
@debug
def __init__(self, obstacles, coeff=None, **kwds):
"""
- @param[in] obstacles : dictionnary or list of
- subsets on which penalization must be applied
- @param[in] coeff : penalization factor or function.
+ Parameters
+ ----------
+ obstacles : dict or list of :class:`~hysop.domain.subsets.Subset`
+ sets of geometries on which penalization must be applied
+ coeff : double, optional
+ penalization factor applied to all geometries.
+ **kwds : extra parameters for parent class
+
+ Notes
+ -----
+ Set::
+
+ obstacles = {obs1: coeff1, obs2: coeff2, ...}
+ coeff = None
- Set :
- - obstacles = {obs1: coeff1, obs2: coeff2, ...} and coeff=None
to apply a different coefficient on each subset.
- - obstacles = [obs1, obs2, ...], coeff=coeff1 to apply the
- same penalization on all subsets.
- obs1, ob2 ... must be some hysop.domain.subsets.Subset
- and coeff1 must be either a scalar or a function of the
- coordinates like
- def coeff(*args):
- return 3 * args[0]
+ Set::
+
+ obstacles = [obs1, obs2, ...]
+ coeff = some_value
+
+ to apply the same penalization on all subsets.
+ obs1, ob2 ... must be some :class:`~hysop.domain.subsets.Subset`
+ and some_value must be either a real scalar or a function of the
+ coordinates like::
+
+ def coeff(*args):
+ return 3 * args[0]
with args[0,1,...] = x,y,...
+ Warning : coeff as a function is not yet implemented!!
"""
super(Penalization, self).__init__(**kwds)
- ## The list of subset on which penalization must be applied
+ # The list of subset on which penalization must be applied
self.obstacles = obstacles
- ## Penalization functions or coef
+ # Penalization functions or coef
self.coeff = coeff
self.input = self.output = self.variables
@@ -57,7 +86,6 @@ class Penalization(Computational):
# all variables must have the same resolution
assert self._single_topo, 'multi-resolution case not allowed.'
topo = self.variables.values()[0]
- from hysop.domain.subsets.subset import Subset
for obs in self.obstacles:
assert isinstance(obs, Subset)
obs.discretize(topo)
@@ -65,8 +93,87 @@ class Penalization(Computational):
@debug
@opsetup
def setup(self, rwork=None, iwork=None):
- self.discrete_op = DiscrPenal(
+ self.discrete_op = DPenalV(
variables=self.discreteFields.values(), obstacles=self.obstacles,
coeff=self.coeff, rwork=rwork, iwork=iwork)
self._is_uptodate = True
+
+
+class PenalizeVorticity(Penalization):
+ """
+ Solve
+ \f{eqnarray*}
+ \frac{\partial w}{\partial t} &=& \lambda\chi_s\nabla\times(v_D - v)
+ \f}
+ using penalization.
+ """
+
+ @debug
+ def __init__(self, velocity, vorticity, **kwds):
+ """
+ Parameters
+ ----------
+ velocity, vorticity: :class:`~hysop.fields.continuous.Field`
+ **kwds : extra parameters for parent class.
+
+ Notes
+ -----
+ velocity is not modified by this operator.
+ vorticity is an in-out parameter.
+ """
+ assert 'variables' not in kwds, 'variables parameter is useless.'
+ super(PenalizeVorticity, self).__init__(
+ variables=[velocity, vorticity], **kwds)
+ # velocity field
+ self.velocity = velocity
+ # vorticity field
+ self.vorticity = vorticity
+ # A method is required to set how the curl will be computed.
+ if self.method is None:
+ self.method = default.DIFFERENTIAL
+ # operator to compute buffer = curl(penalised velocity)
+ self._curl = None
+
+ def discretize(self):
+
+ if self.method[SpaceDiscretisation] is FD_C_4:
+ # Finite differences method
+ # Minimal number of ghost points
+ nb_ghosts = 2
+ elif self.method[SpaceDiscretisation] is FD_C_2:
+ nb_ghosts = 1
+ else:
+ raise ValueError("Unknown method for space discretization of the\
+ differential operator in penalization.")
+ super(PenalizeVorticity, self)._standard_discretize(nb_ghosts)
+ # all variables must have the same resolution
+ assert self._single_topo, 'multi-resolution case not allowed.'
+ topo = self.variables[self.velocity]
+ for obs in self.obstacles:
+ assert isinstance(obs, Subset)
+ obs.discretize(topo)
+ invar = Field(domain=self.velocity.domain,
+ name='curl_in', is_vector=True)
+ dimension = self.domain.dimension
+ outvar = Field(domain=self.velocity.domain,
+ name='curl_out',
+ is_vector=dimension == 3)
+ self._curl = Curl(invar=invar, outvar=outvar,
+ discretization=topo, method=self.method)
+ self._curl.discretize()
+
+ def get_work_properties(self):
+ return self._curl.get_work_properties()
+
+ @debug
+ @opsetup
+ def setup(self, rwork=None, iwork=None):
+ self._curl.setup(rwork, iwork)
+ self.discrete_op = DPenalW(
+ vorticity=self.discreteFields[self.vorticity],
+ velocity=self.discreteFields[self.velocity],
+ curl=self._curl.discrete_op,
+ obstacles=self.obstacles,
+ coeff=self.coeff, rwork=rwork, iwork=iwork)
+ self._is_uptodate = True
diff --git a/HySoP/hysop/operator/penalize_vorticity.py b/HySoP/hysop/operator/penalize_vorticity.py
index 9d1ab278af5d67543173b01968a16b073ba879c2..adaa33a24b39ea70487911c837643cbdd9d09645 100644
--- a/HySoP/hysop/operator/penalize_vorticity.py
+++ b/HySoP/hysop/operator/penalize_vorticity.py
@@ -1,90 +1,3 @@
# -*- coding: utf-8 -*-
+"""Compute the vorticity field from velocity using penalization.
"""
-@file operator/penalization_and_curl.py
-
-Compute the vorticity field from velocity using penalization.
-"""
-from hysop.operator.penalization import Penalization
-from hysop.operator.discrete.penalize_vorticity\
- import PenalizeVorticity as PW
-from hysop.constants import debug
-from hysop.operator.continuous import opsetup
-import hysop.default_methods as default
-from hysop.methods_keys import SpaceDiscretisation
-from hysop.numerics.finite_differences import FD_C_4,\
- FD_C_2
-from hysop.operator.differential import Curl
-from hysop.fields.continuous import Field
-
-
-class PenalizeVorticity(Penalization):
- """
- Solve
- \f{eqnarray*}
- \frac{\partial w}{\partial t} &=& \lambda\chi_s\nabla\times(v_D - v)
- \f}
- using penalization.
- """
-
- @debug
- def __init__(self, velocity, vorticity, **kwds):
- """
- @param[in,out] velocity : velocity field
- @param[in,out] vorticity : vorticity field
- """
- assert 'variables' not in kwds, 'variables parameter is useless.'
- super(PenalizeVorticity, self).__init__(
- variables=[velocity, vorticity], **kwds)
- # velocity field
- self.velocity = velocity
- # vorticity field
- self.vorticity = vorticity
- # A method is required to set how the curl will be computed.
- if self.method is None:
- self.method = default.DIFFERENTIAL
- # operator to compute buffer = curl(penalised velocity)
- self._curl = None
-
- def discretize(self):
-
- if self.method[SpaceDiscretisation] is FD_C_4:
- # Finite differences method
- # Minimal number of ghost points
- nb_ghosts = 2
- elif self.method[SpaceDiscretisation] is FD_C_2:
- nb_ghosts = 1
- else:
- raise ValueError("Unknown method for space discretization of the\
- differential operator in penalization.")
- super(PenalizeVorticity, self)._standard_discretize(nb_ghosts)
- # all variables must have the same resolution
- assert self._single_topo, 'multi-resolution case not allowed.'
- topo = self.variables[self.velocity]
- from hysop.domain.subsets.subset import Subset
- for obs in self.obstacles:
- assert isinstance(obs, Subset)
- obs.discretize(topo)
- invar = Field(domain=self.velocity.domain,
- name='curl_in', is_vector=True)
- dimension = self.domain.dimension
- outvar = Field(domain=self.velocity.domain,
- name='curl_out',
- is_vector=dimension == 3)
- self._curl = Curl(invar=invar, outvar=outvar,
- discretization=topo, method=self.method)
- self._curl.discretize()
-
- def get_work_properties(self):
- return self._curl.get_work_properties()
-
- @debug
- @opsetup
- def setup(self, rwork=None, iwork=None):
- self._curl.setup(rwork, iwork)
- self.discrete_op = PW(
- vorticity=self.discreteFields[self.vorticity],
- velocity=self.discreteFields[self.velocity],
- curl=self._curl.discrete_op,
- obstacles=self.obstacles,
- coeff=self.coeff, rwork=rwork, iwork=iwork)
- self._is_uptodate = True
diff --git a/HySoP/hysop/operator/poisson.py b/HySoP/hysop/operator/poisson.py
index 73d29b5351b5080641d2a63d42d161db72b59b80..5151d2380bb00321dd18aac451c85f0322a21b56 100644
--- a/HySoP/hysop/operator/poisson.py
+++ b/HySoP/hysop/operator/poisson.py
@@ -1,8 +1,5 @@
# -*- coding: utf-8 -*-
-"""
-@file poisson.py
-Poisson problem.
-
+"""Poisson problem.
"""
from hysop.operator.computational import Computational
from hysop.operator.discrete.poisson_fft import PoissonFFT
@@ -16,7 +13,6 @@ import hysop.default_methods as default
class Poisson(Computational):
"""
- Poisson problem,
\f{eqnarray*}
v = Op(\omega)
\f} with :
diff --git a/HySoP/hysop/operator/redistribute.py b/HySoP/hysop/operator/redistribute.py
index c66846c7b230019657a0a2ef2115cd255ebecec2..c37db1da90577dc2fc684637dba68a281f12101e 100644
--- a/HySoP/hysop/operator/redistribute.py
+++ b/HySoP/hysop/operator/redistribute.py
@@ -177,7 +177,7 @@ class Redistribute(Operator):
else:
msg = "the source/target is neither an operator or a topology."
raise AttributeError(msg)
- assert result.task_id() == self.domain.currentTask()
+ assert result.task_id() == self.domain.current_task()
return result
def printComputeTime(self):
diff --git a/HySoP/hysop/operator/redistribute_inter.py b/HySoP/hysop/operator/redistribute_inter.py
index 96fc24521c4edae94d76eacb3ee3b2ea89e16445..f72bdb7c7d76fba564cafcbf151528cbd830f0b2 100644
--- a/HySoP/hysop/operator/redistribute_inter.py
+++ b/HySoP/hysop/operator/redistribute_inter.py
@@ -53,7 +53,7 @@ class RedistributeInter(Redistribute):
self._set_domain_and_tasks()
# Domain is set, we can check if we are on source or target
- current_task = self.domain.currentTask()
+ current_task = self.domain.current_task()
self._is_source = current_task == self._source_id
self._is_target = current_task == self._target_id
assert self._is_target or self._is_source
diff --git a/HySoP/hysop/operator/tests/ref_files/penal_vort_2d_multi_sphere_00000.h5 b/HySoP/hysop/operator/tests/ref_files/penal_vort_2d_multi_sphere_00000.h5
index 6878cadc01a867c59511df698ffc5a9e5eebaec7..1c91c8933ca496944febab6633f8a4678c5c1e49 100644
Binary files a/HySoP/hysop/operator/tests/ref_files/penal_vort_2d_multi_sphere_00000.h5 and b/HySoP/hysop/operator/tests/ref_files/penal_vort_2d_multi_sphere_00000.h5 differ
diff --git a/HySoP/hysop/operator/tests/ref_files/penal_vort_3d_multi_sphere_00000.h5 b/HySoP/hysop/operator/tests/ref_files/penal_vort_3d_multi_sphere_00000.h5
index 6e261313f258d3197c0ab79c5db6ba3215670f19..56b96e09bfa8827e24d40fa5f0df2aa3d5199657 100644
Binary files a/HySoP/hysop/operator/tests/ref_files/penal_vort_3d_multi_sphere_00000.h5 and b/HySoP/hysop/operator/tests/ref_files/penal_vort_3d_multi_sphere_00000.h5 differ
diff --git a/HySoP/hysop/operator/tests/test_adaptive_time_step.py b/HySoP/hysop/operator/tests/test_adaptive_time_step.py
index d31377a96e2d2abcf4021404c8dc8f02006751e0..6792f4b3a49fa79a2010890c5b90b4f106306896 100644
--- a/HySoP/hysop/operator/tests/test_adaptive_time_step.py
+++ b/HySoP/hysop/operator/tests/test_adaptive_time_step.py
@@ -119,18 +119,18 @@ def test_adapt_4():
discretization=d3d, lcfl=0.125, cfl=0.5,
mpi_params=cpu_task)
simu.initialize()
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
op.discretize()
op.setup()
vorti.initialize()
while not simu.isOver:
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
op.apply()
op.wait()
simu.advance()
refval = 0
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
refval = simu.timeStep
refval = main_comm.bcast(refval, root=0)
assert refval == simu.timeStep
diff --git a/HySoP/hysop/operator/tests/test_drag_and_lift.py b/HySoP/hysop/operator/tests/test_drag_and_lift.py
index 50fcb6e69d08122ebef15c2c2985c12680e966e0..27282aa9ff18006a7a8b89f022cc349f8e48b6d7 100644
--- a/HySoP/hysop/operator/tests/test_drag_and_lift.py
+++ b/HySoP/hysop/operator/tests/test_drag_and_lift.py
@@ -1,6 +1,6 @@
# -*- coding: utf-8 -*-
-from hysop.domain.subsets.sphere import Sphere
-from hysop.operator.penalize_vorticity import PenalizeVorticity
+from hysop.domain.subsets import Sphere
+from hysop.operator.penalization import PenalizeVorticity
from hysop.problem.simulation import Simulation
from hysop.tools.parameters import Discretization, IOParams
from hysop.mpi.topology import Cartesian
@@ -162,6 +162,6 @@ def test_noca3():
# This may be useful to run mpi tests
if __name__ == "__main__":
test_momentum_forces_3d()
- test_noca1()
- test_noca2()
+ #test_noca1()
+ #test_noca2()
#test_noca3()
diff --git a/HySoP/hysop/operator/tests/test_energy_enstrophy.py b/HySoP/hysop/operator/tests/test_energy_enstrophy.py
index 17125d6b0050ea9cda3ea52a608b8844a60231de..124466d2694b8c8980e758352d7d991872552859 100644
--- a/HySoP/hysop/operator/tests/test_energy_enstrophy.py
+++ b/HySoP/hysop/operator/tests/test_energy_enstrophy.py
@@ -3,11 +3,22 @@ import hysop as pp
from hysop.operator.energy_enstrophy import EnergyEnstrophy
from hysop.problem.simulation import Simulation
from hysop.tools.parameters import Discretization
-from hysop import VariableParameter, Field
+from hysop import Field
import numpy as np
-import hysop.tools.numpywrappers as npw
-import os
-from scipy.integrate import nquad
+try:
+ from scipy.integrate import nquad
+
+except:
+ def nquad(func, coords_range):
+ coords = []
+ nbsteps = 1000
+ for x in coords_range:
+ coords.append(np.linspace(x[0], x[1], nbsteps))
+ coords = tuple(coords)
+ hstep = coords[0][1] - coords[0][0]
+ ll = [coords_range[i][1] - coords_range[i][0] for i in xrange(1, 3)]
+ return [np.sum(func(*coords)[:-1]) * hstep * np.prod(ll)]
+
sin = np.sin
cos = np.cos
@@ -58,19 +69,19 @@ def test_energy_enstrophy():
op.apply(simu)
intrange = []
box = topo.domain
- invvol = 1./ np.prod(box.length)
+ invvol = 1. / np.prod(box.length)
for i in xrange(dim):
origin = box.origin[i]
end = origin + box.length[i]
intrange.append([origin, end])
intrange = 2 * intrange
- Eref = nquad(energy_ref, intrange[:dim])[0]
- Eref += nquad(energy_ref, intrange[1:dim + 1])[0]
- Eref += nquad(energy_ref, intrange[2:dim + 2])[0]
- Eref *= invvol
+ eref = nquad(energy_ref, intrange[:dim])[0]
+ eref += nquad(energy_ref, intrange[1:dim + 1])[0]
+ eref += nquad(energy_ref, intrange[2:dim + 2])[0]
+ eref *= invvol
tol = (topo.mesh.space_step).max() ** 2
- assert (op.energy() - Eref * 0.5) < tol
- assert (op.enstrophy() - Eref) < tol
+ assert (op.energy() - eref * 0.5) < tol
+ assert (op.enstrophy() - eref) < tol
if __name__ == "__main__":
test_energy_enstrophy()
diff --git a/HySoP/hysop/operator/tests/test_hdf5_io.py b/HySoP/hysop/operator/tests/test_hdf5_io.py
index 458f15f07e0ac299637f68fe1d7c49f10ba78341..c1c0f5ab7bd4ecc68e4bdbd30da19e570288ccbd 100644
--- a/HySoP/hysop/operator/tests/test_hdf5_io.py
+++ b/HySoP/hysop/operator/tests/test_hdf5_io.py
@@ -13,6 +13,7 @@ from hysop.tools.parameters import Discretization, IOParams
from hysop.operator.hdf_io import HDF_Writer, HDF_Reader
from hysop.tools.io_utils import IO
from hysop.mpi import main_rank, main_size
+from hysop.domain.subsets import SubBox
Lx = 2.
nb = 65
@@ -286,7 +287,7 @@ def test_write_read_subset_1():
velo = Field(domain=dom, formula=vec3D, name='velo', is_vector=True)
# A subset of the current domain
- from hysop.domain.subsets.boxes import SubBox
+ from hysop.domain.subsets import SubBox
mybox = SubBox(origin=[-0.5, 2.3, 4.1], length=[Lx / 2, Lx / 3, Lx],
parent=dom)
# Write a vector field, using default for output location
@@ -339,7 +340,6 @@ def test_write_read_subset_2():
velo = Field(domain=dom, formula=vec3D, name='velo', is_vector=True)
# A subset of the current domain
- from hysop.domain.subsets.boxes import SubBox
# a plane ...
mybox = SubBox(origin=[-0.5, 2.3, 4.1], length=[Lx / 2, Lx / 3, 0.0],
parent=dom)
diff --git a/HySoP/hysop/operator/tests/test_operators.py b/HySoP/hysop/operator/tests/test_operators.py
index 2a86d6d45b28f7e03d71093803f603f47e43adf8..1b184f29bd99d3faa2e2df79e48e613ebdbdb16e 100644
--- a/HySoP/hysop/operator/tests/test_operators.py
+++ b/HySoP/hysop/operator/tests/test_operators.py
@@ -25,18 +25,18 @@ def test_tasks():
test proper tasks assignment
"""
dom, topo = create_multitask_context(dim=3, discr=r_ng)
- assert topo.task_id() == dom.currentTask()
+ assert topo.task_id() == dom.current_task()
velo = pp.Field(domain=dom, name='velocity',
is_vector=True, formula=v3d)
op = Analytic(variables={velo: r_ng})
op.discretize()
op.setup()
- assert op.task_id() == dom.currentTask()
+ assert op.task_id() == dom.current_task()
if main_size == 8:
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
assert op.variables[velo].size == 5
- elif dom.isOnTask(GPU):
+ elif dom.is_on_task(GPU):
assert op.variables[velo].size == 2
- elif dom.isOnTask(OTHER):
+ elif dom.is_on_task(OTHER):
assert op.variables[velo].size == 1
diff --git a/HySoP/hysop/operator/tests/test_penalization.py b/HySoP/hysop/operator/tests/test_penalization.py
index 50e5915585f1ea65902153f7f0409721b3c403c5..b2ae2276220efacbd684fce33c4674942098f992 100644
--- a/HySoP/hysop/operator/tests/test_penalization.py
+++ b/HySoP/hysop/operator/tests/test_penalization.py
@@ -1,8 +1,7 @@
# -*- coding: utf-8 -*-
from hysop.domain.subsets import HemiSphere, Sphere, Cylinder
from hysop.domain.porous import Porous
-from hysop.operator.penalization import Penalization
-from hysop.operator.penalize_vorticity import PenalizeVorticity
+from hysop.operator.penalization import Penalization, PenalizeVorticity
from hysop.problem.simulation import Simulation
from hysop.tools.parameters import Discretization, IOParams
from hysop.mpi.topology import Cartesian
@@ -11,6 +10,7 @@ import numpy as np
import os
from hysop import Field, Box
from hysop.operator.hdf_io import HDF_Reader
+from hysop.domain.subsets import SubBox
def v2d(res, x, y, t):
@@ -51,7 +51,6 @@ ldom = npw.asrealarray([math.pi * 2., ] * 3)
xdef = npw.asrealarray(xdom + 0.2)
xpos = npw.asrealarray(ldom * 0.5)
xpos[-1] += 0.1
-from hysop.mpi import main_size
working_dir = os.getcwd() + '/'
@@ -128,7 +127,6 @@ def test_penal_2d_multi():
radius=rd + 0.3)
ll = topo.domain.length.copy()
ll[1] = 0.
- from hysop.domain.subsets.boxes import SubBox
downplane = SubBox(parent=topo.domain, origin=topo.domain.origin,
length=ll)
penal = Penalization(variables=[scal, velo], discretization=topo,
@@ -170,7 +168,6 @@ def test_penal_3d_multi():
radius=rd + 0.3)
ll = topo.domain.length.copy()
ll[1] = 0.
- from hysop.domain.subsets.boxes import SubBox
downplane = SubBox(parent=topo.domain, origin=topo.domain.origin,
length=ll)
penal = Penalization(variables=[scal, velo], discretization=topo,
@@ -192,7 +189,6 @@ def test_penal_3d_porous():
source=Sphere, layers=[0.5, 0.7, 0.3])
ll = topo.domain.length.copy()
ll[1] = 0.
- from hysop.domain.subsets.boxes import SubBox
downplane = SubBox(parent=topo.domain, origin=topo.domain.origin,
length=ll)
penal = Penalization(variables=[scal, velo], discretization=topo,
@@ -214,7 +210,6 @@ def test_penal_3d_porous_cyl():
source=Cylinder, layers=[0.5, 0.7, 0.3])
ll = topo.domain.length.copy()
ll[1] = 0.
- from hysop.domain.subsets.boxes import SubBox
downplane = SubBox(parent=topo.domain, origin=topo.domain.origin,
length=ll)
penal = Penalization(variables=[scal, velo], discretization=topo,
@@ -236,7 +231,6 @@ def test_penal_2d_porous():
source=Sphere, layers=[0.5, 0.7, 0.3])
ll = topo.domain.length.copy()
ll[1] = 0.
- from hysop.domain.subsets.boxes import SubBox
downplane = SubBox(parent=topo.domain, origin=topo.domain.origin,
length=ll)
penal = Penalization(variables=[scal, velo], discretization=topo,
diff --git a/HySoP/hysop/operator/tests/test_poisson.py b/HySoP/hysop/operator/tests/test_poisson.py
index 7676e87dbfb32d91466dae6b7a4b098ff98ea1e2..30bc1106e39a49b4014d08e90d5a07b7b041f4ee 100755
--- a/HySoP/hysop/operator/tests/test_poisson.py
+++ b/HySoP/hysop/operator/tests/test_poisson.py
@@ -6,12 +6,13 @@ from hysop.operator.analytic import Analytic
from hysop.operator.reprojection import Reprojection
from hysop.problem.simulation import Simulation
from hysop.tools.parameters import Discretization
-from hysop.methods_keys import Scales, TimeIntegrator, Interpolation,\
- Remesh, Support, Splitting, dtCrit, SpaceDiscretisation, \
+from hysop.methods_keys import SpaceDiscretisation, \
GhostUpdate, Formulation
import numpy as np
import hysop.tools.numpywrappers as npw
import math
+from hysop.domain.subsets import SubBox
+
pi = math.pi
sin = np.sin
cos = np.cos
@@ -217,7 +218,6 @@ def test_Poisson3D_correction():
refOp.apply(simu)
refD = ref.discretize(topo)
vd = velocity.discretize(topo)
- from hysop.domain.subsets.boxes import SubBox
surf = SubBox(parent=dom, origin=dom.origin,
length=[0., LL[1], LL[2]])
surf.discretize(topo)
diff --git a/HySoP/hysop/operator/tests/test_redistribute.py b/HySoP/hysop/operator/tests/test_redistribute.py
index d5a15c8490dda0378a585d71a1c953e35efb0559..e72ce82681a10cf3564fd68c160b83ed96c8a908 100644
--- a/HySoP/hysop/operator/tests/test_redistribute.py
+++ b/HySoP/hysop/operator/tests/test_redistribute.py
@@ -14,6 +14,7 @@ from hysop.mpi.main_var import main_comm, main_size
from hysop.mpi.tests.utils import create_inter_topos, CPU, GPU, OTHER,\
create_multitask_context
from hysop.fields.tests.func_for_tests import v3d, v2d, v3dbis
+from hysop.operator.poisson import Poisson
dim3 = 3
@@ -413,7 +414,7 @@ def test_distribute_inter():
velocity = fields['velocity']
# Inititialize velocity on CPU task
vd = velocity.discretize(topo=topo1)
- if dom_tasks.isOnTask(CPU):
+ if dom_tasks.is_on_task(CPU):
velocity.initialize(time=simu.time, topo=topo1)
# A field to compute a reference solution, initialized with an analytic
@@ -425,10 +426,10 @@ def test_distribute_inter():
op.apply(simu)
wnorm = reference.norm(topo1)
vnorm = velocity.norm(topo1)
- if dom_tasks.isOnTask(CPU):
+ if dom_tasks.is_on_task(CPU):
assert (vnorm > 0).all()
assert np.allclose(vnorm, wnorm)
- elif dom_tasks.isOnTask(GPU):
+ elif dom_tasks.is_on_task(GPU):
assert (wnorm > 0).all()
assert np.allclose(vnorm, 0)
@@ -440,10 +441,10 @@ def test_distribute_inter():
red.apply(simu)
red.wait()
wd = reference.discretize(topo1)
- if dom.isOnTask(CPU):
+ if dom.is_on_task(CPU):
assert (vnorm > 0).all()
assert np.allclose(vnorm, wnorm)
- elif dom.isOnTask(GPU):
+ elif dom.is_on_task(GPU):
assert (wnorm > 0).all()
assert np.allclose(vnorm, wnorm)
for d in xrange(dom.dimension):
@@ -461,19 +462,19 @@ def test_distribute_inter_2():
if main_size > 2:
proc_tasks[-1] = GPU
proc_tasks[0] = GPU
- domTasks = pp.Box(proc_tasks=proc_tasks)
- fields, simu = init_3d(domTasks)
+ domtasks = pp.Box(proc_tasks=proc_tasks)
+ fields, simu = init_3d(domtasks)
velocity = fields['velocity']
# Inititialize velocity on GPU task
- if domTasks.isOnTask(GPU):
- topo_GPU = domTasks.create_topology(r_ng)
+ if domtasks.is_on_task(GPU):
+ topo_GPU = domtasks.create_topology(r_ng)
vd = velocity.discretize(topo=topo_GPU)
velocity.initialize(time=simu.time, topo=topo_GPU)
vnorm = velocity.norm(topo_GPU)
assert (vnorm > 0).all()
topo_CPU = None
- elif domTasks.isOnTask(CPU):
+ elif domtasks.is_on_task(CPU):
# A field to compute a reference solution, initialized with an analytic
# operator, on both tasks.
reference = fields['vorticity']
@@ -491,7 +492,7 @@ def test_distribute_inter_2():
red.setup()
red.apply(simu)
red.wait()
- if domTasks.isOnTask(CPU):
+ if domtasks.is_on_task(CPU):
vd = velocity.discretize(topo=topo_CPU)
wd = reference.discretize(topo=topo_CPU)
vnorm = velocity.norm(topo_CPU)
@@ -512,7 +513,7 @@ def test_distribute_inter_3():
fields, simu = init_3d(dom_tasks)
velocity = fields['velocity']
# Inititialize velocity on GPU task
- if dom_tasks.isOnTask(GPU):
+ if dom_tasks.is_on_task(GPU):
vd = velocity.discretize(topo=topo1)
velocity.initialize(time=simu.time, topo=topo1)
@@ -524,7 +525,7 @@ def test_distribute_inter_3():
op.setup()
op.apply(simu)
wnorm = reference.norm(topo2)
- if dom_tasks.isOnTask(GPU):
+ if dom_tasks.is_on_task(GPU):
vnorm = velocity.norm(topo1)
assert (vnorm > 0).all()
assert np.allclose(vnorm, wnorm)
@@ -535,7 +536,7 @@ def test_distribute_inter_3():
red.setup()
red.apply(simu)
red.wait()
- if dom_tasks.isOnTask(CPU):
+ if dom_tasks.is_on_task(CPU):
wd = reference.discretize(topo=topo2)
vd = velocity.discretize(topo=topo2)
vnorm = velocity.norm(topo2)
@@ -551,11 +552,11 @@ def test_distribute_inter_4():
"""
if main_size < 4:
return
- dom_tasks, topo= create_multitask_context(3, r_ng)
+ dom_tasks, topo = create_multitask_context(3, r_ng)
fields, simu = init_3d(dom_tasks)
velocity = fields['velocity']
# Inititialize velocity on GPU task
- if dom_tasks.isOnTask(GPU):
+ if dom_tasks.is_on_task(GPU):
vd = velocity.discretize(topo=topo)
velocity.initialize(time=simu.time, topo=topo)
@@ -568,7 +569,7 @@ def test_distribute_inter_4():
op.apply(simu)
# Redistribute from topo on CPU to topo on GPU, ignoring OTHER
- if not dom_tasks.isOnTask(OTHER):
+ if not dom_tasks.is_on_task(OTHER):
red = RedistributeInter(source=topo, target=topo, parent=main_comm,
variables=[velocity],
source_id=GPU, target_id=CPU)
@@ -576,14 +577,14 @@ def test_distribute_inter_4():
red.apply(simu)
red.wait()
- if dom_tasks.isOnTask(CPU):
+ if dom_tasks.is_on_task(CPU):
wd = reference.discretize(topo=topo)
vd = velocity.discretize(topo=topo)
ind = topo.mesh.iCompute
for d in xrange(dom.dimension):
assert np.allclose(wd.data[d][ind], vd.data[d][ind])
- if dom_tasks.isOnTask(OTHER):
+ if dom_tasks.is_on_task(OTHER):
assert topo not in velocity.discreteFields
@@ -591,31 +592,30 @@ def test_distribute_inter_5():
"""
2 tasks, redistribute op to op
"""
- from hysop.operator.poisson import Poisson
if main_size < 4:
return
proc_tasks = [CPU, ] * main_size
proc_tasks[-1] = GPU
proc_tasks[0] = GPU
- domTasks = pp.Box(proc_tasks=proc_tasks)
+ domtasks = pp.Box(proc_tasks=proc_tasks)
- fields, simu = init_3d(domTasks)
+ fields, simu = init_3d(domtasks)
velocity = fields['velocity']
reference = fields['vorticity']
- if domTasks.isOnTask(CPU):
+ if domtasks.is_on_task(CPU):
# initialize velocity on CPU
op = Analytic(variables={velocity: r_ng})
op.discretize()
op.setup()
op.apply(simu)
- elif domTasks.isOnTask(GPU):
+ elif domtasks.is_on_task(GPU):
# initialize reference on CPU
op_init = Analytic(variables={reference: r_ng})
op_init.discretize()
op_init.setup()
op_init.apply(simu)
# An empty operator for velocity
- op = Poisson(velocity=velocity, vorticity=reference,
+ op = Poisson(output_field=velocity, input_field=reference,
discretization=r_ng)
op.discretize()
op.setup()
@@ -628,11 +628,11 @@ def test_distribute_inter_5():
red.apply(simu)
red.wait()
- if domTasks.isOnTask(GPU):
+ if domtasks.is_on_task(GPU):
toporef = op.variables[reference]
vd = velocity.discretize(toporef)
wd = velocity.discretize(toporef)
- for d in xrange(dom.dimension):
+ for d in xrange(domtasks.dimension):
assert np.allclose(wd.data[d], vd.data[d])
@@ -640,34 +640,33 @@ def test_distribute_inter_2d():
"""
2 tasks, redistribute op to op, 2D domain
"""
- from hysop.operator.poisson import Poisson
if main_size < 4:
return
proc_tasks = [CPU, ] * main_size
proc_tasks[-1] = GPU
proc_tasks[0] = GPU
- domTasks = pp.Box(dimension=2, proc_tasks=proc_tasks)
- velocity = pp.Field(domain=domTasks, name='velocity',
+ domtasks = pp.Box(dimension=2, proc_tasks=proc_tasks)
+ velocity = pp.Field(domain=domtasks, name='velocity',
is_vector=True, formula=v2d)
- vort = pp.Field(domain=domTasks, name='vort')
+ vort = pp.Field(domain=domtasks, name='vort')
simu = Simulation(nbIter=3)
- reference = pp.Field(domain=domTasks, name='ref',
+ reference = pp.Field(domain=domtasks, name='ref',
is_vector=True, formula=v2d)
r_2d = Discretization([33, ] * 2)
- if domTasks.isOnTask(CPU):
+ if domtasks.is_on_task(CPU):
# initialize velocity on CPU
op = Analytic(variables={velocity: r_2d})
op.discretize()
op.setup()
op.apply(simu)
- elif domTasks.isOnTask(GPU):
+ elif domtasks.is_on_task(GPU):
# initialize reference on CPU
op_init = Analytic(variables={reference: r_2d})
op_init.discretize()
op_init.setup()
op_init.apply(simu)
# An empty operator for velocity
- op = Poisson(velocity=velocity, vorticity=vort,
+ op = Poisson(output_field=velocity, input_field=vort,
discretization=r_2d)
op.discretize()
op.setup()
@@ -680,7 +679,7 @@ def test_distribute_inter_2d():
red.apply(simu)
red.wait()
- if domTasks.isOnTask(GPU):
+ if domtasks.is_on_task(GPU):
toporef = op.variables[velocity]
vd = velocity.discretize(toporef)
wd = velocity.discretize(toporef)
diff --git a/HySoP/hysop/operator/velocity_correction.py b/HySoP/hysop/operator/velocity_correction.py
index aa1eb885db8bef0c6b1ef35384f920d166827a29..088eabd876bf26f84df9e8304f84fba237feb5a2 100755
--- a/HySoP/hysop/operator/velocity_correction.py
+++ b/HySoP/hysop/operator/velocity_correction.py
@@ -8,7 +8,7 @@ Operator to shift velocity to fit with a required input flowrate.
from hysop.constants import debug
from hysop.operator.discrete.velocity_correction import VelocityCorrection_D
from hysop.operator.computational import Computational
-from hysop.domain.subsets.control_box import ControlBox
+from hysop.domain.control_box import ControlBox
from hysop.operator.continuous import opsetup
diff --git a/HySoP/hysop/tools/io_utils.py b/HySoP/hysop/tools/io_utils.py
index e838752cb6eaf8bea3663c7c68dcc308da11d9c3..f6c77076c4e66184e33834f7000c1ffc6828ae4d 100644
--- a/HySoP/hysop/tools/io_utils.py
+++ b/HySoP/hysop/tools/io_utils.py
@@ -286,7 +286,7 @@ class XMF(object):
Parameters
-----------
- topo : hysop.mpi.topology.Cartesian
+ topo : :class:`hysop.mpi.topology.Cartesian`
used as reference to define local and global meshes in xdmf file.
dataset_names : list
all datasets names
@@ -296,7 +296,7 @@ class XMF(object):
current time
filename : string
name of the hdf file which contains datas for the current process.
- subset : hysop.domain.subsets.subset.Subset, optional
+ subset : :class:`hysop.domain.subsets.Subset`, optional
to define a grid only on this subset.
If None, grid on the whole domain (from topo)
diff --git a/HySoP/hysop/tools/parameters.py b/HySoP/hysop/tools/parameters.py
index 2cc89f445ccbaf120496ecf9b1ca3023cf1bb6d2..860c639c721c61a3267f799d5a42b7f7b4795fc1 100644
--- a/HySoP/hysop/tools/parameters.py
+++ b/HySoP/hysop/tools/parameters.py
@@ -10,7 +10,7 @@ from hysop.constants import DEFAULT_TASK_ID
class MPIParams(namedtuple('MPIParams', ['comm', 'task_id',
- 'rank', 'on_task'])):
+ 'rank', 'on_task'])):
"""
Struct to save mpi parameters :
- comm : parent mpi communicator (default = main_comm)
@@ -27,10 +27,10 @@ class MPIParams(namedtuple('MPIParams', ['comm', 'task_id',
---------
op = SomeOperator(..., task_id=1)
- if op.isOnTask():
+ if op.is_on_task():
...
- 'isOnTask' will return MPIParams.on_task value for op
+ 'is_on_task' will return MPIParams.on_task value for op
and tell if the current operator belongs to the current process
mpi task.
"""
@@ -41,7 +41,7 @@ class MPIParams(namedtuple('MPIParams', ['comm', 'task_id',
else:
rank = MPI.UNDEFINED
return super(MPIParams, cls).__new__(cls, comm, task_id,
- rank, on_task)
+ rank, on_task)
import hysop.tools.numpywrappers as npw
diff --git a/HySoP/hysop/tools/remeshing_formula_parsing.py b/HySoP/hysop/tools/remeshing_formula_parsing.py
index a3616178bb0397ecbfe90a55fabaf722ed75a6f7..cb158d9334b5716fe940b00efb8e2f9adef14e21 100644
--- a/HySoP/hysop/tools/remeshing_formula_parsing.py
+++ b/HySoP/hysop/tools/remeshing_formula_parsing.py
@@ -7,108 +7,115 @@ hysop.numerics.remeshing module or in OpenCL code in
hysop/gpu/cl_src/remeshing/weights*
"""
import re
-import sympy as sp
+try:
+ import sympy as sp
-# Weights names
-weights_names = ['alpha',
- 'beta',
- 'gamma',
- 'delta',
- 'eta',
- 'zeta',
- 'theta',
- 'iota',
- 'kappa',
- 'mu']
+ # Weights names
+ weights_names = ['alpha',
+ 'beta',
+ 'gamma',
+ 'delta',
+ 'eta',
+ 'zeta',
+ 'theta',
+ 'iota',
+ 'kappa',
+ 'mu']
-
-def parse(f, fac=1, vec=False, toOpenCL=True, CLBuiltins=False, keep=False):
- """
- Parsing function.
- @param f : functions to parse as string
- @param fac : numeric factor for all formulas
- @param vec : is OpenCL output is generated to use vector builtin types
- @param toOpenCL : is OpenCL output
- @param CLBuiltins : is OpenCL output uses fma builtin function
- @param keep : low parsing
- """
- assert not (vec and not toOpenCL), "Vector works only in OpenCL parsing"
- assert not (CLBuiltins and not toOpenCL),\
- "CLBuiltins only in OpenCL parsing"
- t = "float__N__" if vec else "float"
- cteEnd = ".0" if toOpenCL else "."
- res = ""
- # Split each line
- fl = f.split('\n')
- # sympy formulas
- y = sp.symbols('y')
- print (y)
- sw = [None] * f.count(';')
- i = 0
- for wl in fl:
- if len(wl) > 2:
- # replace pow
- power = re.search('pow\(y, ([0-9]+)\)', wl)
- if power is not None:
- np = "y" + "*y" * (int(power.group(1)) - 1)
- wl = wl.replace(power.group(0), np)
- sw[i] = '('
- sw[i] += str(sp.horner(eval(wl.split(';')[0].split('=')[1]) * fac))
- sw[i] += ')/' + str(fac)
- i += 1
- for i, s in enumerate(sw):
- if not keep:
- if toOpenCL:
- res += "inline " + t + " "
- res += weights_names[i] + "(" + t + " y){\n"
- res += ' return '
+ def parse(f, fac=1, vec=False, toOpenCL=True,
+ CLBuiltins=False, keep=False):
+ """
+ Parsing function.
+ @param f : functions to parse as string
+ @param fac : numeric factor for all formulas
+ @param vec : is OpenCL output is generated to use vector builtin types
+ @param toOpenCL : is OpenCL output
+ @param CLBuiltins : is OpenCL output uses fma builtin function
+ @param keep : low parsing
+ """
+ msg = 'Vector works only in OpenCL parsing'
+ assert not (vec and not toOpenCL), msg
+ assert not (CLBuiltins and not toOpenCL),\
+ "CLBuiltins only in OpenCL parsing"
+ t = "float__N__" if vec else "float"
+ cteEnd = ".0" if toOpenCL else "."
+ res = ""
+ # Split each line
+ fl = f.split('\n')
+ # sympy formulas
+ y = sp.symbols('y')
+ print (y)
+ sw = [None] * f.count(';')
+ i = 0
+ for wl in fl:
+ if len(wl) > 2:
+ # replace pow
+ power = re.search('pow\(y, ([0-9]+)\)', wl)
+ if power is not None:
+ np = "y" + "*y" * (int(power.group(1)) - 1)
+ wl = wl.replace(power.group(0), np)
+ sw[i] = '('
+ sw[i] += str(sp.horner(eval(wl.split(';')[0].split('=')[1]) * fac))
+ sw[i] += ')/' + str(fac)
+ i += 1
+ for i, s in enumerate(sw):
+ if not keep:
+ if toOpenCL:
+ res += "inline " + t + " "
+ res += weights_names[i] + "(" + t + " y){\n"
+ res += ' return '
+ else:
+ res += 'lambda y, s: s * '
+ res += '('
+ # replace y**n
+ power = re.findall('y\*\*[0-9]+', s)
+ if power is not None:
+ for pw in power:
+ n = int(pw.split('**')[1])
+ npower = 'y' + "*y" * (n - 1)
+ s = s.replace(pw, npower)
+ s = s.replace(' ', '')
+ if CLBuiltins:
+ s = createFMA(s)
+ # From integers to floats
+ s = re.sub(r"(?P<id>[0-9]+)", r"\g<id>" + cteEnd, s)
+ s = s.replace('*', ' * ')
+ s = s.replace('/', ' / ')
+ s = s.replace('+', ' + ')
+ s = s.replace('-', ' - ')
+ s = s.replace('( - ', '(-')
+ s = s.replace(' ', ' ')
+ s = s.replace(", - ", ", -")
+ res += s + ')'
+ if toOpenCL:
+ res += ";}"
+ if i < len(sw) - 1:
+ res += "\n"
else:
- res += 'lambda y, s: s * '
- res += '('
- # replace y**n
- power = re.findall('y\*\*[0-9]+', s)
- if power is not None:
- for pw in power:
- n = int(pw.split('**')[1])
- npower = 'y' + "*y" * (n - 1)
- s = s.replace(pw, npower)
- s = s.replace(' ', '')
- if CLBuiltins:
- s = createFMA(s)
- # From integers to floats
- s = re.sub(r"(?P<id>[0-9]+)", r"\g<id>" + cteEnd, s)
- s = s.replace('*', ' * ')
- s = s.replace('/', ' / ')
- s = s.replace('+', ' + ')
- s = s.replace('-', ' - ')
- s = s.replace('( - ', '(-')
- s = s.replace(' ', ' ')
- s = s.replace(", - ", ", -")
- res += s + ')'
- if toOpenCL:
- res += ";}"
- if i < len(sw) - 1:
- res += "\n"
- else:
- res += "w[{0}] = ".format(i)
- # replace y**n
- power = re.findall('y\*\*[0-9]+', s)
- if power is not None:
- for pw in power:
- n = int(pw.split('**')[1])
- npower = 'y' + "*y" * (n - 1)
- s = s.replace(pw, npower)
- # From integers to floats
- s = re.sub(r"(?P<id>[0-9]+)", r"\g<id>.", s)
- s = s.replace('*', ' * ')
- s = s.replace('/', ' / ')
- s = s.replace('+', ' + ')
- s = s.replace('-', ' - ')
- s = s.replace('( - ', '(-')
- s = s.replace(' ', ' ')
- s = s.replace(", - ", ", -")
- res += s + "\n"
- return res
+ res += "w[{0}] = ".format(i)
+ # replace y**n
+ power = re.findall('y\*\*[0-9]+', s)
+ if power is not None:
+ for pw in power:
+ n = int(pw.split('**')[1])
+ npower = 'y' + "*y" * (n - 1)
+ s = s.replace(pw, npower)
+ # From integers to floats
+ s = re.sub(r"(?P<id>[0-9]+)", r"\g<id>.", s)
+ s = s.replace('*', ' * ')
+ s = s.replace('/', ' / ')
+ s = s.replace('+', ' + ')
+ s = s.replace('-', ' - ')
+ s = s.replace('( - ', '(-')
+ s = s.replace(' ', ' ')
+ s = s.replace(", - ", ", -")
+ res += s + "\n"
+ return res
+
+except:
+ msge = 'Sympy not available - remeshing formula parsing will not work.'
+ msge += 'If you need parsing, try "pip install sympy" and reinstall hysop.'
+ print msge
def createFMA(s):
diff --git a/HySoP/hysop/tools/tests/test_formula_parsing.py b/HySoP/hysop/tools/tests/test_formula_parsing.py
index acd9fdad565c3f7fdd084e414f1e3251c598cabb..1f093ec1e9c68deeadb23fdb6e62f63b7eb8817f 100644
--- a/HySoP/hysop/tools/tests/test_formula_parsing.py
+++ b/HySoP/hysop/tools/tests/test_formula_parsing.py
@@ -1,5 +1,3 @@
-from hysop.tools.remeshing_formula_parsing import parse
-
m4p = """
w[0] = (-1 + (2 - y) * y) * y / 2;
w[1] = 1 + (-5 + 3 * y) * y * y / 2;
@@ -81,15 +79,19 @@ l63 += """* y) * y) * y / 720;
w[7] = (-145 + (367 + (-311 + 89 * y) * y) * y) * pow(y, 4) / 720;
"""
+try:
+ from hysop.tools.remeshing_formula_parsing import parse
-def test_parsing_toPython():
- assert m4p_python == parse(m4p, 2, toOpenCL=False)
+ def test_parsing_toPython():
+ assert m4p_python == parse(m4p, 2, toOpenCL=False)
+ def test_parsing_toOpenCL():
+ assert m4p_cl_novec == parse(m4p, 2, vec=False, toOpenCL=True)
+ assert m4p_cl == parse(m4p, 2, vec=True, toOpenCL=True)
+ assert m4p_cl_builtin == \
+ parse(m4p, 2, vec=True, toOpenCL=True, CLBuiltins=True)
+ res = parse(l63, 720, vec=True, toOpenCL=True, CLBuiltins=True)
+ assert l6_cl_builtin == res
-def test_parsing_toOpenCL():
- assert m4p_cl_novec == parse(m4p, 2, vec=False, toOpenCL=True)
- assert m4p_cl == parse(m4p, 2, vec=True, toOpenCL=True)
- assert m4p_cl_builtin == \
- parse(m4p, 2, vec=True, toOpenCL=True, CLBuiltins=True)
- assert l6_cl_builtin == parse(l63, 720,
- vec=True, toOpenCL=True, CLBuiltins=True)
+except:
+ pass