add velocity penalization. Add OpenCl implementation of Dummy operator

hysop/backend/host/python/operator/penalization.py

@@ -39,11 +39,11 @@ class PythonPenalizeVorticity(HostOperator):
 
     @debug
     def __new__(cls, obstacles, variables,
-                 velocity, vorticity,
-                 dt, coeff=None, ubar=None, formulation=None, **kwds):
+                velocity, vorticity,
+                dt, coeff=None, ubar=None, formulation=None, **kwds):
         return super(PythonPenalizeVorticity, cls).__new__(cls,
-            input_fields=None, output_fields=None,
-            input_params=None, **kwds)
+                                                           input_fields=None, output_fields=None,
+                                                           input_params=None, **kwds)
 
     @debug
     def __init__(self, obstacles, variables,
@@ -252,3 +252,105 @@ class PythonPenalizeVorticity(HostOperator):
         # Y direction
         self.tmp = self.stencil[1](a=self.tmp_v[0], out=self.tmp, axis=0)
         self.W[0][...] += -self.tmp
+
+
+class PythonPenalizeVelocity(HostOperator, CommonPenalization):
+    """
+
+    """
+    __default_method = {
+        SpaceDiscretization: SpaceDiscretization.FDC4
+    }
+
+    __available_methods = {
+        SpaceDiscretization: (InstanceOf(SpaceDiscretization), InstanceOf(int))
+    }
+
+    @classmethod
+    def default_method(cls):
+        dm = super(PythonPenalizeVelocity, cls).default_method()
+        dm.update(cls.__default_method)
+        return dm
+
+    @classmethod
+    def available_methods(cls):
+        am = super(PythonPenalizeVelocity, cls).available_methods()
+        am.update(cls.__available_methods)
+        return am
+
+    @debug
+    def __new__(cls, obstacles, variables, velocity,
+                dt, coeff=None, ubar=None, formulation=None, **kwds):
+        return super(PythonPenalizeVelocity, cls).__new__(
+            cls, input_fields=None, output_fields=None, input_params=None, **kwds)
+
+    @debug
+    def __init__(self, obstacles, variables,
+                 velocity,
+                 dt, coeff=None, ubar=None, formulation=None, **kwds):
+        check_instance(velocity, Field)
+        check_instance(variables, dict, keys=Field,
+                       values=CartesianTopologyDescriptors)
+
+        input_fields = {velocity: variables[velocity], }
+        output_fields = {velocity: variables[velocity]}
+        input_params = {dt.name: dt}
+        for o in (obstacles.values() if isinstance(obstacles, dict) else obstacles):
+            input_fields[o] = variables[o]
+
+        self.velocity = velocity
+        self.dt = dt
+        super(PythonPenalizeVelocity, self).__init__(
+            input_fields=input_fields, output_fields=output_fields,
+            input_params=input_params,
+            obstacles=obstacles, velocity=velocity, dt=dt, coeff=coeff, ubar=ubar, formulation=formulation, **kwds)
+
+    @debug
+    def discretize(self):
+        if self.discretized:
+            return
+        super(PythonPenalizeVelocity, self).discretize()
+        self.dvelocity = self.input_discrete_tensor_fields[self.velocity]
+
+        dv = self.dvelocity
+        self.dobstacles = {}
+        for c, o in self.obstacles.items():
+            o_df = self.input_discrete_fields[o]
+            self.dobstacles[c] = o_df.data[0]
+
+        self.v_ghost_exchanger = dv.build_ghost_exchanger()
+        if self.v_ghost_exchanger is None:
+            def _dummy_func():
+                pass
+            self.v_ghost_exchanger = _dummy_func
+
+    @debug
+    def get_work_properties(self):
+        requests = super(PythonPenalizeVelocity, self).get_work_properties()
+        buffers = MemoryRequest.empty_like(
+            a=self.dvelocity.buffers[0], nb_components=1, backend=self.dvelocity.backend)
+        requests.push_mem_request('Wtmp', buffers)
+        return requests
+
+    def setup(self, work):
+        super(PythonPenalizeVelocity, self).setup(work=work)
+        Wtmp = work.get_buffer(self, 'Wtmp', handle=True)
+        self.tmp = Wtmp[-1]
+
+    @classmethod
+    def supported_dimensions(cls):
+        return (3, 2)
+
+    @classmethod
+    def supports_mpi(cls):
+        return True
+
+    @op_apply
+    def apply(self, **kwds):
+        super(PythonPenalizeVelocity, self).apply(**kwds)
+
+        # Penalize velocity
+        self._compute_penalization()
+        for ubar, v in zip(self._ubar()[::-1], self.dvelocity.buffers):
+            v[...] += (v - ubar) * self.tmp
+        self.v_ghost_exchanger()

hysop/operator/dummy.py

@@ -16,8 +16,8 @@ class PythonDummy(HostOperator):
 
     @debug
     def __new__(cls, variables, **kwds):
-        return super(PythonDummy, cls).__new__(cls,
-            input_fields=None, output_fields=None, **kwds)
+        return super(PythonDummy, cls).__new__(
+            cls, input_fields=None, output_fields=None, **kwds)
 
     @debug
     def __init__(self, variables, **kwds):
@@ -49,10 +49,49 @@ class PythonDummy(HostOperator):
         return True
 
 
+class OpenClDummy(OpenClOperator):
+
+    @debug
+    def __new__(cls, variables, **kwds):
+        return super(OpenClDummy, cls).__new__(
+            cls,
+            input_fields=None, output_fields=None, **kwds)
+
+    @debug
+    def __init__(self, variables, **kwds):
+        check_instance(variables, dict, keys=Field,
+                       values=CartesianTopologyDescriptors)
+        input_fields = variables
+        output_fields = variables
+        super(OpenClDummy, self).__init__(
+            input_fields=input_fields, output_fields=output_fields,
+            **kwds)
+
+    @debug
+    def get_field_requirements(self):
+        requirements = super(OpenClDummy, self).get_field_requirements()
+        for (is_input, reqs) in requirements.iter_requirements():
+            if (reqs is None):
+                continue
+            (field, td, req) = reqs
+            req.axes = ((0, 1, 2), )
+        return requirements
+
+    @op_apply
+    def apply(self, **kwds):
+        super(OpenClDummy, self).apply(**kwds)
+        # Here doing nothing
+
+    @classmethod
+    def supports_mpi(cls):
+        return True
+
+
 class Dummy(ComputationalGraphNodeFrontend):
 
     __implementations = {
         Implementation.PYTHON: PythonDummy
+        Implementation.OPENCL: OpenClDummy
     }
 
     @classmethod
@@ -62,4 +101,3 @@ class Dummy(ComputationalGraphNodeFrontend):
     @classmethod
     def default_implementation(cls):
         return Implementation.PYTHON
-