diff --git a/examples/inviscid-flow/CMakeLists.txt b/examples/inviscid-flow/CMakeLists.txt
index 634cce16141410c692254c2b7385b2cac5acd5a7..19f8d03003c876c4c510a7d5ed5c616a3e2e0888 100644
--- a/examples/inviscid-flow/CMakeLists.txt
+++ b/examples/inviscid-flow/CMakeLists.txt
@@ -6,19 +6,19 @@ set( tnl_inviscid_flow_SOURCES
      euler.cu )
                
 IF( BUILD_CUDA )
-   CUDA_ADD_EXECUTABLE(tnl-euler-2d${debugExt} euler.cu)
-   target_link_libraries (tnl-euler-2d${debugExt} tnl${debugExt}-${tnlVersion}  ${CUSPARSE_LIBRARY} )
+   CUDA_ADD_EXECUTABLE(tnl-euler${debugExt} euler.cu)
+   target_link_libraries (tnl-euler${debugExt} tnl${debugExt}-${tnlVersion}  ${CUSPARSE_LIBRARY} )
 ELSE(  BUILD_CUDA )               
-   ADD_EXECUTABLE(tnl-euler-2d${debugExt} euler.cpp)     
-   target_link_libraries (tnl-euler-2d${debugExt} tnl${debugExt}-${tnlVersion} )
+   ADD_EXECUTABLE(tnl-euler${debugExt} euler.cpp)     
+   target_link_libraries (tnl-euler${debugExt} tnl${debugExt}-${tnlVersion} )
 ENDIF( BUILD_CUDA )
 
 
-INSTALL( TARGETS tnl-euler-2d${debugExt}
+INSTALL( TARGETS tnl-euler${debugExt}
          RUNTIME DESTINATION bin
          PERMISSIONS OWNER_READ OWNER_WRITE OWNER_EXECUTE GROUP_READ GROUP_EXECUTE WORLD_READ WORLD_EXECUTE )
         
 INSTALL( FILES run-euler
                ${tnl_inviscid_flow_SOURCES}
-         DESTINATION share/tnl-${tnlVersion}/examples/inviscid-flow-2d )
+         DESTINATION share/tnl-${tnlVersion}/examples/inviscid-flow )
 
diff --git a/examples/narrow-band/main.h b/examples/narrow-band/main.h
index 7f30a9da9083f651de35bf09974456785137f33f..e9f55656f1128eb268214eb85a4bc7fe26ec5773 100644
--- a/examples/narrow-band/main.h
+++ b/examples/narrow-band/main.h
@@ -47,7 +47,7 @@ int main( int argc, char* argv[] )
 
    if(dim == 2)
    {
-		tnlNarrowBand<tnlGrid<2,double,tnlHost, int>, double, int> solver;
+		tnlNarrowBand<tnlGrid<2,double,TNL::Devices::Host, int>, double, int> solver;
 		if(!solver.init(parameters))
 	   {
 			cerr << "Solver failed to initialize." << endl;
@@ -60,7 +60,7 @@ int main( int argc, char* argv[] )
    }
 //   else if(dim == 3)
 //   {
-//		tnlNarrowBand<tnlGrid<3,double,tnlHost, int>, double, int> solver;
+//		tnlNarrowBand<tnlGrid<3,double,TNL::Devices::Host, int>, double, int> solver;
 //		if(!solver.init(parameters))
 //	   {
 //			cerr << "Solver failed to initialize." << endl;
diff --git a/examples/narrow-band/tnlNarrowBand.h b/examples/narrow-band/tnlNarrowBand.h
index e2817424a1feb0f24bcbf8a759d37e8c9697977d..7d3d19bc03b43247f735cf04c5c82368360a748d 100644
--- a/examples/narrow-band/tnlNarrowBand.h
+++ b/examples/narrow-band/tnlNarrowBand.h
@@ -17,10 +17,10 @@
 #define TNLNARROWBAND_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
 #include <functions/tnlMeshFunction.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 #include <limits.h>
@@ -55,7 +55,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -135,7 +135,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 	tnlNarrowBand();
diff --git a/examples/narrow-band/tnlNarrowBand2D_CUDA_v4_impl.h b/examples/narrow-band/tnlNarrowBand2D_CUDA_v4_impl.h
index 3aab9ede171d869b3a20c89a0d20c61c1543b82d..310d0fb239260fcc740a0f2bc56b37a36662a2f6 100644
--- a/examples/narrow-band/tnlNarrowBand2D_CUDA_v4_impl.h
+++ b/examples/narrow-band/tnlNarrowBand2D_CUDA_v4_impl.h
@@ -277,7 +277,7 @@ void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: u
 	int subgridID = i/NARROWBAND_SUBGRID_SIZE + (j/NARROWBAND_SUBGRID_SIZE) * statusGridSize;
 	if(cudaStatusVector[subgridID] != 0 && i<Mesh.getDimensions().x() && j < Mesh.getDimensions().y())
 	{
-		tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+		tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 		Entity.setCoordinates(CoordinatesType(i,j));
 		Entity.refresh();
 		tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -317,7 +317,7 @@ void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: u
 }
 
 
-__global__ void initCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 
 
@@ -346,7 +346,7 @@ bool tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: i
 	int i = threadIdx.x + blockDim.x*blockIdx.x;
 	int j = blockDim.y*blockIdx.y + threadIdx.y;
 
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
@@ -497,7 +497,7 @@ Real tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: f
 
 
 
-__global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 
@@ -603,7 +603,7 @@ __global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double
 
 
 
-__global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	__shared__ double u0;
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
@@ -641,7 +641,7 @@ __global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int
 
 
 // run this with one thread per block
-__global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 //	printf("Hello\n");
 	if(solver->cudaStatusVector[blockIdx.x + gridDim.x*blockIdx.y] == 1)
@@ -670,7 +670,7 @@ __global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, in
 
 
 
-__global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, double tau)
+__global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, double tau)
 {
 	int gid = (blockDim.y*blockIdx.y + threadIdx.y)*solver->Mesh.getDimensions().x()+ threadIdx.x;
 	int i = threadIdx.x + blockIdx.x*blockDim.x;
@@ -688,10 +688,10 @@ __global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int
 
 		if(status != 0)
 		{
-			tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(solver->Mesh);
+			tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(solver->Mesh);
 			Entity.setCoordinates(Containers::StaticVector<2,double>(i,j));
 			Entity.refresh();
-			tnlNeighborGridEntityGetter<tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
+			tnlNeighborGridEntityGetter<tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
 			double value = solver->cudaDofVector2[Entity.getIndex()];
 			double xf,xb,yf,yb, grad, fu, a,b;
 			a = b = 0.0;
@@ -837,7 +837,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -856,7 +856,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -875,7 +875,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -908,7 +908,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -941,7 +941,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -974,7 +974,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1008,7 +1008,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1041,7 +1041,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1074,7 +1074,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1107,7 +1107,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1144,7 +1144,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1177,7 +1177,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1210,7 +1210,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1234,7 +1234,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1267,7 +1267,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1300,7 +1300,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
diff --git a/examples/narrow-band/tnlNarrowBand2D_CUDA_v5_impl.h b/examples/narrow-band/tnlNarrowBand2D_CUDA_v5_impl.h
index 92d5e1e982e8d6cefe62ea088c3039d122254fef..cb41d572674bd468cfa0286e274e630dd56fdfa6 100644
--- a/examples/narrow-band/tnlNarrowBand2D_CUDA_v5_impl.h
+++ b/examples/narrow-band/tnlNarrowBand2D_CUDA_v5_impl.h
@@ -276,7 +276,7 @@ void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: u
 	int subgridID = i/NARROWBAND_SUBGRID_SIZE + (j/NARROWBAND_SUBGRID_SIZE) * ((Mesh.getDimensions().x() + NARROWBAND_SUBGRID_SIZE-1 ) / NARROWBAND_SUBGRID_SIZE);
 	if(/*cudaStatusVector[subgridID] != 0 &&*/ i<Mesh.getDimensions().x() && Mesh.getDimensions().y())
 	{
-		tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+		tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 		Entity.setCoordinates(CoordinatesType(i,j));
 		Entity.refresh();
 		tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -316,7 +316,7 @@ void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: u
 }
 
 
-__global__ void initCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 
 
@@ -345,7 +345,7 @@ bool tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: i
 	int i = threadIdx.x + blockDim.x*blockIdx.x;
 	int j = blockDim.y*blockIdx.y + threadIdx.y;
 
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
@@ -495,7 +495,7 @@ Real tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: f
 
 
 
-__global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 
@@ -601,7 +601,7 @@ __global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double
 
 
 
-__global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	__shared__ double u0;
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
@@ -639,7 +639,7 @@ __global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int
 
 
 // run this with one thread per block
-__global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 //	printf("Hello\n");
 	if(solver->cudaStatusVector[blockIdx.x + gridDim.x*blockIdx.y] == 1)
@@ -666,7 +666,7 @@ __global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, in
 
 
 
-__global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, double tau)
+__global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, double tau)
 {
 	int gid = (blockDim.y*blockIdx.y + threadIdx.y)*solver->Mesh.getDimensions().x()+ threadIdx.x;
 	int i = threadIdx.x + blockIdx.x*blockDim.x;
@@ -684,10 +684,10 @@ __global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int
 
 //		if(status != 0)
 		{
-			tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(solver->Mesh);
+			tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(solver->Mesh);
 			Entity.setCoordinates(Containers::StaticVector<2,double>(i,j));
 			Entity.refresh();
-			tnlNeighborGridEntityGetter<tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
+			tnlNeighborGridEntityGetter<tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
 			double value = solver->cudaDofVector2[Entity.getIndex()];
 			double xf,xb,yf,yb, grad, fu, a,b;
 			a = b = 0.0;
@@ -833,7 +833,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -852,7 +852,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -871,7 +871,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -904,7 +904,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -937,7 +937,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -970,7 +970,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1004,7 +1004,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1037,7 +1037,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1070,7 +1070,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1103,7 +1103,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1140,7 +1140,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1173,7 +1173,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1206,7 +1206,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1230,7 +1230,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1263,7 +1263,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1296,7 +1296,7 @@ template< typename MeshReal,
           typename Index >
 void tnlNarrowBand< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
diff --git a/examples/narrow-band/tnlNarrowBand3D_CUDA_impl.h b/examples/narrow-band/tnlNarrowBand3D_CUDA_impl.h
index 4a972db1d06a7f66561e591df51970e5114be93a..2b96ad58205d0578776b2b2d645eb71292c0fad4 100644
--- a/examples/narrow-band/tnlNarrowBand3D_CUDA_impl.h
+++ b/examples/narrow-band/tnlNarrowBand3D_CUDA_impl.h
@@ -165,7 +165,7 @@ template< typename MeshReal,
 __device__
 void tnlNarrowBand< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > :: updateValue( Index i, Index j, Index k)
 {
-	tnlGridEntity< tnlGrid< 3,double, tnlHost, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 3,double, TNL::Devices::Host, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j,k));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 3, tnlGridEntityNoStencilStorage >,3> neighborEntities(Entity);
@@ -222,7 +222,7 @@ template< typename MeshReal,
 __device__
 bool tnlNarrowBand< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > :: initGrid(int i, int j, int k)
 {
-	tnlGridEntity< tnlGrid< 3,double, tnlHost, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 3,double, TNL::Devices::Host, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j,k));
 	Entity.refresh();
 	int gid = Entity.getIndex();
