diff --git a/hysop/codegen/functions/stretching_rhs.py b/hysop/codegen/functions/stretching_rhs.py
index b7d8503491344dee69845bff290e2700ceb9571f..d3482594dc07c000d23d71b1812f9842730285e6 100644
--- a/hysop/codegen/functions/stretching_rhs.py
+++ b/hysop/codegen/functions/stretching_rhs.py
@@ -102,8 +102,9 @@ class DirectionalStretchingRhsFunction(OpenClFunctionCodeGenerator):
                 
         args[W]  = CodegenVectorClBuiltin(W, ftype, dim, typegen, add_impl_const=True)
         args['inv_dx'] = CodegenVariable('inv_dx', ftype, typegen, add_impl_const=True, nl=True)
-
-        args['rk_step'] = CodegenVariable('rk_step', itype, typegen, add_impl_const=True)
+        
+        if is_conservative:
+            args['rk_step'] = CodegenVariable('rk_step', itype, typegen, add_impl_const=True)
         if is_periodic and (not cached):
             args['base']   = CodegenVariable('base',  itype, typegen, add_impl_const=True,nl=True)
             args['offset'] = CodegenVariable('offset', itype, typegen, add_impl_const=True,nl=True)
@@ -218,8 +219,10 @@ class DirectionalStretchingRhsFunction(OpenClFunctionCodeGenerator):
         W = s.args[W]
         lidx = s.args['lidx']
         Lx = s.args['Lx']
-        rk_step = s.args['rk_step']
+        if is_conservative:
+            rk_step = s.args['rk_step']
 
+        offset = s.args['offset']
         
         dw_dt  = CodegenVectorClBuiltin('dW_dt', ftype, dim, tg)
         ghosts = CodegenVariable('ghosts', itype, tg, const=True, value=order/2, symbolic_mode=True)
@@ -232,19 +235,23 @@ class DirectionalStretchingRhsFunction(OpenClFunctionCodeGenerator):
             if is_conservative:
                 s.jumpline()
                 s.comment('Synchronize required vorticiy component across workitems')
-                s.mem_fence(read=True,_local=cached,_global=not cached)
+                s.barrier(_local=cached,_global=not cached)
                 with s._block_():
                     code = '{} = {};'.format(Wd[offset()], W[direction])
                     s.append(code)
-                s.mem_fence(write=True,_local=cached,_global=not cached)
+                s.barrier(_local=cached,_global=not cached)
                 s.jumpline(2)
 
             dw_dt.declare(s,init=0)
         
-            cond = '({lid}>={step}*{ghosts}) && ({lid}<{L}-{step}*{ghosts})'.format(
-                    lid=lidx(), L=Lx(), ghosts=ghosts(), step='{step}')
+            if is_conservative:
+                cond = '({lid}>={step}*{ghosts}) && ({lid}<{L}-{step}*{ghosts})'.format(
+                        lid=lidx(), L=Lx(), ghosts=ghosts(), step='({}+1)'.format(rk_step()))
+            else:
+                cond = '({lid}>={ghosts}) && ({lid}<{L}-{ghosts})'.format(
+                        lid=lidx(), L=Lx(), ghosts=ghosts())
             
-            with s._if_(cond.format(step='({}+1)'.format(rk_step()))):
+            with s._if_(cond):
 
                 if is_conservative:
                     s.comment('compute d(U*W{0})/d{0} using finite difference stencil'.format(xyz[direction]))
diff --git a/hysop/codegen/kernels/directional_stretching.py b/hysop/codegen/kernels/directional_stretching.py
index 16b852ca4d6b778b39219195ef1a512b71f4de61..a7949439e28510d680ee9aa72d6ba5b8854394a8 100644
--- a/hysop/codegen/kernels/directional_stretching.py
+++ b/hysop/codegen/kernels/directional_stretching.py
@@ -40,14 +40,14 @@ from hysop.gpu.gpu_discrete import GPUDiscreteField
 class DirectionalStretchingKernel(KernelCodeGenerator):
 
     @staticmethod
