diff --git a/Examples/testPenalization.py b/Examples/testPenalization.py
index 192becbff17e797cce928260a0b70c7570ec9608..a060db6aa48175529449cc2a7dcff1b67c454fe4 100644
--- a/Examples/testPenalization.py
+++ b/Examples/testPenalization.py
@@ -12,11 +12,11 @@ from parmepy.problem.simulation import Simulation
 
 
 def vitesse(x, y, z):
-    return x, y, z
+    return x
 
 
 def vitesse2(x, y):
-    return x, y
+    return x
 
 nb = 33
 
@@ -31,19 +31,21 @@ scal = Field(domain=dom, name='Scalar')
 op = Analytic(scal, formula=vitesse, resolutions={scal: resol3D})
 scal2 = Field(domain=dom2, name='Scalar')
 op2 = Analytic(scal2, formula=vitesse2, resolutions={scal2: resol2D})
+velo = Field(domain=dom2, name='Velo', isVector=True)
 
 op.setUp()
 op2.setUp()
 topo = dom.topologies[0]
 coords = topo.mesh.coords
 
-ff = scal.discreteFields[topo].data
+ff = scal.discreteFields[topo].data[0]
+
 
 ff[...] = 128
 
 topo2 = dom2.topologies.values()[0]
 
-ff2 = scal2.discreteFields[topo2].data
+ff2 = scal2.discreteFields[topo2].data[0]
 ff2[...] = 128
 sphere = Sphere(dom, position=[0., 0., 0.], radius=0.5, porousLayers=[0.13])
 
@@ -59,13 +61,20 @@ hcyl = SemiCylinder2D(dom2, position=[0., 0.], radius=0.5, porousLayers=[0.13])
 ind = hsphere.discretize(topo)
 
 penal = Penalization(scal, [hsphere, plates], coeff=[1e6, 10],
-                     resolutions={scal: resol3D})
-penal2 = Penalization(scal2, [hcyl, plates2], coeff=[1e6, 10],
+                     resolutions={velo: resol3D})
+penal2 = Penalization(velo, [hcyl, plates2], coeff=[1e6, 10],
                       resolutions={scal2: resol2D})
 
 penal.setUp()
 penal2.setUp()
-printer = Printer(fields=[scal2], frequency=1)
+
+velod = velo.discreteFields[topo2]
+velod[0] = 128
+velod[1] = 12
+#velod[2] = 4.3
+
+
+printer = Printer(fields=[velo], frequency=1)
 printer.setUp()
 print scal.norm()
 simulation = Simulation()
@@ -74,3 +83,5 @@ penal2.apply(simulation)
 printer.apply(simulation)
 print "print ..."
 print scal.norm()
+
+velo.dump('sc')
diff --git a/Examples/testPoisson.py b/Examples/testPoisson.py
index 2f3fb5cc01be137f132f9ee4da882995b8a424ff..728769663df9bc751e6f20386ec79b7133f83b5a 100755
--- a/Examples/testPoisson.py
+++ b/Examples/testPoisson.py
@@ -2,6 +2,10 @@
 import parmepy as pp
 from parmepy.operator.poisson import Poisson
 from parmepy.operator.analytic import Analytic
+from parmepy.problem.problem import Problem
+from parmepy.problem.navier_stokes import NSProblem
+from parmepy.problem.simulation import Simulation
+from parmepy.operator.monitors.printer import Printer
 
 import numpy as np
 import math