@@ -255,7 +255,7 @@ Real tnlNarrowBand< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > :: f
 
 
 
-__global__ void runCUDA(tnlNarrowBand< tnlGrid< 3,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlNarrowBand< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 	int gx = 0;
@@ -474,7 +474,7 @@ __global__ void runCUDA(tnlNarrowBand< tnlGrid< 3,double, tnlHost, int >, double
 }
 
 
-__global__ void initCUDA(tnlNarrowBand< tnlGrid< 3,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlNarrowBand< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
 	int gy = blockDim.y*blockIdx.y + threadIdx.y;
diff --git a/examples/narrow-band/tnlNarrowBand_CUDA.h b/examples/narrow-band/tnlNarrowBand_CUDA.h
index 8da92f5fc570275c85ebba69113f0e59a4be7988..ca9b1da2cc6e26b14bc003532b6eea75e89d907d 100644
--- a/examples/narrow-band/tnlNarrowBand_CUDA.h
+++ b/examples/narrow-band/tnlNarrowBand_CUDA.h
@@ -17,9 +17,9 @@
 #define TNLNARROWBAND_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 
@@ -54,7 +54,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 	tnlNarrowBand();
@@ -138,7 +138,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -183,16 +183,16 @@ protected:
 
 #ifdef HAVE_CUDA
 //template<int sweep_t>
-__global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i);
-//__global__ void runCUDA(tnlNarrowBand< tnlGrid< 3,double, tnlHost, int >, double, int >* solver, int sweep, int i);
+__global__ void runCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i);
+//__global__ void runCUDA(tnlNarrowBand< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i);
 
-__global__ void initCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver);
+__global__ void initCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver);
 
-__global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver);
-__global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver);
-__global__ void initSetupGrid1_2CUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver);
-__global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, double tau);
-//__global__ void initCUDA(tnlNarrowBand< tnlGrid< 3,double, tnlHost, int >, double, int >* solver);
+__global__ void initSetupGridCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver);
+__global__ void initSetupGrid2CUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver);
+__global__ void initSetupGrid1_2CUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver);
+__global__ void runNarrowBandCUDA(tnlNarrowBand< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, double tau);
+//__global__ void initCUDA(tnlNarrowBand< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver);
 #endif
 
 