-    def codegen_name(ftype,cached,direction,sformulation):
-        cache = 'cached_' if cached else ''
+    def codegen_name(ftype,is_cached,direction,sformulation):
+        cache = 'cached_' if is_cached else ''
         sformulation = sformulation.lower()
         return 'directional_{}stretching_{}_{}{}'.format(cache,sformulation,
                 ftype[0],S_DIR[direction])
    
     def __init__(self, typegen, dim,  ftype, order, direction,
-                       cached, boundary, 
+                       is_cached, boundary, 
                        formulation, rk_scheme,
                        symbolic_mode = False,
                        known_vars = None):
@@ -72,18 +72,21 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
         is_conservative  = (formulation==StretchingFormulation.CONSERVATIVE)
         is_periodic      = (boundary==BoundaryCondition.PERIODIC)
         
-        if cached:
+        if is_cached:
             storage = OpenClCodeGenerator.default_keywords['local']
         else:
             storage = OpenClCodeGenerator.default_keywords['global']
+       
+        if is_conservative and not is_cached:
+            raise ValueError('Conservetive stretching requires caching.')
 
-        name = DirectionalStretchingKernel.codegen_name(ftype,cached,direction,sformulation)
+        name = DirectionalStretchingKernel.codegen_name(ftype,is_cached,direction,sformulation)
 
