diff --git a/Solver/TESTCASES/CylinderEddies.lua b/Solver/TESTCASES/CylinderEddies.lua
index dc80c39967e3bb61df8391f03fe48e7a726578aa..42b5d087cec8962fada6e86500c6f1b1fe2c99e1 100644
--- a/Solver/TESTCASES/CylinderEddies.lua
+++ b/Solver/TESTCASES/CylinderEddies.lua
@@ -32,8 +32,8 @@ kappa:set(0,0,0.01)
 
 print'*** Loading the mesh and the model ***'
 myModel   = GModel  ()
-myModel:load ('cyl.geo')	
-myModel:load ('cyl.msh')
+myModel:load ('cyl2.geo')	
+myModel:load ('cyl2.msh')
 
 print'*** Create a dg solver ***'
 DG = dgSystemOfEquations (myModel)
@@ -59,7 +59,6 @@ DG:exportSolution('output/cyl_0')
 print'*** solve ***'
 
 CFL = 20.1;
-CFL = 5;
 dt = CFL * DG:computeInvSpectralRadius();
 print('DT = ',dt)
 T = 0;
diff --git a/Solver/TESTCASES/cyl.geo b/Solver/TESTCASES/cyl.geo
index 586107725df39284cc2a735dca4748b45ef54f38..5779ec056692987e17dcc56ab56dfa994f7268d8 100644
--- a/Solver/TESTCASES/cyl.geo
+++ b/Solver/TESTCASES/cyl.geo
@@ -62,7 +62,7 @@ Transfinite Surface(1) = {2,12,14,4};
 Transfinite Surface(2) = {5,15,12,2};
 Transfinite Surface(3) = {3,13,15,5};
 Transfinite Surface(4) = {4,14,13,3};
-//Recombine Surface {1,2,3,4};
+Recombine Surface {1,2,3,4};
 
 Line Loop(9) = {6, 7, 8, 5};
 Line Loop(10) = {11,12,13,14};
diff --git a/Solver/TESTCASES/cyl2.geo b/Solver/TESTCASES/cyl2.geo
new file mode 100644
index 0000000000000000000000000000000000000000..08735cab6f3ab49ad96dd77e6057661e5c4c835b
--- /dev/null
+++ b/Solver/TESTCASES/cyl2.geo
@@ -0,0 +1,76 @@
+a = .1;
+b = 1;
+radius = 0.5;
+radiusBoundaryLayer = 0.7;
+Point(1) = {0, 0, 0, a};
+
+nlayers=3;
+np=5;
+np2=2;
+
+
+For n In {0:nlayers}
+  r= radius + (radiusBoundaryLayer-radius)*n/nlayers;
+  Point(n*10+2) = {r,0,0,a};
+  Point(n*10+3) = {0,-r,0,a};
+  Point(n*10+4) = {-r,0,0,a};
+  Point(n*10+5) = {0,r,0,a};
+
+  Circle(n*10+2) = {n*10+2,1,n*10+3};
+  Circle(n*10+3) = {n*10+3,1,n*10+4};
+  Circle(n*10+4) = {n*10+4,1,n*10+5};
+  Circle(n*10+5) = {n*10+5,1,n*10+2};
+
+  Transfinite Line (n*10+2) = np;
+  Transfinite Line (n*10+3) = np;
+  Transfinite Line (n*10+4) = np;
+  Transfinite Line (n*10+5) = np;
+
+  If (n>0)
+    Line (n*10+7) = {(n-1)*10+2, n*10+2};
+    Line (n*10+8) = {(n-1)*10+3, n*10+3};
+    Line (n*10+9) = {(n-1)*10+4, n*10+4};
+    Line (n*10+10) = {(n-1)*10+5, n*10+5};
+
+    Transfinite Line (n*10+7) = np2;
+    Transfinite Line (n*10+8) = np2;
+    Transfinite Line (n*10+9) = np2;
+    Transfinite Line (n*10+10) = np2;
+
+    Line Loop (n*10+2) = {n*10+2, -n*10-8, -(n-1)*10-2, n*10+7};
+    Plane Surface(n*10+2) = {n*10+2};
+    Transfinite Surface (n*10+2) = {n*10+2, n*10+3, (n-1)*10+3, (n-1)*10+2};
+    Line Loop (n*10+3) = {n*10+3, -n*10-9, -(n-1)*10-3, n*10+8};
+    Plane Surface(n*10+3) = {n*10+3};
+    Transfinite Surface (n*10+3) = {n*10+3, n*10+4, (n-1)*10+4, (n-1)*10+3};
+    Line Loop (n*10+4) = {n*10+4, -n*10-10, -(n-1)*10-4, n*10+9};
+    Plane Surface(n*10+4) = {n*10+4};
+    Transfinite Surface (n*10+4) = {n*10+4, n*10+5, (n-1)*10+5, (n-1)*10+4};
+    Line Loop (n*10+5) = {n*10+5, -n*10-7, -(n-1)*10-5, n*10+10};
+    Plane Surface(n*10+5) = {n*10+5};
+    Transfinite Surface (n*10+5) = {n*10+5, n*10+2, (n-1)*10+2, (n-1)*10+5};
+
+    Recombine Surface {n*10+2, n*10+3, n*10+4,n*10+5};
+  EndIf
+EndFor
+
+Point(6) = {-5, -5, 0, b};
+Point(7) = {-5, 5, 0, b};
+Point(8) = {15, 5, 0, b};
+Point(9) = {15, -5, 0, b};
+
+Line(10001) = {7, 8};
+Line(10002) = {8, 9};
+Line(10003) = {9, 6};
+Line(10004) = {6, 7};
+
+Line Loop(1)={10001,10002,10003,10004};
+Line Loop(2)={-nlayers*10-5,-nlayers*10-4,-nlayers*10-3,-nlayers*10-2};
+Plane Surface(1)={1,2};
+Mesh.CharacteristicLengthExtendFromBoundary=1;
+
+Physical Line("Box") = {10004, 10001, 10002, 10003};
+Physical Line("Cylinder") = {2,3,4,5};
+
+Physical Surface("BoundaryLayer") = {2:10000};
+Physical Surface("Air") = {1};