diff --git a/src/TNL/Containers/Algorithms/VectorOperationsCuda_impl.h b/src/TNL/Containers/Algorithms/VectorOperationsCuda_impl.h
index e1d5ce50b7c7cb8d4e74c9e02628bd61ccc06d18..705369fe0ab0d2a3baca0fb7e7f1964da2372c10 100644
--- a/src/TNL/Containers/Algorithms/VectorOperationsCuda_impl.h
+++ b/src/TNL/Containers/Algorithms/VectorOperationsCuda_impl.h
@@ -413,12 +413,10 @@ VectorOperations< Devices::Cuda >::
 vectorScalarMultiplication( Vector& v,
                             const typename Vector::RealType& alpha )
 {
-   typedef typename Vector::RealType Real;
-   typedef typename Vector::IndexType Index;
-
    TNL_ASSERT_GT( v.getSize(), 0, "Vector size must be positive." );
 
 #ifdef HAVE_CUDA
+   typedef typename Vector::IndexType Index;   
    dim3 blockSize( 0 ), gridSize( 0 );
    const Index& size = v.getSize();
    blockSize.x = 256;
@@ -497,13 +495,12 @@ addVector( Vector1& y,
            const typename Vector2::RealType& alpha,
            const typename Vector1::RealType& thisMultiplicator )
 {
-   typedef typename Vector1::RealType Real;
-   typedef typename Vector1::IndexType Index;
-
    TNL_ASSERT_GT( x.getSize(), 0, "Vector size must be positive." );
    TNL_ASSERT_EQ( x.getSize(), y.getSize(), "The vector sizes must be the same." );
 
 #ifdef HAVE_CUDA
+   typedef typename Vector1::IndexType Index;
+   
    dim3 blockSize( 0 ), gridSize( 0 );
 
    const Index& size = x.getSize();
@@ -566,14 +563,12 @@ addVectors( Vector1& v,
             const typename Vector3::RealType& multiplicator2,
             const typename Vector1::RealType& thisMultiplicator )
 {
-   typedef typename Vector1::RealType Real;
-   typedef typename Vector1::IndexType Index;
-
    TNL_ASSERT_GT( v.getSize(), 0, "Vector size must be positive." );
    TNL_ASSERT_EQ( v.getSize(), v1.getSize(), "The vector sizes must be the same." );
    TNL_ASSERT_EQ( v.getSize(), v2.getSize(), "The vector sizes must be the same." );
 
 #ifdef HAVE_CUDA
+   typedef typename Vector1::IndexType Index;   
    dim3 blockSize( 0 ), gridSize( 0 );
 
    const Index& size = v.getSize();
diff --git a/src/TNL/Containers/Algorithms/VectorOperationsHost_impl.h b/src/TNL/Containers/Algorithms/VectorOperationsHost_impl.h
index 492153c729e23a6022f441aee75ddcb0f4a5dabf..b103340c26ee1d5c4b8135a3e53224ca823e6fe4 100644
--- a/src/TNL/Containers/Algorithms/VectorOperationsHost_impl.h
+++ b/src/TNL/Containers/Algorithms/VectorOperationsHost_impl.h
@@ -444,7 +444,6 @@ VectorOperations< Devices::Host >::
 vectorScalarMultiplication( Vector& v,
                             const typename Vector::RealType& alpha )
 {
-   typedef typename Vector::RealType Real;
    typedef typename Vector::IndexType Index;
 
    TNL_ASSERT_GT( v.getSize(), 0, "Vector size must be positive." );
@@ -528,7 +527,6 @@ addVector( Vector1& y,
            const typename Vector2::RealType& alpha,
            const typename Vector1::RealType& thisMultiplicator )
 {
-   typedef typename Vector1::RealType Real;
    typedef typename Vector1::IndexType Index;
 
    TNL_ASSERT_GT( x.getSize(), 0, "Vector size must be positive." );
@@ -540,6 +538,7 @@ addVector( Vector1& y,
 #ifdef __GNUC__
    // We need to get the address of the first element to avoid
    // bounds checking in TNL::Array::operator[]
+   typedef typename Vector1::RealType Real;   
          Real* Y = y.getData();
    const Real* X = x.getData();
 #endif
@@ -598,7 +597,6 @@ addVectors( Vector1& v,
             const typename Vector3::RealType& multiplicator2,
             const typename Vector1::RealType& thisMultiplicator )
 {
-   typedef typename Vector1::RealType Real;
    typedef typename Vector1::IndexType Index;
 
    TNL_ASSERT_GT( v.getSize(), 0, "Vector size must be positive." );
diff --git a/src/TNL/Containers/Array_impl.h b/src/TNL/Containers/Array_impl.h
index 2612bcc25bdfdb7354490f429840537aed4b11a2..52504e47a31f057adf878c86e32891a1210daa7c 100644
--- a/src/TNL/Containers/Array_impl.h
+++ b/src/TNL/Containers/Array_impl.h
@@ -354,7 +354,9 @@ Array< Element, Device, Index >&
 Array< Element, Device, Index >::
 operator = ( const Array< Element, Device, Index >& array )
 {
-   TNL_ASSERT_EQ( array.getSize(), this->getSize(), "Array sizes must be the same." );
+   //TNL_ASSERT_EQ( array.getSize(), this->getSize(), "Array sizes must be the same." );
+   if( this->getSize() != array.getSize() )
+      this->setLike( array );
    if( this->getSize() > 0 )
       Algorithms::ArrayOperations< Device >::
          template copyMemory< Element,
@@ -374,7 +376,9 @@ Array< Element, Device, Index >&
 Array< Element, Device, Index >::
 operator = ( const ArrayT& array )
 {
-   TNL_ASSERT_EQ( array.getSize(), this->getSize(), "Array sizes must be the same." );
+   //TNL_ASSERT_EQ( array.getSize(), this->getSize(), "Array sizes must be the same." );
+   if( this->getSize() != array.getSize() )
+      this->setLike( array );   
    if( this->getSize() > 0 )
       Algorithms::ArrayOperations< Device, typename ArrayT::DeviceType >::
          template copyMemory< Element,
@@ -491,7 +495,10 @@ load( File& file )
       return false;
    Index _size;
    if( ! file.read( &_size ) )
+   {
+      std::cerr << "Unable to read the array size." << std::endl;
       return false;
+   }
    if( _size < 0 )
    {
       std::cerr << "Error: The size " << _size << " of the file is not a positive number or zero." << std::endl;
diff --git a/src/TNL/Containers/List_impl.h b/src/TNL/Containers/List_impl.h
index 95c73269aa6ff04d9cf4f282748f743f1e803d2c..230719b1dfb95f073e4ea817a0b8fb37a3aefb77 100644
--- a/src/TNL/Containers/List_impl.h
+++ b/src/TNL/Containers/List_impl.h
@@ -1,5 +1,5 @@
 /***************************************************************************
-                          tnlList_impl.h  -  description
+                          TNL::Containers::List_impl.h  -  description
                              -------------------
     begin                : Mar, 5 Apr 2016 12:46 PM
     copyright            : (C) 2016 by Tomas Oberhuber
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Eikonal/godunovEikonal.h b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Eikonal/godunovEikonal.h
index 4e79b48f67944f970dcbab911711bbd533825ebd..9cebf0ea9a2608526cff1dcf1290e41fa4f7b083 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Eikonal/godunovEikonal.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Eikonal/godunovEikonal.h
@@ -37,7 +37,7 @@ class godunovEikonalScheme< tnlGrid< 1,MeshReal, Device, MeshIndex >, Real, Inde
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 1, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
       static String getType();
@@ -99,7 +99,7 @@ class godunovEikonalScheme< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Inde
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 2, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -180,7 +180,7 @@ class godunovEikonalScheme< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Inde
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 3, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;   
 
       static String getType();
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/godunovEikonal.h b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/godunovEikonal.h
index d4941ba38a992fd76047a7ead2f538b79af2f079..16dc991dd615244efb991dcb5d04c869ac7f6477 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/godunovEikonal.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/godunovEikonal.h
@@ -20,7 +20,7 @@
 #include <solvers/preconditioners/tnlDummyPreconditioner.h>
 #include <solvers/tnlSolverMonitor.h>
 #include <core/tnlLogger.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <core/vectors/tnlSharedVector.h>
 #include <core/mfilename.h>
 #include <mesh/tnlGrid.h>
@@ -49,7 +49,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 1, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -98,7 +98,7 @@ class godunovEikonalScheme< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Inde
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 2, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
       static String getType();
@@ -143,7 +143,7 @@ class godunovEikonalScheme< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Inde
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 3, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovEikonal.h b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovEikonal.h
index 1ad2a4dea10013ee608abe2e5a147b0317c69500..69596e0f246f3b80ff8c8f9aaaa86a3776c9fece 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovEikonal.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovEikonal.h
@@ -21,7 +21,7 @@
 #include <solvers/preconditioners/tnlDummyPreconditioner.h>
 #include <solvers/tnlSolverMonitor.h>
 #include <core/tnlLogger.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <core/vectors/tnlSharedVector.h>
 #include <core/mfilename.h>
 #include <mesh/tnlGrid.h>
@@ -52,7 +52,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 1, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -121,7 +121,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 #ifdef HAVE_CUDA
    __device__ __host__
@@ -197,7 +197,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 #ifdef HAVE_CUDA
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovMap.h b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovMap.h
index d70d360231eb60edc47c4c1037646851a0f46dc8..95276e7998571f0119357cbba027599a97fca276 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovMap.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/Godunov-Eikonal/parallelGodunovMap.h
@@ -21,7 +21,7 @@
 #include <solvers/preconditioners/tnlDummyPreconditioner.h>
 #include <solvers/tnlSolverMonitor.h>
 #include <core/tnlLogger.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <core/vectors/tnlSharedVector.h>
 #include <core/mfilename.h>
 #include <mesh/tnlGrid.h>
@@ -52,7 +52,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 1, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -121,7 +121,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 #ifdef HAVE_CUDA
    __device__ __host__
@@ -199,7 +199,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 #ifdef HAVE_CUDA
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/upwindEikonal.h b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/upwindEikonal.h
index 009c21c4cd546951c2724d24e1f9573d3cb7cdd8..7467f66715e6dcc5fcf9a04ac0ecd66d1f026a11 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/upwindEikonal.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Operators/Hamilton-Jacobi/upwindEikonal.h
@@ -41,7 +41,7 @@ class upwindEikonalScheme< tnlGrid< 1,MeshReal, Device, MeshIndex >, Real, Index
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 1, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -94,7 +94,7 @@ class upwindEikonalScheme< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 2, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
       static String getType();
@@ -151,7 +151,7 @@ class upwindEikonalScheme< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index
       typedef Device DeviceType;
       typedef Index IndexType;
       typedef tnlGrid< 3, Real, Device, Index > MeshType;
-      typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+      typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
       typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/main.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/main.h
index c7a3d34cb0ddffcf31548f70ae6cdd74dcce4a7c..813e99c72cd6de9b597c10a943ab666220453018 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/main.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/main.h
@@ -47,7 +47,7 @@ int main( int argc, char* argv[] )
 
    if(dim == 2)
    {
-		tnlFastSweepingMap<tnlGrid<2,double,tnlHost, int>, double, int> solver;
+		tnlFastSweepingMap<tnlGrid<2,double,TNL::Devices::Host, int>, double, int> solver;
 		if(!solver.init(parameters))
 	   {
 			cerr << "Solver failed to initialize." << endl;
@@ -60,7 +60,7 @@ int main( int argc, char* argv[] )
    }
 //   else if(dim == 3)
 //   {
-//		tnlFastSweepingMap<tnlGrid<3,double,tnlHost, int>, double, int> solver;
+//		tnlFastSweepingMap<tnlGrid<3,double,TNL::Devices::Host, int>, double, int> solver;
 //		if(!solver.init(parameters))
 //	   {
 //			cerr << "Solver failed to initialize." << endl;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap.h
index f6f2150a4df3ae6f0fe85efeeb6ea9abac0007d3..c568329ba2aa5fdb8fed303d43b25e73f210c014 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap.h
@@ -17,10 +17,10 @@
 #define TNLFASTSWEEPING_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
 #include <functions/tnlMeshFunction.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 #include <limits.h>
@@ -55,7 +55,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -137,7 +137,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 	tnlFastSweepingMap();
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap2D_CUDA_v4_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap2D_CUDA_v4_impl.h
index 7b4024b26432bd7c8ae474e84f509660d64fcee6..b865bef001413f10d2862e0c16675e5c33309470 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap2D_CUDA_v4_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap2D_CUDA_v4_impl.h
@@ -207,7 +207,7 @@ template< typename MeshReal,
 __device__
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: updateValue( Index i, Index j, Index* something_changed)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 
@@ -269,7 +269,7 @@ bool tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index >
 	int i = threadIdx.x + blockDim.x*blockIdx.x;
 	int j = blockDim.y*blockIdx.y + threadIdx.y;
 
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
@@ -397,7 +397,7 @@ Real tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index >
 
 
 
-__global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i, int* changed)
+__global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i, int* changed)
 {
 
 	__shared__ int something_changed;
@@ -512,7 +512,7 @@ __global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 2,double, tnlHost, int >, d
 }
 
 
-__global__ void initCUDA(tnlFastSweepingMap< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweepingMap< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 
 
@@ -571,7 +571,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -590,7 +590,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -609,7 +609,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -642,7 +642,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -675,7 +675,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -708,7 +708,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -742,7 +742,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -775,7 +775,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -808,7 +808,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -841,7 +841,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -878,7 +878,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -911,7 +911,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -944,7 +944,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -968,7 +968,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1001,7 +1001,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -1034,7 +1034,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweepingMap< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap_CUDA.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap_CUDA.h
index aa606ea47999ce739fb23d7ad6c38f937cf23f26..a23057e78c745e74467db4c4190d6f217024bc5a 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap_CUDA.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping-map/tnlFastSweepingMap_CUDA.h
@@ -17,9 +17,9 @@
 #define TNLFASTSWEEPING_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 
@@ -54,7 +54,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 	tnlFastSweepingMap();
@@ -131,7 +131,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -176,11 +176,11 @@ protected:
 
 #ifdef HAVE_CUDA
 //template<int sweep_t>
-__global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i, int* changed);
-//__global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 3,double, tnlHost, int >, double, int >* solver, int sweep, int i);
+__global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i, int* changed);
+//__global__ void runCUDA(tnlFastSweepingMap< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i);
 
-__global__ void initCUDA(tnlFastSweepingMap< tnlGrid< 2,double, tnlHost, int >, double, int >* solver);
-//__global__ void initCUDA(tnlFastSweepingMap< tnlGrid< 3,double, tnlHost, int >, double, int >* solver);
+__global__ void initCUDA(tnlFastSweepingMap< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver);
+//__global__ void initCUDA(tnlFastSweepingMap< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver);
 #endif
 
 /*various implementtions.... choose one*/
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/main.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/main.h
index b5a27952f3d138ab66fc58849f335722b3181333..8aca8f1b81db3a4a462d231e4724f3eef50e4723 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/main.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/main.h
@@ -47,7 +47,7 @@ int main( int argc, char* argv[] )
 
    if(dim == 2)
    {
-		tnlFastSweeping<tnlGrid<2,double,tnlHost, int>, double, int> solver;
+		tnlFastSweeping<tnlGrid<2,double,TNL::Devices::Host, int>, double, int> solver;
 		if(!solver.init(parameters))
 	   {
 			cerr << "Solver failed to initialize." << endl;
@@ -60,7 +60,7 @@ int main( int argc, char* argv[] )
    }
    else if(dim == 3)
    {
-		tnlFastSweeping<tnlGrid<3,double,tnlHost, int>, double, int> solver;
+		tnlFastSweeping<tnlGrid<3,double,TNL::Devices::Host, int>, double, int> solver;
 		if(!solver.init(parameters))
 	   {
 			cerr << "Solver failed to initialize." << endl;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping.h
index 55f145af95deba9752fb5c4f564bbbfa331bf153..96d26db7b5a2077d8e2199292f0e888b0171a5c2 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping.h
@@ -17,10 +17,10 @@
 #define TNLFASTSWEEPING_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
 #include <functions/tnlMeshFunction.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 #include <limits.h>
@@ -55,7 +55,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -135,7 +135,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 	tnlFastSweeping();
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_impl.h
index f5a0088cf430e701d23b49b3042c649923da881a..21e45020924453bd5676650028f6fdc6b210bc42 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_impl.h
@@ -439,7 +439,7 @@ Real tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > ::
 
 
 
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
@@ -502,7 +502,7 @@ __global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, doub
 }
 
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
 	int gy = blockDim.y*blockIdx.y + threadIdx.y;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v2_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v2_impl.h
index 1320037227b995c493f6f6f5a94c54a8fd1f081b..1d4ee11b0856a618f7b86452828a2bd2e71bf24e 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v2_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v2_impl.h
@@ -438,7 +438,7 @@ Real tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > ::
 
 
 
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 {
 
 	//int gx = threadIdx.x;
@@ -568,7 +568,7 @@ __global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, doub
 }
 
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
 	int gy = blockDim.y*blockIdx.y + threadIdx.y;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v3_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v3_impl.h
index 916593edac041c9221acaedf585ff18865da2855..fad3bc293cc452c1fe077e892b13b51750a78bb6 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v3_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v3_impl.h
@@ -492,7 +492,7 @@ Real tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > ::
 
 
 template<>
-__global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+__global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 {
 
 	if(blockIdx.x+blockIdx.y == k)
@@ -516,7 +516,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 			/*---------------------------------------------------------------------------------------------------------------------------*/
 }
 	template<>
-	__global__ void runCUDA<2>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+	__global__ void runCUDA<2>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 	{
 	if((gridDim.x - blockIdx.x - 1)+blockIdx.y == k)
 	{
@@ -538,7 +538,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 	}
 	}			/*---------------------------------------------------------------------------------------------------------------------------*/
 	template<>
-	__global__ void runCUDA<4>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+	__global__ void runCUDA<4>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 	{
 	if(blockIdx.x+blockIdx.y == gridDim.x+gridDim.y-k-2)
 		{
@@ -564,7 +564,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 	}
 
 	template<>
-	__global__ void runCUDA<8>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+	__global__ void runCUDA<8>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 	{
 	if((gridDim.x - blockIdx.x - 1)+blockIdx.y == gridDim.x+gridDim.y-k-2)
 		{
@@ -595,7 +595,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 
 
 	template<>
-		__global__ void runCUDA<5>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+		__global__ void runCUDA<5>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 		{
 
 			if(blockIdx.x+blockIdx.y == k)
@@ -640,7 +640,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 
 
 	template<>
-		__global__ void runCUDA<10>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+		__global__ void runCUDA<10>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 		{
 			if((gridDim.x - blockIdx.x - 1)+blockIdx.y == k)
 			{
@@ -687,7 +687,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 
 
 	template<>
-	__global__ void runCUDA<15>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+	__global__ void runCUDA<15>(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 	{
 
 		if(blockIdx.x+blockIdx.y == k)
@@ -806,7 +806,7 @@ __global__ void runCUDA<1>(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, d
 
 
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
 	int gy = blockDim.y*blockIdx.y + threadIdx.y;
@@ -825,7 +825,7 @@ __global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, dou
 
 
 
-//__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int k)
+//__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int k)
 //{
 //
 //	if(sweep==1 && blockIdx.x+blockIdx.y == k)
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v4_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v4_impl.h
index d52845d1f656a6062cb059be0a21851d2ace86cf..dd878d4953d46ac014622a43e03ad0e8b6b9a492 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v4_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v4_impl.h
@@ -176,7 +176,7 @@ template< typename MeshReal,
 __device__
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: updateValue( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -227,7 +227,7 @@ bool tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > ::
 	int i = threadIdx.x + blockDim.x*blockIdx.x;
 	int j = blockDim.y*blockIdx.y + threadIdx.y;
 
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
@@ -357,7 +357,7 @@ Real tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > ::
 
 
 
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 	int gx = 0;
@@ -464,7 +464,7 @@ __global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, doub
 }
 
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 
 
@@ -523,7 +523,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1111( Index i, Index j)
 {
-//	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+//	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 //	Entity.setCoordinates(CoordinatesType(i,j));
 //	Entity.refresh();
 //	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -542,7 +542,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0000( Index i, Index j)
 {
-//	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+//	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 //	Entity.setCoordinates(CoordinatesType(i,j));
 //	Entity.refresh();
 //	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -561,7 +561,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -594,7 +594,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -627,7 +627,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -660,7 +660,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0111( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -694,7 +694,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -727,7 +727,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -760,7 +760,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -793,7 +793,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1000( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -830,7 +830,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1100( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -863,7 +863,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1010( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -896,7 +896,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare1001( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -920,7 +920,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0011( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -953,7 +953,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0101( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
@@ -986,7 +986,7 @@ template< typename MeshReal,
           typename Index >
 void tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > :: setupSquare0110( Index i, Index j)
 {
-	tnlGridEntity< tnlGrid< 2,double, tnlHost, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 2,double, TNL::Devices::Host, int >, 2, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 2, tnlGridEntityNoStencilStorage >,2> neighborEntities(Entity);
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v5_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v5_impl.h
index 7c28c805e697fc8954564a0fc4ff69eeecf285b1..4f8c90df818b7121663e8053b0577d37fd119055 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v5_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping2D_CUDA_v5_impl.h
@@ -491,7 +491,7 @@ Real tnlFastSweeping< tnlGrid< 2,MeshReal, Device, MeshIndex >, Real, Index > ::
 
 
 
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 	extern __shared__ double u[];
@@ -677,7 +677,7 @@ __global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, doub
 }
 
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
 	int gy = blockDim.y*blockIdx.y + threadIdx.y;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping3D_CUDA_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping3D_CUDA_impl.h
index 62aae3112b4ecd0ccb9b157bed957bf14b07f0dd..51bb61716d1b49e2921d974e93692cdf0a34ddf1 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping3D_CUDA_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping3D_CUDA_impl.h
@@ -165,7 +165,7 @@ template< typename MeshReal,
 __device__
 void tnlFastSweeping< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > :: updateValue( Index i, Index j, Index k)
 {
-	tnlGridEntity< tnlGrid< 3,double, tnlHost, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 3,double, TNL::Devices::Host, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j,k));
 	Entity.refresh();
 	tnlNeighborGridEntityGetter<tnlGridEntity< MeshType, 3, tnlGridEntityNoStencilStorage >,3> neighborEntities(Entity);
@@ -222,7 +222,7 @@ template< typename MeshReal,
 __device__
 bool tnlFastSweeping< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > :: initGrid(int i, int j, int k)
 {
-	tnlGridEntity< tnlGrid< 3,double, tnlHost, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
+	tnlGridEntity< tnlGrid< 3,double, TNL::Devices::Host, int >, 3, tnlGridEntityNoStencilStorage > Entity(Mesh);
 	Entity.setCoordinates(CoordinatesType(i,j,k));
 	Entity.refresh();
 	int gid = Entity.getIndex();
@@ -255,7 +255,7 @@ Real tnlFastSweeping< tnlGrid< 3,MeshReal, Device, MeshIndex >, Real, Index > ::
 
 
 
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 3,double, tnlHost, int >, double, int >* solver, int sweep, int i)
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i)
 {
 
 	int gx = 0;
@@ -474,7 +474,7 @@ __global__ void runCUDA(tnlFastSweeping< tnlGrid< 3,double, tnlHost, int >, doub
 }
 
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 3,double, tnlHost, int >, double, int >* solver)
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver)
 {
 	int gx = threadIdx.x + blockDim.x*blockIdx.x;
 	int gy = blockDim.y*blockIdx.y + threadIdx.y;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping_CUDA.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping_CUDA.h
index 310cdf3f3028eb68e71c31d821ffea1538615724..f531da431bfec5d16da8ea7deabe6595031a0873 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping_CUDA.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/fast-sweeping/tnlFastSweeping_CUDA.h
@@ -17,9 +17,9 @@
 #define TNLFASTSWEEPING_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 
@@ -54,7 +54,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 2, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 	tnlFastSweeping();
@@ -129,7 +129,7 @@ public:
 	typedef Device DeviceType;
 	typedef Index IndexType;
 	typedef tnlGrid< 3, Real, Device, Index > MeshType;
-	typedef tnlVector< RealType, DeviceType, IndexType> DofVectorType;
+	typedef TNL::Containers::Vector< RealType, DeviceType, IndexType> DofVectorType;
 	typedef typename MeshType::CoordinatesType CoordinatesType;
 
 
@@ -174,11 +174,11 @@ protected:
 
 #ifdef HAVE_CUDA
 //template<int sweep_t>
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver, int sweep, int i);
-__global__ void runCUDA(tnlFastSweeping< tnlGrid< 3,double, tnlHost, int >, double, int >* solver, int sweep, int i);
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i);
+__global__ void runCUDA(tnlFastSweeping< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver, int sweep, int i);
 
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, tnlHost, int >, double, int >* solver);
-__global__ void initCUDA(tnlFastSweeping< tnlGrid< 3,double, tnlHost, int >, double, int >* solver);
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 2,double, TNL::Devices::Host, int >, double, int >* solver);
+__global__ void initCUDA(tnlFastSweeping< tnlGrid< 3,double, TNL::Devices::Host, int >, double, int >* solver);
 #endif
 
 /*various implementtions.... choose one*/
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/main.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/main.h
index 55a0942f81721c1e90e670aaae18941471f6e13d..911f0a29d9f0a82aa73767cfcee407093bc2eb66 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/main.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/main.h
@@ -41,24 +41,24 @@ int main( int argc, char* argv[] )
 
 
 	tnlDeviceEnum device;
-	device = tnlHostDevice;
+	device = TNL::Devices::HostDevice;
 
 	const int& dim = parameters.getParameter< int >( "dim" );
 
 	if(dim == 2)
 	{
 
-	   typedef parallelGodunovMapScheme< tnlGrid<2,double,tnlHost, int>, double, int > SchemeTypeHost;
+	   typedef parallelGodunovMapScheme< tnlGrid<2,double,TNL::Devices::Host, int>, double, int > SchemeTypeHost;
 /*#ifdef HAVE_CUDA
 		   typedef parallelGodunovMapScheme< tnlGrid<2,double,tnlCuda, int>, double, int > SchemeTypeDevice;
 #endif
 #ifndef HAVE_CUDA*/
-	   typedef parallelGodunovMapScheme< tnlGrid<2,double,tnlHost, int>, double, int > SchemeTypeDevice;
+	   typedef parallelGodunovMapScheme< tnlGrid<2,double,TNL::Devices::Host, int>, double, int > SchemeTypeDevice;
 /*#endif*/
 
-	   if(device==tnlHostDevice)
+	   if(device==TNL::Devices::HostDevice)
 	   {
-		   typedef tnlHost Device;
+		   typedef TNL::Devices::Host Device;
 
 
 		   tnlParallelMapSolver<2,SchemeTypeHost,SchemeTypeDevice, Device> solver;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver.h
index f18fbef8c668870bcf7ac15947c92947e8a35eff..400e163c9dcc8d536a478a0952aabf8ccbb1a2d8 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver.h
@@ -18,10 +18,10 @@
 #define TNLPARALLELMAPSOLVER_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
 #include <functions/tnlMeshFunction.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 #include <limits.h>
@@ -52,13 +52,13 @@ public:
 	typedef SchemeDevice SchemeTypeDevice;
 	typedef SchemeHost SchemeTypeHost;
 	typedef Device DeviceType;
-	typedef tnlVector< double, tnlHost, int > VectorType;
-	typedef tnlVector< int, tnlHost, int > IntVectorType;
-	typedef tnlGrid< 2, double, tnlHost, int > MeshType;
+	typedef TNL::Containers::Vector< double, TNL::Devices::Host, int > VectorType;
+	typedef TNL::Containers::Vector< int, TNL::Devices::Host, int > IntVectorType;
+	typedef tnlGrid< 2, double, TNL::Devices::Host, int > MeshType;
 #ifdef HAVE_CUDA
-	typedef tnlVector< double, tnlHost, int > VectorTypeCUDA;
-	typedef tnlVector< int, tnlHost, int > IntVectorTypeCUDA;
-	typedef tnlGrid< 2, double, tnlHost, int > MeshTypeCUDA;
+	typedef TNL::Containers::Vector< double, TNL::Devices::Host, int > VectorTypeCUDA;
+	typedef TNL::Containers::Vector< int, TNL::Devices::Host, int > IntVectorTypeCUDA;
+	typedef tnlGrid< 2, double, TNL::Devices::Host, int > MeshTypeCUDA;
 #endif
 	tnlParallelMapSolver();
 	bool init( const Config::ParameterContainer& parameters );
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver2D_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver2D_impl.h
index 3604f31f96abffae680554495731bb14992abda1..2925df8b236f13b0d3c460b5ca9b2970ff730e80 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver2D_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel-map/tnlParallelMapSolver2D_impl.h
@@ -31,7 +31,7 @@
 template< typename SchemeHost, typename SchemeDevice, typename Device>
 tnlParallelMapSolver<2,SchemeHost, SchemeDevice, Device, double, int>::tnlParallelMapSolver()
 {
-	this->device = tnlHostDevice;  /////////////// tnlCuda Device --- vypocet na GPU, tnlHostDevice   ---    vypocet na CPU
+	this->device = TNL::Devices::HostDevice;  /////////////// tnlCuda Device --- vypocet na GPU, TNL::Devices::HostDevice   ---    vypocet na CPU
 
 #ifdef HAVE_CUDA
 	if(this->device == tnlCudaDevice)
@@ -143,7 +143,7 @@ bool tnlParallelMapSolver<2,SchemeHost, SchemeDevice, Device, double, int>::init
 	}
 #endif
 
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 	{
 		VectorType tmp_map;
 		tmp_map.setSize(this->n * this->n);
@@ -201,7 +201,7 @@ bool tnlParallelMapSolver<2,SchemeHost, SchemeDevice, Device, double, int>::init
 
 
 	this->currentStep = 1;
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 		synchronize();
 #ifdef HAVE_CUDA
 	else if(this->device == tnlCudaDevice)
@@ -226,7 +226,7 @@ bool tnlParallelMapSolver<2,SchemeHost, SchemeDevice, Device, double, int>::init
 template< typename SchemeHost, typename SchemeDevice, typename Device>
 void tnlParallelMapSolver<2,SchemeHost, SchemeDevice, Device, double, int>::run()
 {
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 	{
 		while ((this->boundaryConditions.max() > 0 )/* || !end*/)
 		{
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/main.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/main.h
index 3a976dae469a77dc33559cf562f722913de0ef28..178e816aff7c1fe2c46b16529b1b263b24796879 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/main.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/main.h
@@ -42,24 +42,24 @@ int main( int argc, char* argv[] )
    //if (parameters.GetParameter <String>("scheme") == "godunov")
    //{
    tnlDeviceEnum device;
-   device = tnlHostDevice;
+   device = TNL::Devices::HostDevice;
 
    const int& dim = parameters.getParameter< int >( "dim" );
 
   if(dim == 2)
   {
 
-	   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlHost, int>, double, int > SchemeTypeHost;
+	   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,TNL::Devices::Host, int>, double, int > SchemeTypeHost;
 		/*#ifdef HAVE_CUDA
 		   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlCuda, int>, double, int > SchemeTypeDevice;
 		#endif
 		#ifndef HAVE_CUDA*/
-	   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlHost, int>, double, int > SchemeTypeDevice;
+	   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,TNL::Devices::Host, int>, double, int > SchemeTypeDevice;
 		/*#endif*/
 
-	   if(device==tnlHostDevice)
+	   if(device==TNL::Devices::HostDevice)
 	   {
-		   typedef tnlHost Device;
+		   typedef TNL::Devices::Host Device;
 
 
 		   tnlParallelEikonalSolver<2,SchemeTypeHost,SchemeTypeDevice, Device> solver;
@@ -92,17 +92,17 @@ int main( int argc, char* argv[] )
   else if(dim == 3)
   {
 
-	   typedef parallelGodunovEikonalScheme< tnlGrid<3,double,tnlHost, int>, double, int > SchemeTypeHost;
+	   typedef parallelGodunovEikonalScheme< tnlGrid<3,double,TNL::Devices::Host, int>, double, int > SchemeTypeHost;
 		/*#ifdef HAVE_CUDA
 		   typedef parallelGodunovEikonalScheme< tnlGrid<2,double,tnlCuda, int>, double, int > SchemeTypeDevice;
 		#endif
 		#ifndef HAVE_CUDA*/
-	   typedef parallelGodunovEikonalScheme< tnlGrid<3,double,tnlHost, int>, double, int > SchemeTypeDevice;
+	   typedef parallelGodunovEikonalScheme< tnlGrid<3,double,TNL::Devices::Host, int>, double, int > SchemeTypeDevice;
 		/*#endif*/
 
-	   if(device==tnlHostDevice)
+	   if(device==TNL::Devices::HostDevice)
 	   {
-		   typedef tnlHost Device;
+		   typedef TNL::Devices::Host Device;
 
 
 		   tnlParallelEikonalSolver<3,SchemeTypeHost,SchemeTypeDevice, Device> solver;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver.h
index 435b04323107902469893885ed806877a87db320..19cdd949359d4349172af820def49169146c8717 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver.h
@@ -18,10 +18,10 @@
 #define TNLPARALLELEIKONALSOLVER_H_
 
 #include <TNL/Config/ParameterContainer.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <TNL/Containers/StaticVector.h>
 #include <functions/tnlMeshFunction.h>
-#include <core/tnlHost.h>
+#include <TNL/Devices/Host.h>
 #include <mesh/tnlGrid.h>
 #include <mesh/grids/tnlGridEntity.h>
 #include <limits.h>
@@ -53,13 +53,13 @@ public:
 	typedef SchemeDevice SchemeTypeDevice;
 	typedef SchemeHost SchemeTypeHost;
 	typedef Device DeviceType;
-	typedef tnlVector< double, tnlHost, int > VectorType;
-	typedef tnlVector< int, tnlHost, int > IntVectorType;
-	typedef tnlGrid< 2, double, tnlHost, int > MeshType;
+	typedef TNL::Containers::Vector< double, TNL::Devices::Host, int > VectorType;
+	typedef TNL::Containers::Vector< int, TNL::Devices::Host, int > IntVectorType;
+	typedef tnlGrid< 2, double, TNL::Devices::Host, int > MeshType;
 #ifdef HAVE_CUDA
-	typedef tnlVector< double, tnlHost, int > VectorTypeCUDA;
-	typedef tnlVector< int, tnlHost, int > IntVectorTypeCUDA;
-	typedef tnlGrid< 2, double, tnlHost, int > MeshTypeCUDA;
+	typedef TNL::Containers::Vector< double, TNL::Devices::Host, int > VectorTypeCUDA;
+	typedef TNL::Containers::Vector< int, TNL::Devices::Host, int > IntVectorTypeCUDA;
+	typedef tnlGrid< 2, double, TNL::Devices::Host, int > MeshTypeCUDA;
 #endif
 	tnlParallelEikonalSolver();
 	bool init( const Config::ParameterContainer& parameters );
@@ -160,7 +160,7 @@ public:
 
 	//__device__ bool initCUDA( tnlParallelEikonalSolver<SchemeHost, SchemeDevice, Device, double, int >* cudaSolver);
 
-	/*__global__ void initRunCUDA(tnlParallelEikonalSolver<Scheme, double, tnlHost, int >* caller);*/
+	/*__global__ void initRunCUDA(tnlParallelEikonalSolver<Scheme, double, TNL::Devices::Host, int >* caller);*/
 
 #endif
 
@@ -180,13 +180,13 @@ public:
 		typedef SchemeDevice SchemeTypeDevice;
 		typedef SchemeHost SchemeTypeHost;
 		typedef Device DeviceType;
-		typedef tnlVector< double, tnlHost, int > VectorType;
-		typedef tnlVector< int, tnlHost, int > IntVectorType;
-		typedef tnlGrid< 3, double, tnlHost, int > MeshType;
+		typedef TNL::Containers::Vector< double, TNL::Devices::Host, int > VectorType;
+		typedef TNL::Containers::Vector< int, TNL::Devices::Host, int > IntVectorType;
+		typedef tnlGrid< 3, double, TNL::Devices::Host, int > MeshType;
 	#ifdef HAVE_CUDA
-		typedef tnlVector< double, tnlHost, int > VectorTypeCUDA;
-		typedef tnlVector< int, tnlHost, int > IntVectorTypeCUDA;
-		typedef tnlGrid< 3, double, tnlHost, int > MeshTypeCUDA;
+		typedef TNL::Containers::Vector< double, TNL::Devices::Host, int > VectorTypeCUDA;
+		typedef TNL::Containers::Vector< int, TNL::Devices::Host, int > IntVectorTypeCUDA;
+		typedef tnlGrid< 3, double, TNL::Devices::Host, int > MeshTypeCUDA;
 	#endif
 		tnlParallelEikonalSolver();
 		bool init( const Config::ParameterContainer& parameters );
@@ -284,7 +284,7 @@ public:
 
 	//__device__ bool initCUDA( tnlParallelEikonalSolver<SchemeHost, SchemeDevice, Device, double, int >* cudaSolver);
 
-	/*__global__ void initRunCUDA(tnlParallelEikonalSolver<Scheme, double, tnlHost, int >* caller);*/
+	/*__global__ void initRunCUDA(tnlParallelEikonalSolver<Scheme, double, TNL::Devices::Host, int >* caller);*/
 
 #endif
 
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver2D_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver2D_impl.h
index a547af90a22d08dedcc89ffa1178bb3b9fd757bc..8370b069d5c96bfd7921390e46e4979131f15b36 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver2D_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver2D_impl.h
@@ -25,7 +25,7 @@ template< typename SchemeHost, typename SchemeDevice, typename Device>
 tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::tnlParallelEikonalSolver()
 {
 	cout << "a" << endl;
-	this->device = tnlCudaDevice;  /////////////// tnlCuda Device --- vypocet na GPU, tnlHostDevice   ---    vypocet na CPU
+	this->device = tnlCudaDevice;  /////////////// tnlCuda Device --- vypocet na GPU, TNL::Devices::HostDevice   ---    vypocet na CPU
 
 #ifdef HAVE_CUDA
 	if(this->device == tnlCudaDevice)
@@ -146,7 +146,7 @@ bool tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::
 	}
 #endif
 
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 	{
 	for(int i = 0; i < this->subgridValues.getSize(); i++)
 	{
@@ -196,7 +196,7 @@ bool tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::
 
 
 	this->currentStep = 1;
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 		synchronize();
 #ifdef HAVE_CUDA
 	else if(this->device == tnlCudaDevice)
@@ -233,7 +233,7 @@ bool tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::
 template< typename SchemeHost, typename SchemeDevice, typename Device>
 void tnlParallelEikonalSolver<2,SchemeHost, SchemeDevice, Device, double, int>::run()
 {
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 	{
 
 	bool end = false;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver3D_impl.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver3D_impl.h
index 1c943a2b4fd1a9ca8caf3ecab65d6316157a0d21..b0871824cbf67f7ef18d276903ec84d8eeb0a109 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver3D_impl.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi-parallel/tnlParallelEikonalSolver3D_impl.h
@@ -25,7 +25,7 @@ template< typename SchemeHost, typename SchemeDevice, typename Device>
 tnlParallelEikonalSolver<3,SchemeHost, SchemeDevice, Device, double, int>::tnlParallelEikonalSolver()
 {
 	cout << "a" << endl;
-	this->device = tnlHostDevice;  /////////////// tnlCuda Device --- vypocet na GPU, tnlHostDevice   ---    vypocet na CPU
+	this->device = TNL::Devices::HostDevice;  /////////////// tnlCuda Device --- vypocet na GPU, TNL::Devices::HostDevice   ---    vypocet na CPU
 
 #ifdef HAVE_CUDA
 	if(this->device == tnlCudaDevice)
@@ -146,7 +146,7 @@ bool tnlParallelEikonalSolver<3,SchemeHost, SchemeDevice, Device, double, int>::
 	}
 #endif
 
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 	{
 #ifdef HAVE_OPENMP
 #pragma omp parallel for num_threads(4) schedule(dynamic)
@@ -209,7 +209,7 @@ bool tnlParallelEikonalSolver<3,SchemeHost, SchemeDevice, Device, double, int>::
 
 
 	this->currentStep = 1;
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 		synchronize();
 #ifdef HAVE_CUDA
 	else if(this->device == tnlCudaDevice)
@@ -246,7 +246,7 @@ bool tnlParallelEikonalSolver<3,SchemeHost, SchemeDevice, Device, double, int>::
 template< typename SchemeHost, typename SchemeDevice, typename Device>
 void tnlParallelEikonalSolver<3,SchemeHost, SchemeDevice, Device, double, int>::run()
 {
-	if(this->device == tnlHostDevice)
+	if(this->device == TNL::Devices::HostDevice)
 	{
 
 	bool end = false;
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/HamiltonJacobiProblem.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/HamiltonJacobiProblem.h
index a4a187cee17f9d8f8c43537347208d1caf10bc56..0523f7e90e1ccd0cef6683ffcff34b324305eee7 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/HamiltonJacobiProblem.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/HamiltonJacobiProblem.h
@@ -20,7 +20,7 @@
 #include <solvers/preconditioners/tnlDummyPreconditioner.h>
 #include <solvers/tnlSolverMonitor.h>
 #include <core/tnlLogger.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 #include <solvers/pde/tnlExplicitUpdater.h>
 #include <solvers/pde/tnlLinearSystemAssembler.h>
 #include <functions/tnlMeshFunction.h>
diff --git a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlFastSweepingMethod.h b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlFastSweepingMethod.h
index 671448e5bfb282d3f9b84a6c3c71e99e8c403432..59ddc2a690b5e6aecb62d5cf4c0997b7c567ff82 100644
--- a/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlFastSweepingMethod.h
+++ b/src/TNL/Experimental/Hamilton-Jacobi/Solvers/hamilton-jacobi/tnlFastSweepingMethod.h
@@ -30,13 +30,13 @@ template< typename Real,
 class tnlFastSweepingMethod< tnlGrid< 1, Real, Device, Index >, Anisotropy >
    : public tnlDirectEikonalMethodsBase< tnlGrid< 1, Real, Device, Index > >
 {
-   static_assert(  std::is_same< Device, tnlHost >::value, "The fast sweeping method works only on CPU." );
+   static_assert(  std::is_same< Device, TNL::Devices::Host >::value, "The fast sweeping method works only on CPU." );
    
    public:
       
       typedef tnlGrid< 1, Real, Device, Index > MeshType;
       typedef Real RealType;
-      typedef tnlHost DeviceType;
+      typedef TNL::Devices::Host DeviceType;
       typedef Index IndexType;
       typedef Anisotropy AnisotropyType;
       typedef tnlDirectEikonalMethodsBase< tnlGrid< 1, Real, Device, Index > > BaseType;
@@ -67,13 +67,13 @@ template< typename Real,
 class tnlFastSweepingMethod< tnlGrid< 2, Real, Device, Index >, Anisotropy >
    : public tnlDirectEikonalMethodsBase< tnlGrid< 2, Real, Device, Index > >
 {
-   static_assert(  std::is_same< Device, tnlHost >::value, "The fast sweeping method works only on CPU." );
+   static_assert(  std::is_same< Device, TNL::Devices::Host >::value, "The fast sweeping method works only on CPU." );
    
    public:
       
       typedef tnlGrid< 2, Real, Device, Index > MeshType;
       typedef Real RealType;
-      typedef tnlHost DeviceType;
+      typedef TNL::Devices::Host DeviceType;
       typedef Index IndexType;
       typedef Anisotropy AnisotropyType;
       typedef tnlDirectEikonalMethodsBase< tnlGrid< 2, Real, Device, Index > > BaseType;
@@ -104,13 +104,13 @@ template< typename Real,
 class tnlFastSweepingMethod< tnlGrid< 3, Real, Device, Index >, Anisotropy >
    : public tnlDirectEikonalMethodsBase< tnlGrid< 3, Real, Device, Index > >
 {
-   static_assert(  std::is_same< Device, tnlHost >::value, "The fast sweeping method works only on CPU." );
+   static_assert(  std::is_same< Device, TNL::Devices::Host >::value, "The fast sweeping method works only on CPU." );
    
    public:
       
       typedef tnlGrid< 3, Real, Device, Index > MeshType;
       typedef Real RealType;
-      typedef tnlHost DeviceType;
+      typedef TNL::Devices::Host DeviceType;
       typedef Index IndexType;
       typedef Anisotropy AnisotropyType;
       typedef tnlDirectEikonalMethodsBase< tnlGrid< 3, Real, Device, Index > > BaseType;
diff --git a/src/TNL/Object.cpp b/src/TNL/Object.cpp
index e9b13bd9d0761e85560e9306227497e30d748f83..b4f364139de556ce410526f1bf94e1e12b4d8e41 100644
--- a/src/TNL/Object.cpp
+++ b/src/TNL/Object.cpp
@@ -121,8 +121,16 @@ bool getObjectType( File& file, String& type )
       return false;
    }
    if( strncmp( mn, magic_number, 5 ) != 0 &&
-       strncmp( mn, "SIM33", 5 ) != 0 ) return false;
-   if( ! type. load( file ) ) return false;
+       strncmp( mn, "SIM33", 5 ) != 0 )
+   {
+       std::cout << "Not a TNL file (wrong magic number)." << std::endl;
+       return false;
+   }
+   if( ! type. load( file ) )
+   {
+       std::cerr << "Cannot load the object type." << std::endl;
+       return false;
+   }
    return true;
 }
 
diff --git a/src/TNL/Solvers/ODE/Merson_impl.h b/src/TNL/Solvers/ODE/Merson_impl.h
index 3d0b98a78184924f6b247a26a1aeb4adf67e46d2..d206087f1cc3d5ba78a888ef916b06f40e76235a 100644
--- a/src/TNL/Solvers/ODE/Merson_impl.h
+++ b/src/TNL/Solvers/ODE/Merson_impl.h
@@ -133,6 +133,11 @@ bool Merson< Problem > :: solve( DofVectorPointer& u )
       std::cerr << "No problem was set for the Merson ODE solver." << std::endl;
       return false;
    }
+   if( this->getTau() == 0.0 )
+   {
+      std::cerr << "The time step for the Merson ODE solver is zero." << std::endl;
+      return false;
+   }
    /****
     * First setup the supporting meshes k1...k5 and kAux.
     */
diff --git a/src/TNL/String.cpp b/src/TNL/String.cpp
index 4b43be36c0e53f7d9ee35a18293a45340c4c34fc..2645801139e3e464037dbe005063d29fb6c676be 100644
--- a/src/TNL/String.cpp
+++ b/src/TNL/String.cpp
@@ -70,6 +70,7 @@ void String::setSize( int size )
 {
    TNL_ASSERT_GE( size, 0, "string size must be non-negative" );
    const int _length = STRING_PAGE * ( size / STRING_PAGE + 1 );
+   TNL_ASSERT_GE( _length, 0, "_length size must be non-negative" );
    if( length != _length ) {
       if( string ) {
          delete[] string;
diff --git a/src/TNL/param-types.h b/src/TNL/param-types.h
index 790cf627abc0c58094e3113b664019130d53e825..8a13e09565dd3ede7c2414c9504b0e12ab142a5d 100644
--- a/src/TNL/param-types.h
+++ b/src/TNL/param-types.h
@@ -16,7 +16,7 @@
 namespace TNL {
 
 template< typename T >
-String getType() { return T :: getType(); };
+String getType() { return T::getType(); };
 
 template<> inline String getType< bool >() { return String( "bool" ); };
 template<> inline String getType< short int >() { return String( "short int" ); };
diff --git a/src/Tools/tnl-diff.h b/src/Tools/tnl-diff.h
index 66f707593154e63bcb9b2362c83fe40646fbb28a..dec5cf7ebab48b6f3d2de373a55f737b56725eb5 100644
--- a/src/Tools/tnl-diff.h
+++ b/src/Tools/tnl-diff.h
@@ -298,7 +298,7 @@ bool computeDifference( const MeshPointer& meshPointer, const String& objectType
        objectType == "tnlMeshFunction" )  // TODO: remove deprecated type name
       return computeDifferenceOfMeshFunctions< MeshPointer, Element, Real, Index >( meshPointer, parameters );
    if( objectType == "Containers::Vector" ||
-       objectType == "tnlVector" || objectType == "tnlSharedVector" )   // TODO: remove deprecated type name
+       objectType == "TNL::Containers::Vector" || objectType == "tnlSharedVector" )   // TODO: remove deprecated type name
       return computeDifferenceOfVectors< MeshPointer, Element, Real, Index >( meshPointer, parameters );
    return false;
 }
@@ -312,12 +312,12 @@ bool setIndexType( const MeshPointer& meshPointer,
 {
    String indexType;
    if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
-       parsedObjectType[ 0 ] == "tnlMultiVector" ||                       // TODO: remove deprecated type names
+       parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||                       // TODO: remove deprecated type names
        parsedObjectType[ 0 ] == "tnlSharedMultiVector"   )                //
       indexType = parsedObjectType[ 4 ];
    if( parsedObjectType[ 0 ] == "Containers::Vector" ||
        parsedObjectType[ 0 ] == "tnlSharedVector" ||                     // TODO: remove deprecated type names
-       parsedObjectType[ 0 ] == "tnlVector" )                            //
+       parsedObjectType[ 0 ] == "TNL::Containers::Vector" )                            //
       indexType = parsedObjectType[ 3 ];
 
    if( parsedObjectType[ 0 ] == "Functions::MeshFunction" ||
@@ -392,7 +392,7 @@ bool setElementType( const MeshPointer& meshPointer,
    String elementType;
 
    if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
-       parsedObjectType[ 0 ] == "tnlMultiVector" ||                         // TODO: remove deprecated type names
+       parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||                         // TODO: remove deprecated type names
        parsedObjectType[ 0 ] == "tnlSharedMultiVector" )                    //
       elementType = parsedObjectType[ 2 ];
    if( parsedObjectType[ 0 ] == "Functions::MeshFunction" ||
@@ -400,7 +400,7 @@ bool setElementType( const MeshPointer& meshPointer,
       elementType = parsedObjectType[ 3 ];
    if( parsedObjectType[ 0 ] == "Containers::Vector" ||
        parsedObjectType[ 0 ] == "tnlSharedVector" ||                        // TODO: remove deprecated type names
-       parsedObjectType[ 0 ] == "tnlVector" )                               //
+       parsedObjectType[ 0 ] == "TNL::Containers::Vector" )                               //
       elementType = parsedObjectType[ 1 ];
 
 
diff --git a/src/Tools/tnl-view.h b/src/Tools/tnl-view.h
index c4d25f1d5520d2886be390bc41763c1f338e45aa..98227f5b8fd02564a8fbd76dc4553bb78d393b2d 100644
--- a/src/Tools/tnl-view.h
+++ b/src/Tools/tnl-view.h
@@ -279,7 +279,7 @@ bool convertObject( const MeshPointer& meshPointer,
 
    if( parsedObjectType[ 0 ] == "Containers::Vector" ||
        parsedObjectType[ 0 ] == "tnlSharedVector" ||   // TODO: remove deprecated type names
-       parsedObjectType[ 0 ] == "tnlVector" )          //
+       parsedObjectType[ 0 ] == "TNL::Containers::Vector" )          //
    {
       using MeshType = typename MeshPointer::ObjectType;
       // FIXME: why is MeshType::GlobalIndexType not the same as Index?
@@ -294,7 +294,7 @@ bool convertObject( const MeshPointer& meshPointer,
    }
 
    if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
-       parsedObjectType[ 0 ] == "tnlMultiVector" ||      // TODO: remove deprecated type names  
+       parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||      // TODO: remove deprecated type names  
        parsedObjectType[ 0 ] == "tnlSharedMultiVector" ) //
    {
       Containers::MultiVector< Dimension, Element, Devices::Host, Index > multiVector;
@@ -321,11 +321,11 @@ bool setDimensions( const MeshPointer& meshPointer,
 {
    int dimensions( 0 );
    if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
-       parsedObjectType[ 0 ] == "tnlMultiVector" ||                     // TODO: remove deprecated type names
+       parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||                     // TODO: remove deprecated type names
        parsedObjectType[ 0 ] == "tnlSharedMultiVector" )                //
       dimensions = atoi( parsedObjectType[ 1 ]. getString() );
    if( parsedObjectType[ 0 ] == "Containers::Vector" ||
-       parsedObjectType[ 0 ] == "tnlVector" ||                          // TODO: remove deprecated type names
+       parsedObjectType[ 0 ] == "TNL::Containers::Vector" ||                          // TODO: remove deprecated type names
        parsedObjectType[ 0 ] == "tnlSharedVector" )                     //
       dimensions = 1;
    switch( dimensions )
@@ -349,12 +349,12 @@ bool setIndexType( const MeshPointer& meshPointer,
 {
    String indexType;
    if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
-       parsedObjectType[ 0 ] == "tnlMultiVector" ||                        // TODO: remove deprecated type names
+       parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||                        // TODO: remove deprecated type names
        parsedObjectType[ 0 ] == "tnlSharedMultiVector" )                   //
       indexType = parsedObjectType[ 4 ];
    if( parsedObjectType[ 0 ] == "Containers::Vector" || 
        parsedObjectType[ 0 ] == "tnlSharedVector" ||                       // TODO: remove deprecated type names
-       parsedObjectType[ 0 ] == "tnlVector" )                              //
+       parsedObjectType[ 0 ] == "TNL::Containers::Vector" )                              //
       indexType = parsedObjectType[ 3 ];
 
    if( indexType == "int" )
@@ -426,12 +426,12 @@ bool setElementType( const MeshPointer& meshPointer,
 
    // TODO: Fix this even for arrays
    if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
-       parsedObjectType[ 0 ] == "tnlMultiVector" ||                            // TODO: remove deprecated type names
+       parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||                            // TODO: remove deprecated type names
        parsedObjectType[ 0 ] == "tnlSharedMultiVector" )                       //
       elementType = parsedObjectType[ 2 ];
    if( parsedObjectType[ 0 ] == "Containers::Vector" ||
        parsedObjectType[ 0 ] == "tnlSharedVector" ||                           // TODO: remove deprecated type names
-       parsedObjectType[ 0 ] == "tnlVector" )                                  //
+       parsedObjectType[ 0 ] == "TNL::Containers::Vector" )                                  //
       elementType = parsedObjectType[ 1 ];
 
 
@@ -516,10 +516,10 @@ bool processFiles( const Config::ParameterContainer& parameters )
          }
          if( parsedObjectType[ 0 ] == "Containers::MultiVector" ||
              parsedObjectType[ 0 ] == "Containers::Vector" ||
-             parsedObjectType[ 0 ] == "tnlMultiVector" ||                     // TODO: remove deprecated type names
+             parsedObjectType[ 0 ] == "TNL::Containers::MultiVector" ||                     // TODO: remove deprecated type names
              parsedObjectType[ 0 ] == "tnlSharedMultiVector" ||               // 
              parsedObjectType[ 0 ] == "tnlSharedVector" ||                    //
-             parsedObjectType[ 0 ] == "tnlVector" )                           //
+             parsedObjectType[ 0 ] == "TNL::Containers::Vector" )                           //
             setElementType< MeshPointer >( meshPointer, inputFiles[ i ], parsedObjectType, parameters );
          if( parsedObjectType[ 0 ] == "Functions::MeshFunction" ||
              parsedObjectType[ 0 ] == "tnlMeshFunction" )                     // TODO: remove deprecated type names
diff --git a/tests/benchmarks/heat-equation-benchmark/tnlTestGrid2D.h b/tests/benchmarks/heat-equation-benchmark/tnlTestGrid2D.h
index 5a6c6b9a9ce3adda72e3832d9c853c34a7396318..c10cec6ee9fa2a2f5762bcac0c80607c2a3d8326 100644
--- a/tests/benchmarks/heat-equation-benchmark/tnlTestGrid2D.h
+++ b/tests/benchmarks/heat-equation-benchmark/tnlTestGrid2D.h
@@ -15,7 +15,7 @@
 #include <core/tnlObject.h>
 #include <core/Devices::Host.h>
 #include <TNL/Containers/StaticVector.h>
-#include <core/vectors/tnlVector.h>
+#include <TNL/Containers/Vector.h>
 
 template< int Dimension,
           typename Real = double,
@@ -214,11 +214,11 @@ template< typename Real,
           typename Index >
 String Meshes::Grid< 2, Real, Device, Index > :: getType()
 {
-   return tnlString( "Meshes::Grid< " ) +
-          tnlString( getMeshDimension() ) + ", " +
-          tnlString( ::getType< RealType >() ) + ", " +
-          tnlString( Device :: getDeviceType() ) + ", " +
-          tnlString( ::getType< IndexType >() ) + " >";
+   return TNL::String( "Meshes::Grid< " ) +
+          TNL::String( getMeshDimension() ) + ", " +
+          TNL::String( ::getType< RealType >() ) + ", " +
+          TNL::String( Device :: getDeviceType() ) + ", " +
+          TNL::String( ::getType< IndexType >() ) + " >";
 }
 
 template< typename Real,
@@ -979,11 +979,11 @@ template< typename Real,
           typename Index >
 String Meshes::Grid< 2, Real, Device, Index > :: getType()
 {
-   return tnlString( "Meshes::Grid< " ) +
-          tnlString( getMeshDimension() ) + ", " +
-          tnlString( ::getType< RealType >() ) + ", " +
-          tnlString( Device :: getDeviceType() ) + ", " +
-          tnlString( ::getType< IndexType >() ) + " >";
+   return TNL::String( "Meshes::Grid< " ) +
+          TNL::String( getMeshDimension() ) + ", " +
+          TNL::String( ::getType< RealType >() ) + ", " +
+          TNL::String( Device :: getDeviceType() ) + ", " +
+          TNL::String( ::getType< IndexType >() ) + " >";
 }
 
 template< typename Real,
diff --git a/tests/benchmarks/heat-equation-benchmark/tnlTestGridEntity.h b/tests/benchmarks/heat-equation-benchmark/tnlTestGridEntity.h
index d798bc374af4bb1fb4826ebbbbbb0056ebd09941..aa8bd8d057309b1cd48fe38d71eab7886ccc0d7c 100644
--- a/tests/benchmarks/heat-equation-benchmark/tnlTestGridEntity.h
+++ b/tests/benchmarks/heat-equation-benchmark/tnlTestGridEntity.h
@@ -53,8 +53,8 @@ class tnlTestGridEntity< Meshes::Grid< Dimension, Real, Device, Index >, Dimensi
       constexpr static int getDimension() { return meshDimension; };
       
       
-      typedef tnlStaticVector< meshDimension, IndexType > EntityOrientationType;
-      typedef tnlStaticVector< meshDimension, IndexType > EntityBasisType;
+      typedef TNL::Containers::StaticVector< meshDimension, IndexType > EntityOrientationType;
+      typedef TNL::Containers::StaticVector< meshDimension, IndexType > EntityBasisType;
       typedef tnlTestGridEntity< GridType, entityDimension, Config > ThisType;
       //typedef tnlTestNeighborGridEntitiesStorage< ThisType > NeighborGridEntitiesStorageType;