-        kernel_reqs = self.build_requirements(typegen, dim, ftype, order, cached, 
+        kernel_reqs = self.build_requirements(typegen, dim, ftype, order, is_cached, 
                 rk_scheme, direction, boundary, symbolic_mode, formulation, storage,
                 is_periodic)
 
-        kernel_args = self.gen_kernel_arguments(typegen, dim, ftype, kernel_reqs, cached,
+        kernel_args = self.gen_kernel_arguments(typegen, dim, ftype, kernel_reqs, is_cached,
                 local_size_known)
         
         super(DirectionalStretchingKernel,self).__init__(
@@ -99,7 +102,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
 
         self.order            = order
         self.ftype            = ftype
-        self.cached           = cached
+        self.is_cached           = is_cached
         self.direction        = direction
         self.dim              = dim
         self.boundary         = boundary
@@ -144,7 +147,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
     def required_workgroup_cache_size(self, local_work_size):
         dim             = self.work_dim
         ftype           = self.ftype
-        cached          = self.cached
+        is_cached          = self.is_cached
         direction       = self.direction
         cache_ghosts    = self.cache_ghosts()
         is_periodic     = self.is_periodic
@@ -152,7 +155,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
         flt_bytes       = self.typegen.FLT_BYTES[ftype]
        
         sc,dc = 0,0
-        if cached: 
+        if is_cached: 
             count = dim*local_work_size[0]
             if is_conservative:
                 count += local_work_size[0]
@@ -172,7 +175,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
         
         return (sc,dc)
 
-    def build_requirements(self,typegen,work_dim,ftype,order,cached,rk_scheme,direction,
+    def build_requirements(self,typegen,work_dim,ftype,order,is_cached,rk_scheme,direction,
             boundary,force_symbolic,formulation,storage,is_periodic):
         tg=typegen 
         reqs = WriteOnceDict()
@@ -189,7 +192,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
         
          
         stretching_rhs = DirectionalStretchingRhsFunction(typegen=typegen, dim=work_dim, 
-                ftype=ftype, cached=cached,
+                ftype=ftype, cached=is_cached,
                 order=order, direction=direction, boundary=boundary,
                 formulation=formulation,
                 restrict=True, vectorize_u=False,
@@ -207,7 +210,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
             
         return reqs
     
-    def gen_kernel_arguments(self, typegen, work_dim, ftype, requirements,cached,
+    def gen_kernel_arguments(self, typegen, work_dim, ftype, requirements,is_cached,
             local_size_known):
         
         xyz = 'xyz'
@@ -231,7 +234,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
 
         kargs['mesh_info'] = requirements['MeshInfoStruct'].build_codegen_variable(const=True,name='mesh_info')
 
-        if cached and not local_size_known:
+        if is_cached and not local_size_known:
              kargs['buffer'] = CodegenVariable(storage=_local,ctype=ftype, add_impl_const=True,
                      name='buffer', ptr=True, restrict=True, typegen=typegen, nl=False)
 
@@ -257,7 +260,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
         dim       = s.dim
         ftype     = s.ftype
         boundary  = s.boundary
-        cached    = s.cached
+        is_cached    = s.is_cached
         storage   = s.storage
         symbolic_mode = s.symbolic_mode
 
@@ -306,7 +309,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
         local_work = CodegenVariable('lwork','int',tg,const=True)
         
         cached_vars = ArgDict()
-        if cached:
+        if is_cached:
             for i in xrange(work_dim):
                 Vi = self.svelocity+self.xyz[i]
                 if local_size_known:
@@ -388,7 +391,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
                         init='{} - 2*{}'.format(local_size[0],cache_ghosts()))
             s.jumpline()
             
-            if cached:
+            if is_cached:
                 with s._align_() as al:
                     for varname,var in cached_vars.iteritems():
                         var.declare(al,align=True)
@@ -484,7 +487,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
 
                 s.jumpline()
 
-                if self.cached:
+                if self.is_cached:
                     for i in xrange(work_dim):
                         Ui = self.svelocity+self.xyz[i]
                         Uic = cached_vars[Ui]
@@ -504,7 +507,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
                 rk_args={'dt': dt,
                          'inv_dx': inv_dx,
                          'W': W}
-                if is_periodic and (not cached):
+                if is_periodic and (not is_cached):
                     base = CodegenVariable('base','int',typegen=tg,const=True)
                     base.declare(s,init='({}/{}) * {}'.format(global_index(),grid_size[0],grid_size[0]))
                     offset = CodegenVariable('offset','int',typegen=tg,const=True)
@@ -513,19 +516,19 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
                     rk_args['offset'] = offset
                     rk_args['width']  = grid_size[0]
                 else:
-                    rk_args['offset'] = local_id[0] if cached else global_index
+                    rk_args['offset'] = local_id[0] if is_cached else global_index
                 rk_args['Lx']  = local_size[0]
                 rk_args['lidx'] = local_id[0]
                 for i in xrange(work_dim):
                     Ui_name  = self.svelocity+xyz[i]
-                    if cached:
+                    if is_cached:
                         Ui = cached_vars[Ui_name]
                     else:
                         Ui  = s.vars[Ui_name]
                     rk_args[Ui_name] = Ui
                 if is_conservative:
                     Wd_name  = self.svorticity+xyz[direction]
-                    if cached:
+                    if is_cached:
                         Wd = cached_vars[Wd_name]
                     else:
                         Wd  = s.vars[Wd_name]
@@ -548,7 +551,7 @@ class DirectionalStretchingKernel(KernelCodeGenerator):
                     cond='{}=={} && {}=={}'.format(global_id[2],compute_grid_ghosts[2],global_id[1],compute_grid_ghosts[1])
                     #cond='{0} == 0 || {0} == 49'.format(global_index[0])
                     with s._if_(cond):
-                        if is_periodic and not cached:
+                        if is_periodic and not is_cached:
                             code = 'printf(\"write k=%i, gidx=%i, lidx=%i, offset=%i, GID=%i, base=%i, offset=%i, width=%i\\n\", k,gid.x,lid.x, gid.x-compute_grid_size.x,GID,base,offset,grid_size.x);'
                         else:
                             code = 'printf(\"write k=%i, gidx=%i, lidx=%i, GID=%i, Ux=%f, U.x=%f, Wx=%f, dW=%f \\n\", k,gid.x,lid.x, GID, Ux[GID], U.x, Wx[GID], W.x-Wx[GID]);'
@@ -583,7 +586,7 @@ if __name__ == '__main__':
         order=4, dim=dim, direction=0, 
         formulation=StretchingFormulation.GRAD_UW,
         rk_scheme=ExplicitRungeKutta('RK2'),
-        cached=True,
+        is_cached=True,
         symbolic_mode=True,
         boundary=BoundaryCondition.NONE,
         known_vars=dict(
diff --git a/hysop/codegen/kernels/tests/test_directional_stretching.py b/hysop/codegen/kernels/tests/test_directional_stretching.py
index 6b83beb3f5e8068bae512d7b2621678296f5af18..ebc05b4fd6419ba71c116c9cd2e12a9d6ca6fdcb 100644
--- a/hysop/codegen/kernels/tests/test_directional_stretching.py
+++ b/hysop/codegen/kernels/tests/test_directional_stretching.py
@@ -105,10 +105,9 @@ class TestDirectionalStretching(object):
         cls.host_buffers_gpu       = host_buffers_gpu
         cls.device_buffers         = device_buffers
 
-        cls.local_work_size = np.asarray([16,1,1])
+        cls.local_work_size = np.asarray([27,1,1])
         cls.inv_dx = inv_dx
-        #cls.dt = dtype(0.5)
-        cls.dt = 2./inv_dx[0]
+        cls.dt = dtype(0.5)
 
         cls.enable_error_plots = enable_error_plots
 
@@ -159,15 +158,121 @@ class TestDirectionalStretching(object):
         velocity  = ['ux','uy','uz']
         U = [host_init_buffers[name].copy() for name in velocity]
         W = [host_init_buffers[name].copy() for name in vorticity]
-        dUdx = [deriv(ui) for ui in U]
+        dUdx  = [deriv(ui) for ui in U]
         
         if rk_scheme.name() == 'Euler':
+            Wc = [ w.copy() for w in W ]
             if formulation==StretchingFormulation.GRAD_UW:
                 for i in xrange(3):
-                    if i!=direction:
-                        W[i] += dt*dUdx[i]*W[direction]
-                W[direction] += dt*dUdx[direction]*W[direction]
-
+                    W[i] += dt*dUdx[i]*Wc[direction]
+            elif formulation==StretchingFormulation.GRAD_UW_T:
+                for i in xrange(3):
+                    W[direction] += dt*dUdx[i]*Wc[i]
+            elif formulation==StretchingFormulation.MIXED_GRAD_UW:
+                for i in xrange(3):
+                    W[i]         += 0.5*dt*dUdx[i]*Wc[direction]
+                    W[direction] += 0.5*dt*dUdx[i]*Wc[i] 
+            elif formulation==StretchingFormulation.CONSERVATIVE:
+                W0 = Wc
+                K0 = [deriv(ui*W0[direction]) for ui in U]
+                W1 = [W0[i] + dt*K0[i]  for i in xrange(3)]
+                W  = W1
+            else:
+                msg = 'Unknown stretching formulation scheme {}.'.format(
+                        StretchingFormulation.svalue(formulation))
+                raise ValueError(msg)
+        elif rk_scheme.name() == 'RK2':
+            Wc = [ w.copy() for w in W ]
+            if formulation==StretchingFormulation.GRAD_UW:
+                W0 = Wc
+                K0 = [dUdx[i]*W0[direction] for i in xrange(3)]
+                W1 = [W0[i] + 0.5*dt*K0[i]  for i in xrange(3)]
+                K1 = [dUdx[i]*W1[direction] for i in xrange(3)]
+                W2 = [W0[i] + 1.0*dt*K1[i]  for i in xrange(3)]
+                W  = W2
+            elif formulation==StretchingFormulation.GRAD_UW_T:
+                W0 = Wc
+                K0 = sum([dUdx[i]*W0[i] for i in xrange(3)])
+                W1 = W0[direction] + 0.5*dt*K0
+                K1 = sum([dUdx[i]*W0[i] for i in xrange(3) if i!=direction]) + W1*dUdx[direction]
+                W2 = W0
+                W2[direction] += dt*K1
+                W  = W2
+            elif formulation==StretchingFormulation.MIXED_GRAD_UW:
+                W0 = Wc
+                K0             =     [0.5*dUdx[i]*W0[direction] for i in xrange(3)]
+                K0[direction] += sum([0.5*dUdx[i]*W0[i]         for i in xrange(3)])
+                W1 = [W0[i] + 0.5*dt*K0[i]  for i in xrange(3)]
+                K1             =     [0.5*dUdx[i]*W1[direction] for i in xrange(3)]
+                K1[direction] += sum([0.5*dUdx[i]*W1[i]         for i in xrange(3)])
+                W2 = [W0[i] + 1.0*dt*K1[i]  for i in xrange(3)]
+                W  = W2
+            elif formulation==StretchingFormulation.CONSERVATIVE:
+                W0 = Wc
+                K0 = [deriv(ui*W0[direction]) for ui in U]
+                W1 = [W0[i] + 0.5*dt*K0[i]  for i in xrange(3)]
+                K1 = [deriv(ui*W1[direction]) for ui in U]
+                W2 = [W0[i] + 1.0*dt*K1[i]  for i in xrange(3)]
+                W  = W2
+            else:
+                msg = 'Unknown stretching formulation scheme {}.'.format(
+                        StretchingFormulation.svalue(formulation))
+                raise ValueError(msg)
+        elif rk_scheme.name() == 'RK4':
+            Wc = [ w.copy() for w in W ]
+            if formulation==StretchingFormulation.GRAD_UW:
+                W0 = Wc
+                K0 = [dUdx[i]*W0[direction] for i in xrange(3)]
+                W1 = [W0[i] + 0.5*dt*K0[i]  for i in xrange(3)]
+                K1 = [dUdx[i]*W1[direction] for i in xrange(3)]
+                W2 = [W0[i] + 0.5*dt*K1[i]  for i in xrange(3)]
+                K2 = [dUdx[i]*W2[direction] for i in xrange(3)]
+                W3 = [W0[i] + 1.0*dt*K2[i]  for i in xrange(3)]
+                K3 = [dUdx[i]*W3[direction] for i in xrange(3)]
+                K =  [1./6*K0[i] + 1./3*K1[i] + 1./3*K2[i] + 1./6*K3[i] for i in xrange(3)]
+                W4 = [W0[i] + dt*K[i]  for i in xrange(3)]
+                W  = W4
+            elif formulation==StretchingFormulation.GRAD_UW_T:
+                W0 = Wc
+                K0 = sum([dUdx[i]*W0[i] for i in xrange(3)])
+                W1 = W0[direction] + 0.5*dt*K0
+                K1 = sum([dUdx[i]*W0[i] for i in xrange(3) if i!=direction]) + W1*dUdx[direction]
+                W2 = W0[direction] + 0.5*dt*K1
+                K2 = sum([dUdx[i]*W0[i] for i in xrange(3) if i!=direction]) + W2*dUdx[direction]
+                W3 = W0[direction] + 1.0*dt*K2
+                K3 = sum([dUdx[i]*W0[i] for i in xrange(3) if i!=direction]) + W3*dUdx[direction]
+                K =  1./6*K0 + 1./3*K1 + 1./3*K2 + 1./6*K3
+                W4 = W0
+                W4[direction] += dt*K
+                W  = W4
+            elif formulation==StretchingFormulation.MIXED_GRAD_UW:
+                W0 = Wc
+                K0             =     [0.5*dUdx[i]*W0[direction] for i in xrange(3)]
+                K0[direction] += sum([0.5*dUdx[i]*W0[i]         for i in xrange(3)])
+                W1 = [W0[i] + 0.5*dt*K0[i]  for i in xrange(3)]
+                K1             =     [0.5*dUdx[i]*W1[direction] for i in xrange(3)]
+                K1[direction] += sum([0.5*dUdx[i]*W1[i]         for i in xrange(3)])
+                W2 = [W0[i] + 0.5*dt*K1[i]  for i in xrange(3)]
+                K2             =     [0.5*dUdx[i]*W2[direction] for i in xrange(3)]
+                K2[direction] += sum([0.5*dUdx[i]*W2[i]         for i in xrange(3)])
+                W3 = [W0[i] + 1.0*dt*K2[i]  for i in xrange(3)]
+                K3             =     [0.5*dUdx[i]*W3[direction] for i in xrange(3)]
+                K3[direction] += sum([0.5*dUdx[i]*W3[i]         for i in xrange(3)])
+                K =  [1./6*K0[i] + 1./3*K1[i] + 1./3*K2[i] + 1./6*K3[i] for i in xrange(3)]
+                W4 = [W0[i] + dt*K[i]  for i in xrange(3)]
+                W  = W4
+            elif formulation==StretchingFormulation.CONSERVATIVE:
+                W0 = Wc
+                K0 = [deriv(ui*W0[direction]) for ui in U]
+                W1 = [W0[i] + 0.5*dt*K0[i]  for i in xrange(3)]
+                K1 = [deriv(ui*W1[direction]) for ui in U]
+                W2 = [W0[i] + 0.5*dt*K1[i]  for i in xrange(3)]
+                K2 = [deriv(ui*W2[direction]) for ui in U]
+                W3 = [W0[i] + 1.0*dt*K2[i]  for i in xrange(3)]
+                K3 = [deriv(ui*W3[direction]) for ui in U]
+                K =  [1./6*K0[i] + 1./3*K1[i] + 1./3*K2[i] + 1./6*K3[i] for i in xrange(3)]
+                W4 = [W0[i] + dt*K[i]  for i in xrange(3)]
+                W  = W4
             else:
                 msg = 'Unknown stretching formulation scheme {}.'.format(
                         StretchingFormulation.svalue(formulation))
@@ -234,7 +339,7 @@ class TestDirectionalStretching(object):
             direction=direction, 
             formulation=formulation,
             rk_scheme=rk_scheme,
-            cached=cached,
+            is_cached=cached,
             symbolic_mode=False,
             boundary=boundary,
             known_vars=known_vars)
@@ -253,7 +358,7 @@ class TestDirectionalStretching(object):
     
         print '\tGenerating and compiling Kernel...'
         source = dsk.__str__()
-        prg    = cl.Program(self.typegen.context, source)
+        prg = cl.Program(self.typegen.context, source)
         prg.build(devices=[self.typegen.device])
         kernel = prg.all_kernels()[0]
         kernel.set_args(*kernel_args)
@@ -290,7 +395,7 @@ class TestDirectionalStretching(object):
         for (name,host,dev) in buffers:
             (l1,l2,linf) = self._distances(host,dev,view)
             print '\t{} -> l1={}  l2={}  linf={}'.format(name,l1,l2,linf)
-            if linf>1e-12:
+            if linf>1e-10:
                 err_buffers.append(name)
                 good = False
         if not good:
@@ -335,59 +440,89 @@ class TestDirectionalStretching(object):
         assert formulation in StretchingFormulation.entries()
         assert isinstance(rk_scheme, ExplicitRungeKutta)
 
-        direction=0
         cached=[False,True]
-        orders=[2,4,6]
         boundaries=[BoundaryCondition.NONE,BoundaryCondition.PERIODIC]
+        directions=[0,1,2]
+        orders=[2,4,6]
         
         for cache in cached:
+            if (formulation==StretchingFormulation.CONSERVATIVE) and not cache: 
+                continue
             for boundary in boundaries:
-                for order in orders:
-                    self._do_compute_cpu(order=order, direction=direction, boundary=boundary,
-                            formulation=formulation, rk_scheme=rk_scheme)
-                    
-                    self._do_compute_gpu_and_check(order=order, direction=direction, boundary=boundary,
-                            formulation=formulation, rk_scheme=rk_scheme, cached=cache)
+                for direction in directions:
+                    for order in orders:
+                        self._do_compute_cpu(order=order, direction=direction, boundary=boundary,
+                                formulation=formulation, rk_scheme=rk_scheme)
+                        
+                        self._do_compute_gpu_and_check(order=order, direction=direction, boundary=boundary,
+                                formulation=formulation, rk_scheme=rk_scheme, cached=cache)
    
 
-    #########
-    # tests #   
-    #########
-    def _test_gencode(self):
-        dim = 3
-        directions = range(dim)
-        orders = [4,6]
-        cached = [False,True]
-        symbolic_mode = True
-        formulations = StretchingFormulation.fields()
-        rk_schemes = ['Euler', 'RK4']
-        boundaries = [BoundaryCondition.NONE, BoundaryCondition.PERIODIC]
-        
-        typegen=self.typegen
 
-        for symbolic in symbolic_mode:
-            for boundary in boundaries:
-                for formulation in formulations:
-                    for rk_scheme in rk_schemes:
-                        for direction in directions:
-                            for cache in cached:
-                                for order in orders:
-                                    dsk = DirectionalStretchingKernel(
-                                        typegen=typegen, 
-                                        ftype=typegen.fbtype,
-                                        order=order, dim=dim, direction=direction, 
-                                        formulation=StretchingFormulation.GRAD_UW,
-                                        rk_scheme=ExplicitRungeKutta(rk_scheme),
-                                        cached=cache,
-                                        symbolic_mode=symbolic,
-                                        boundary=boundary
-                                    )
-                                    dsk.test_compile()
 
     def test_gradUW_Euler(self):
         formulation=StretchingFormulation.GRAD_UW
         rk_scheme=ExplicitRungeKutta('Euler')
         self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_gradUW_T_Euler(self):
+        formulation=StretchingFormulation.GRAD_UW_T
+        rk_scheme=ExplicitRungeKutta('Euler')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_mixed_gradUW_Euler(self):
+        formulation=StretchingFormulation.MIXED_GRAD_UW
+        rk_scheme=ExplicitRungeKutta('Euler')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_conservative_Euler(self):
+        formulation=StretchingFormulation.CONSERVATIVE
+        rk_scheme=ExplicitRungeKutta('Euler')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+   
+
+
+    def test_gradUW_RK2(self):
+        formulation=StretchingFormulation.GRAD_UW
+        rk_scheme=ExplicitRungeKutta('RK2')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_gradUW_T_RK2(self):
+        formulation=StretchingFormulation.GRAD_UW_T
+        rk_scheme=ExplicitRungeKutta('RK2')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_mixed_gradUW_RK2(self):
+        formulation=StretchingFormulation.MIXED_GRAD_UW
+        rk_scheme=ExplicitRungeKutta('RK2')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_conservative_RK2(self):
+        formulation=StretchingFormulation.CONSERVATIVE
+        rk_scheme=ExplicitRungeKutta('RK2')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+   
+
+    def test_gradUW_RK4(self):
+        formulation=StretchingFormulation.GRAD_UW
+        rk_scheme=ExplicitRungeKutta('RK4')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_gradUW_T_RK4(self):
+        formulation=StretchingFormulation.GRAD_UW_T
+        rk_scheme=ExplicitRungeKutta('RK4')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_mixed_gradUW_RK4(self):
+        formulation=StretchingFormulation.MIXED_GRAD_UW
+        rk_scheme=ExplicitRungeKutta('RK4')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
+    
+    def test_conservative_RK4(self):
+        formulation=StretchingFormulation.CONSERVATIVE
+        rk_scheme=ExplicitRungeKutta('RK4')
+        self._check_kernels(formulation=formulation, rk_scheme=rk_scheme)
 
 
 
@@ -396,7 +531,20 @@ if __name__ == '__main__':
     TestDirectionalStretching.setup_class(enable_error_plots=True)
     test = TestDirectionalStretching()
     
-    test.test_gradUW_Euler()
+    #test.test_gradUW_Euler()
+    #test.test_gradUW_T_Euler()
+    #test.test_mixed_gradUW_Euler()
+    #test.test_conservative_Euler()
+    
+    #test.test_gradUW_RK2()
+    #test.test_gradUW_T_RK2()
+    #test.test_mixed_gradUW_RK2()
+    #test.test_conservative_RK2()
+    
+    test.test_gradUW_RK4()
+    test.test_gradUW_T_RK4()
+    test.test_mixed_gradUW_RK4()
+    test.test_conservative_RK4()
 
     TestDirectionalStretching.teardown_class()