@@ -37,7 +41,7 @@ vorticity = pp.Field(domain=dom, formula=computeVort,
 
 
 # ref. field
-def computeRef(x, y, z, t, dt, ite):
+def computeRef(x, y, z, t):
     refx = -2. * pi / LL[1] * \
         (cos(x * cc[0]) * sin(y * cc[1]) * sin(z * cc[2])) \
         - 2. * pi / LL[2] * (cos(x * cc[0]) * sin(y * cc[1]) * cos(z * cc[2]))
@@ -62,26 +66,39 @@ ref = pp.Field(domain=dom, name='reference', isVector=True)
 refOp = Analytic(ref, formula=computeRef,
                  resolutions={ref: resol})
 
-refOp.setUp()
-refOp.apply(0, 0, 0)
+timeStep = 0.25
+simu = Simulation()
+#simu.iterMax = simu.iterMax - 2
 
+refOp.setUp()
 poisson.setUp()
+refOp.apply(simu)
 
-vorticity.initialize()
+#pb = NSProblem([refOp, poisson], simu,
+#               monitors=[Printer(fields=[velocity],
+#                                 frequency=1)],
+#               dumpFreq=2, name='toto')
 
-print vorticity.norm()
+#pb.setUp()
 
-poisson.apply()
+#pb.solve()
 
-print velocity.norm()
-print ref.norm()
+## refOp.apply(0, 0, 0)
 
-assert np.allclose(ref.norm(), velocity.norm())
+vorticity.initialize()
+## print vorticity.norm()
+poisson.apply(simu)
+## print velocity.norm()
+## print ref.norm()
 
-poisson.finalize()
 
+refD = ref.discreteFields.values()[0]
+vd = velocity.discreteFields.values()[0]
 
-poisson.printComputeTime()
+assert np.allclose(ref.norm(), velocity.norm())
+for i in range(dom.dimension):
+    assert np.allclose(vd[i], refD[i])
 
+poisson.finalize()
 
-print poisson
+#print poisson
diff --git a/HySoP/hysop/fields/continuous.py b/HySoP/hysop/fields/continuous.py
index 1888638bda9b64ce92fc448487e2e59ed6ecf8d4..c1ffa113a36dac3d3fb331b518192f210fd5ffef 100644
--- a/HySoP/hysop/fields/continuous.py
+++ b/HySoP/hysop/fields/continuous.py
@@ -38,7 +38,8 @@ class Field(object):
         return object.__new__(cls, *args, **kw)
 
     @debug
-    def __init__(self, domain, formula=None, name=None, isVector=False):
+    def __init__(self, domain, formula=None, name=None, isVector=False,
+                 nbComponents=None):
         """
         Create a scalar of vector field which dimension is determined
         by its domain of definition.
@@ -69,6 +70,11 @@ class Field(object):
         self._formula = formula
         ## A list of extra parameters, used in _formula
         self.extraParameters = tuple([])
+        ## Number of components of the field
+        if nbComponents is None:
+            self.nbComponents = self.dimension if self.isVector else 1
+        else:
+            self.nbComponents = nbComponents
 
     @debug
     def discretize(self, topo):
@@ -85,9 +91,12 @@ class Field(object):
             return self.discreteFields[topo]
         else:
             nameD = self.name + '_' + str(topo.getId())
-            self.discreteFields[topo] = DiscreteField(topo,
-                                                      isVector=self.isVector,
-                                                      name=nameD)
+            self.discreteFields[topo] = DiscreteField(
+                topo,
+                isVector=self.isVector,
+                name=nameD,
+                nbComponents=self.nbComponents)
+
         return self.discreteFields[topo]
 
     def setFormula(self, formula):
diff --git a/HySoP/hysop/fields/discrete.py b/HySoP/hysop/fields/discrete.py
index 5d0d338f3bad9c70eb1daa1dca45df5636b0731c..63c76c5dcc1f7d82726b0173b8813d5b7a1bc96e 100644
--- a/HySoP/hysop/fields/discrete.py
+++ b/HySoP/hysop/fields/discrete.py
@@ -69,7 +69,7 @@ class DiscreteField(object):
         return object.__new__(cls, *args, **kw)
 
     @debug
-    def __init__(self, topology, isVector=False, name=None):
+    def __init__(self, topology, isVector=False, name=None, nbComponents=None):
         """ Discretize a field on a given topology.
 
         @param topology : a topology (mpi topo + mesh, parmepy.mpi.topology)
@@ -98,7 +98,10 @@ class DiscreteField(object):
         ## Object to store computational times
         self.timer = Timer(self)
         ## Number of components of the field
-        self.nbComponents = self.dimension if self.isVector else 1
+        if nbComponents is None:
+            self.nbComponents = self.dimension if self.isVector else 1
+        else:
+            self.nbComponents = nbComponents
         ## The memory space for data ...
         self.data = [np.zeros((self.resolution),
                               dtype=PARMES_REAL, order=ORDER)
diff --git a/HySoP/hysop/mpi/topology.py b/HySoP/hysop/mpi/topology.py
index 4e6b33a872bcc8c99efa358a70bca535be55a202..1f586eb54037d8cf7c21e8af72551ed258733c29 100644
--- a/HySoP/hysop/mpi/topology.py
+++ b/HySoP/hysop/mpi/topology.py
@@ -488,12 +488,9 @@ class Bridge(object):
         ## --> no mpi messages
         self.ito = []
 
-        ## Flag for useless bridges when topologies are identical
-        self.isUseless = False
+        self.hasLocalInter = True
+
         if self.topoFrom == self.topoTo:
-            # Nothing to do ...
-            print "Warning, Building a useles bridge."
-            self.isUseless = True
             return
 
         # Both topologies must have the
@@ -550,7 +547,6 @@ class Bridge(object):
 
         localFrom = np.zeros((2 * dom.dimension), dtype=np.int32)
 
-        self.hasLocalInter = True
         for dim in range(dom.dimension):
             indexTo = range(iglob2[main_rank, dim * 2],
                             iglob2[main_rank, dim * 2 + 1] + 1)
diff --git a/HySoP/hysop/numerics/differential_operations.py b/HySoP/hysop/numerics/differential_operations.py
index d9df8082685cbaf84480433bd340f0fdd714a5c3..a1ad50a4ee0287b977057db5f5b6fec578be960a 100755
--- a/HySoP/hysop/numerics/differential_operations.py
+++ b/HySoP/hysop/numerics/differential_operations.py
@@ -21,7 +21,7 @@ class DifferentialOperation(object):
     @abstractmethod
     def __init__(self, topo, method=None):
         pass
-    
+
 
 class Curl(DifferentialOperation):
     """
diff --git a/HySoP/hysop/operator/discrete/curl.py b/HySoP/hysop/operator/discrete/curl.py
index 84160f232faca9557d8db228132a24f7a82a10ab..78e2e991b709c130dba6289faf1ef342544446c3 100644
--- a/HySoP/hysop/operator/discrete/curl.py
+++ b/HySoP/hysop/operator/discrete/curl.py
@@ -1,3 +1,4 @@
+# -*- coding: utf-8 -*-
 """
 @file discrete/curl.py
 
diff --git a/HySoP/hysop/operator/discrete/scales_advection.py b/HySoP/hysop/operator/discrete/scales_advection.py
index 9a552589510200b8ee0d0a76189f1b641b721684..038e30618f7256f11d3688a7b8a08d085e4a8b76 100644
--- a/HySoP/hysop/operator/discrete/scales_advection.py
+++ b/HySoP/hysop/operator/discrete/scales_advection.py
@@ -47,91 +47,13 @@ class ScalesAdvection(DiscreteOperator):
 
         dt = simulation.timeStep
 
-#        for adF in self.advectedFields:
-#            if isinstance(adF, VectorField):
-#                for d in xrange(adF.dimension):
-#                    adF[d][...] = scales2py.solve_advection(
-#                        dt, self.velocity.data[0],
-#                        self.velocity.data[1],
-#                        self.velocity.data[2],
-#                        adF[d][...])
-#            else:
-#                adF[...] = scales2py.solve_advection(
-#                    dt, self.velocity.data[0],
-#                    self.velocity.data[1],
-#                    self.velocity.data[2],
-#                    adF.data)
-
-
-        ## Space discretization of grad(u)
-        mesh_size = self.velocity.topology.mesh.space_step
-        ind=np.arange(self.velocity.topology.mesh.resolution[0])
-
-###      First order scheme
-###       X components of temp and result
-        temp1 = (self.velocity[0][...] -
-                 self.velocity[0][np.roll(ind,+1),...]) / (mesh_size[0]) 
-
-        temp2 = (self.velocity[0][...] -
-                 self.velocity[0][:,np.roll(ind,+1),:]) / (mesh_size[1]) 
-
-        temp3 = (self.velocity[0][...] -
-                 self.velocity[0][...,np.roll(ind,+1)]) / (mesh_size[2]) 
-
-        maxstr1= np.max(abs(temp1)+abs(temp2)+abs(temp3))
-        maxadv1= np.max(abs(temp1))
-
-#       Y components of temp and result
-        temp4 = (self.velocity[1][...] -
-                 self.velocity[1][np.roll(ind,+1),...]) / (mesh_size[0]) 
-
-        temp5 = (self.velocity[1][...] -
-                 self.velocity[1][:,np.roll(ind,+1),:]) / (mesh_size[1]) 
-
-        temp6 = (self.velocity[1][...] -
-                 self.velocity[1][...,np.roll(ind,+1)]) / (mesh_size[2]) 
-
-        maxstr2= np.max(abs(temp4)+abs(temp5)+abs(temp6))
-        maxadv2= np.max(abs(temp5))
-
-#       Z components of temp and result
-        temp7 = (self.velocity[2][...] -
-                 self.velocity[2][np.roll(ind,+1),...]) / (mesh_size[0]) 
-
-        temp8 = (self.velocity[2][...] -
-                 self.velocity[2][:,np.roll(ind,+1),:]) / (mesh_size[1]) 
-
-        temp9 = (self.velocity[2][...] -
-                 self.velocity[2][...,np.roll(ind,+1)]) / (mesh_size[2]) 
-
-        maxstr3= np.max(abs(temp7)+abs(temp8)+abs(temp9))
-        maxadv3= np.max(abs(temp9))
-        maxAdv = max(maxadv1,maxadv2,maxadv3)
-
-        stabcoeff = 1.0 / 2.0
-
-        print "dtAdvec", stabcoeff / maxAdv
-
-        ## Time Subcycles
-        dtadv = dt
-        dtadv = min(dtadv, stabcoeff / maxAdv)
-        if dt == dtadv:
-            print "dtadv, ndt, stab/max ", dtadv, int(dt / dtadv),\
-                  stabcoeff / maxAdv
-            ndt = int(dt / dtadv)
-        else:
-            print "dtadv, ndt", dtadv, int(dt / dtadv) + 1, stabcoeff / maxAdv
-            ndt = int(dt / dtadv) + 1
-        dtadv = dt / float(ndt)
-
-        for i in range(ndt):
-            for adF in self.advectedFields:
-                for d in xrange(adF.nbComponents):
-                    adF[d][...] = scales2py.solve_advection(
-                        dtadv, self.velocity.data[0],
-                        self.velocity.data[1],
-                        self.velocity.data[2],
-                        adF[d][...])
+        for adF in self.advectedFields:
+            for d in xrange(adF.nbComponents):
+                adF[d][...] = scales2py.solve_advection(
+                    dt, self.velocity.data[0],
+                    self.velocity.data[1],
+                    self.velocity.data[2],
+                    adF[d][...])
 
     def finalize(self):
         """
diff --git a/HySoP/hysop/operator/energy_enstrophy.py b/HySoP/hysop/operator/energy_enstrophy.py
index b118150495693b07426e916b036ff7c27f8996c9..ec1a6238a9d43f66d17e8936c0918f157fd73967 100644
--- a/HySoP/hysop/operator/energy_enstrophy.py
+++ b/HySoP/hysop/operator/energy_enstrophy.py
@@ -23,55 +23,58 @@ class Energy_enstrophy(Monitoring):
     """
 
     def __init__(self, velocity, vorticity,
-                 resolutions=None, viscosity=None,
-                 frequency=None, prefix=None):
+                 topo, viscosity, frequency, prefix=None):
         """
         Constructor.
         @param fields : List of fields
         @param resolutions :  list of resolutions (one per variable)
         @param viscosity :  kinematic viscosity
+        @param topo : the topology on which we want to monitor the fields
         @param frequency : output file producing frequency.
         @param prefix : file name prefix, contains relative path.
         """
-        Monitoring.__init__(self, [velocity, vorticity], frequency,
+        Monitoring.__init__(self, [velocity, vorticity], frequency, topo=topo,
                             name="Energy Enstrophy")
         if prefix is None:
             self.prefix = './Energy'
         else:
             self.prefix = prefix
-
+        ## velocity field
         self.velocity = velocity
+        ## vorticity field
         self.vorticity = vorticity
-        self.resolutions = resolutions
+        ## viscosity (scalar)
         self.viscosity = viscosity
-        self.buffer_1 = 0.
-        self.buffer_2 = 0.
+        self._buffer_1 = 0.
+        self._buffer_2 = 0.
         if (main_rank == 0):
             self.f = open(self.prefix, 'w')
         self.input = [velocity, vorticity]
         self.output = []
+        self.topo = self._predefinedTopo
+        self.topovel = self.topo
+        self.topovort = self.topo
+        # \todo : rewrite for multiresolution case
 
     def setUp(self):
         Monitoring.setUp(self)
-        discreteFields = {}
-        self.spaceStep = {}
-        topodims = np.ones((self.domain.dimension))
-        topodims[-1] = main_size
-
-        # Variables discretization
-        for v in self.variables:
-            # the topology for v ...
-            topo = self.domain.getOrCreateTopology(self.domain.dimension,
-                                                   self.resolutions[v],
-                                                   topodims)
-            # ... and the corresponding discrete field
-            discreteFields[v] = v.discretize(topo)
-            self.spaceStep[v] = topo.mesh.space_step
-        self.velo = discreteFields[self.velocity]
-        self.vort = discreteFields[self.vorticity]
+
+        self.discreteFields[self.velocity] = \
+            self.velocity.discretize(self.topovel)
+        self.discreteFields[self.vorticity] =\
+            self.vorticity.discretize(self.topovort)
 
         self._isUpToDate = True
 
+        spaceStep = [self.topovel.mesh.spaceStep,
+                     self.topovort.mesh.spaceStep]
+
+        length = self.topo.domain.length
+        self.coeffEnergy = 0.5 * (np.prod(spaceStep[0]) /
+                                  np.prod(length))
+        self.coeffEnstrophy = (np.prod(spaceStep[1]) /
+                               np.prod(length))
+
     @debug
     @timed_function
     def apply(self, simulation=None):
@@ -79,38 +82,25 @@ class Energy_enstrophy(Monitoring):
         Computation of kinetic energy, enstrophy &
         Checking energy and enstrophy decay
         """
+        if simulation is None:
+            raise ValueError("Missing simulation value for computation.")
+
         t = simulation.time
         dt = simulation.timeStep
         ite = simulation.currentIteration
         # Kinetic energy computation
-        dvol = np.prod(self.spaceStep[self.velocity])
-
-        if (self.velo.dimension == 2):
-            squareNormVel = np.sum(self.velo.data[0] ** 2 +
-                                   self.velo.data[1] ** 2) * dvol
-        elif (self.velo.dimension == 3):
-            squareNormVel = np.sum(self.velo.data[0] ** 2 +
-                                   self.velo.data[1] ** 2 +
-                                   self.velo.data[2] ** 2) * dvol
-        else:
-            raise ValueError("invalid problem dimension")
-
-        energy = (0.5 * squareNormVel) / (2.0 * np.pi) ** 3
+        vd = self.discreteFields[self.velocity]
+        # get the list of computation points (no ghosts)
+        ind = self.topo.mesh.iCompute
+        nbc = vd.nbComponents
+        energy = self.coeffEnergy * np.sum([vd[i][ind] ** 2
+                                            for i in xrange(nbc)])
 
         # Enstrophy computation
-        dvol = np.prod(self.spaceStep[self.vorticity])
-
-        if (self.vort.dimension == 2):
-            enstrophy = np.sum(self.vort.data ** 2) * dvol
-        elif (self.vort.dimension == 3):
-            enstrophy = np.sum(self.vort.data[0] ** 2 +
-                               self.vort.data[1] ** 2 +
-                               self.vort.data[2] ** 2) * dvol
-        else:
-            raise ValueError("invalid problem dimension")
-
-        enstrophy = enstrophy / (2.0 * np.pi) ** 3
-
+        vortd = self.discreteFields[self.vorticity]
+        nbc = vortd.nbComponents
+        enstrophy = self.coeffEnstrophy * np.sum([vortd[i][ind] ** 2
+                                                  for i in xrange(nbc)])
         # Collective communications
 
         topovelo = self.velo.topology
diff --git a/HySoP/hysop/operator/monitors/monitoring.py b/HySoP/hysop/operator/monitors/monitoring.py
index 5c4c0344363537fa29dc2f9fcb069284140399fe..7118f88899da07a61c05cd231714ddf2875c733a 100644
--- a/HySoP/hysop/operator/monitors/monitoring.py
+++ b/HySoP/hysop/operator/monitors/monitoring.py
@@ -13,13 +13,14 @@ class Monitoring(Operator):
     __metaclass__ = ABCMeta
 
     @abstractmethod
-    def __init__(self, variables, frequency, name=None):
+    def __init__(self, variables, frequency, topo=None,
+                 ghosts=None, name=None):
         """ Constructor
         @param variables : list of fields to monitor.
         @param frequency : output rate.
         @param name : Monitor name.
         """
-        Operator.__init__(self, variables, name=name)
+        Operator.__init__(self, variables, topo=topo, ghosts=ghosts, name=name)
         ## Monitor frequency
         self.freq = frequency
         ## Object to store computational times of lower level functions
diff --git a/HySoP/hysop/operator/redistribute.py b/HySoP/hysop/operator/redistribute.py
index da6a258e6906deb939da305b4b89813072a68c61..aa385c297f009ec426864f2f6664eec766d2da1a 100644
--- a/HySoP/hysop/operator/redistribute.py
+++ b/HySoP/hysop/operator/redistribute.py
@@ -45,9 +45,6 @@ class Redistribute(Operator):
         @param opFrom : source operator
         @param opTo : target (i.e.) the operator that handles the topology on
         which data must be redistributed.
-        @param method : in {'FromDiscreteOperators'} :
-          - 'FromDiscreteOperators' : create the operator from discrete
-          operators that are using variable's discretisations.
         """
 
         Operator.__init__(self, variables, method, name="Redistribute")
@@ -57,28 +54,11 @@ class Redistribute(Operator):
         ## Targeted operator
         self.opTo = opTo
 
-        ## Variables discretisations used un opFrom
-        self.vd_from = {}
-        ## Variables discretisations used un opTo
-        self.vd_to = {}
         # Then check if variables belong to both operators
         for v in self.variables:
-            if method == 'FromDiscreteOperators':
-                # find the discretisation used by operators
-                for vd in v.discreteFields.values():
-                    if vd in opFrom.variables:
-                        self.vd_from[v] = vd
-                    if vd in opTo.variables:
-                        self.vd_to[v] = vd
-                assert v in self.vd_from.keys() and v in self.vd_to.keys(), \
-                    'Redistribute error : one of the variable is not present\
- in both source and target operator.'
-            else:
-                assert v in opFrom.variables and v in opTo.variables, \
-                    'Redistribute error : one of the variable is not present\
- in both source and target operator.'
-                self.vd_from[v] = self.opFrom.discreteFields[v]
-                self.vd_to[v] = self.opTo.discreteFields[v]
+            assert v in opFrom.variables and v in opTo.variables, \
+                'Redistribute error : one of the variable is not present\
+                in both source and target operator.'
         self.input = self.output = self.variables
 
     @debug
@@ -93,8 +73,9 @@ class Redistribute(Operator):
 
         self.bridges = {}
         for v in self.variables:
-            self.bridges[v] = Bridge(self.vd_from[v].topology,
-                                     self.vd_to[v].topology)
+            vd_from = self.opFrom.discreteFields[v]
+            vd_to = self.opTo.discreteFields[v]
+            self.bridges[v] = Bridge(vd_from.topology, vd_to.topology)
 
         self._isUpToDate = True
 
@@ -112,75 +93,47 @@ class Redistribute(Operator):
             br = self.bridges[v]
             self.r_request = {}
             self.s_request = {}
-            if not br.isUseless:
-                if v.isVector:
-                    # Apply for each component of the data
-                    for d in range(dim):
-                        vTo = self.vd_to[v].data[d]
-                        vFrom = self.vd_from[v].data[d]
-                        if br.hasLocalInter:
-                            vTo[br.ito] = vFrom[br.ifrom]
-
-                        for rk in br.recvFrom.keys():
-                            ind = rk * (d + 1)
-                            subvshape = tuple([br.recvFrom[rk][i].stop -
-                                               br.recvFrom[rk][i].start
-                                               for i in range(dim)])
-                            substart = tuple([br.recvFrom[rk][i].start
-                                              for i in range(dim)])
-                            recvTYPE = PARMES_MPI_REAL.Create_subarray(
-                                vTo.shape, subvshape, substart,
-                                order=ORDERMPI)
-                            recvTYPE.Commit()
-                            self.r_request[ind] = main_comm.Irecv(
-                                [vTo, 1, recvTYPE], source=rk, tag=rk)
-                        for rk in br.sendTo.keys():
-                            ind = rk * (d + 1)
-                            subvshape = tuple([br.sendTo[rk][i].stop -
-                                               br.sendTo[rk][i].start
-                                               for i in range(dim)])
-                            substart = tuple([br.sendTo[rk][i].start
-                                              for i in range(dim)])
-                            sendTYPE = PARMES_MPI_REAL.Create_subarray(
-                                vFrom.shape, subvshape, substart,
-                                order=ORDERMPI)
-                            sendTYPE.Commit()
-                            self.s_request[ind] = main_comm.Issend(
-                                [vFrom, 1, sendTYPE], dest=rk, tag=main_rank)
-
-                # scalar
-                else:
-                    vTo = self.vd_to[v].data
-                    vFrom = self.vd_from[v].data
+            # Apply for each component of the data
+            for d in range(v.nbComponents):
+                vTo = self.opTo.discreteFields[v].data[d]
+                vFrom = self.opFrom.discreteFields[v].data[d]
+                if br.hasLocalInter:
                     vTo[br.ito] = vFrom[br.ifrom]
-                    for rk in br.recvFrom.keys():
-                        ind = rk * (d + 1)
-                        subvshape = tuple([br.recvFrom[rk][i].stop -
-                                           br.recvFrom[rk][i].start
-                                           for i in range(dim)])
-                        substart = tuple([br.recvFrom[rk][i].start
-                                          for i in range(dim)])
-                        recvTYPE = PARMES_MPI_REAL.Create_subarray(
-                            vTo.shape,  subvshape, substart,
-                            order=ORDERMPI)
-                        recvTYPE.Commit()
-
-                        self.r_request[ind] = main_comm.Irecv(
-                            [vTo, 1, recvTYPE], source=rk, tag=rk)
-                    for rk in br.sendTo.keys():
-                        ind = rk * (d + 1)
-                        subvshape = tuple([br.sendTo[rk][i].stop -
-                                           br.sendTo[rk][i].start
-                                           for i in range(dim)])
-                        substart = tuple([br.sendTo[rk][i].start
-                                          for i in range(dim)])
-                        sendTYPE = PARMES_MPI_REAL.Create_subarray(
-                            vFrom.shape, subvshape, substart,
-                            order=ORDERMPI)
-                        sendTYPE.Commit()
-
-                        self.s_request[ind] = main_comm.Issend(
-                            [vFrom, 1, sendTYPE], dest=rk, tag=main_rank)
+
+                for rk in br.recvFrom.keys():
+                    ind = rk * (d + 1)
+                    subvshape = tuple([br.recvFrom[rk][i].stop -
+                                       br.recvFrom[rk][i].start
+                                       for i in range(dim)])
+                    substart = tuple([br.recvFrom[rk][i].start
+                                      for i in range(dim)])
+                    recvTYPE = PARMES_MPI_REAL.Create_subarray(vTo.shape,
+                                                               subvshape,
+                                                               substart,
+                                                               order=
+                                                               ORDERMPI)
+                    recvTYPE.Commit()
+                    self.r_request[ind] = main_comm.Irecv([vTo, 1,
+                                                           recvTYPE],
+                                                          source=rk,
+                                                          tag=rk)
+                for rk in br.sendTo.keys():
+                    ind = rk * (d + 1)
+                    subvshape = tuple([br.sendTo[rk][i].stop -
+                                       br.sendTo[rk][i].start
+                                       for i in range(dim)])
+                    substart = tuple([br.sendTo[rk][i].start
+                                      for i in range(dim)])
+                    sendTYPE = PARMES_MPI_REAL.Create_subarray(vFrom.shape,
+                                                               subvshape,
+                                                               substart,
+                                                               order=
+                                                               ORDERMPI)
+                    sendTYPE.Commit()
+                    self.s_request[ind] = main_comm.Issend([vFrom, 1,
+                                                            sendTYPE],
+                                                           dest=rk,
+                                                           tag=main_rank)
 
     def wait(self):
         """