diff --git a/examples/bubble/periodic_bubble_levelset.py b/examples/bubble/periodic_bubble_levelset.py
index d4b0559858cb5be5198326848f0b9c0da6789239..c8a5afdbf457a6174b679ca42bdb35739b87b604 100644
--- a/examples/bubble/periodic_bubble_levelset.py
+++ b/examples/bubble/periodic_bubble_levelset.py
@@ -182,7 +182,7 @@ def compute(args):
     
     #> External force rot(-rho*g)
     Fext = np.zeros(shape=(dim,), dtype=object).view(SymbolicTensor)
-    Fext[dim-1] = -1.0 #-9.8196
+    Fext[1] = -1.0 #-9.8196
     Fext *= rhos
     lhs = Ws.diff(frame.time)
     rhs = curl(Fext, frame)
@@ -236,7 +236,7 @@ def compute(args):
     
 
     ### Adaptive timestep operator
-    adapt_dt = AdaptiveTimeStep(dt, equivalent_CFL=True, max_dt=0.01,
+    adapt_dt = AdaptiveTimeStep(dt, equivalent_CFL=True, max_dt=1e-2,
                                     name='merge_dt', pretty_name='dt', )
     dt_cfl   = adapt_dt.push_cfl_criteria(cfl=args.cfl, Finf=min_max_U.Finf,
                                             equivalent_CFL=True, 
@@ -401,6 +401,6 @@ if __name__=='__main__':
                         dump_freq=0, 
                         dump_times=(0.0, 0.1, 0.20, 0.30, 0.325, 0.4, 0.45, 0.50),
                         rho1=1.0, rho2=10.0, mu1=0.00025, mu2=0.00050,
-                        Bc = ((0.5,0.2,0.5),(0.5,0.5,0.5),), Br = (0.1,0.15,))
+                        Bc = ((0.5,0.1,0.5),(0.5,0.4,0.5),), Br = (0.1,0.15,))
 
     parser.run(compute)
diff --git a/hysop/operator/directional/stretching_diffusion_dir.py b/hysop/operator/directional/stretching_diffusion_dir.py
index 9ccac9b2a64f6a05e6334fda488719c7255cb321..598027c76c33fd59e74920c8e3d759dae35f3ead 100644
--- a/hysop/operator/directional/stretching_diffusion_dir.py
+++ b/hysop/operator/directional/stretching_diffusion_dir.py
@@ -49,13 +49,15 @@ class DirectionalStretchingDiffusion(DirectionalStretching):
             Base class keywords arguments.
             See hysop.operator.directional.stretching_dir.DirectionalStretching.
         """
-        check_instance(viscosity, float)
+        check_instance(viscosity, (float, Field))
         self.viscosity = viscosity
         super(DirectionalStretchingDiffusion, self).__init__(**kwds)
         
     def _gen_expressions(self, formulation, velocity, vorticity, C, A):
         from hysop.symbolic.field import laplacian
         viscosity = self.viscosity
+        if isinstance(viscosity, Field):
+            viscosity = viscosity.s(*viscosity.domain.frame.vars)
         _exprs = super(DirectionalStretchingDiffusion, self)._gen_expressions(formulation=formulation,
                 velocity=velocity, vorticity=vorticity, C=C, A=A)
         exprs = ()
diff --git a/hysop/operator/directional/symbolic_dir.py b/hysop/operator/directional/symbolic_dir.py
index ac9a659550744e77269bd87b01f470149ddfe3a4..03a6c39111e96fa0570aac12db23fe3802244c53 100644
--- a/hysop/operator/directional/symbolic_dir.py
+++ b/hysop/operator/directional/symbolic_dir.py
@@ -111,6 +111,7 @@ class DirectionalSymbolic(DirectionalOperatorFrontend):
                         directional_exprs[k] += (v,)
                     else:
                         directional_exprs[k]  = (v,)
+                print
         self._directional_exprs = directional_exprs
         self._no_split = no_split
 
diff --git a/hysop/symbolic/directional.py b/hysop/symbolic/directional.py
index 9f11c279cca3a6ef43dd89ec114b852f1a64e6f1..00ee06ff2b32df5a63a60cdff4c64c9b60c4f0da 100644
--- a/hysop/symbolic/directional.py
+++ b/hysop/symbolic/directional.py
@@ -40,8 +40,9 @@ def split(F, coords=None):
         res[i] = collect_direction(F, xi)
     residue = (F - sum(res.values())).simplify()
     count = sum((r != 0) for r in res.values())
-    for i in res.keys():
-        res[i] += (residue/count)
+    if (count>0):
+        for i in res.keys():
+            res[i] += (residue/count)
     return res
 
 def split_assignement(expr, coords=None):