From 55858c39f655e24e222336346628154acd517d7c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Tom=C3=A1=C5=A1=20Oberhuber?= <oberhuber.tomas@gmail.com>
Date: Wed, 3 Jul 2019 22:20:40 +0200
Subject: [PATCH] Implementing distributed expression templates.

---
 src/TNL/Containers/DistributedArray.h         |   2 +-
 .../Containers/DistributedVectorExpressions.h |   3 +-
 .../DistributedVectorViewExpressions.h        |   2 +-
 src/TNL/Containers/Expressions/Comparison.h   | 261 +++-----
 .../DistributedExpressionTemplates.h          | 615 ++++++++++++++++--
 .../Expressions/ExpressionTemplates.h         | 468 +++++++++++--
 .../Containers/Expressions/StaticComparison.h | 192 ++++++
 .../Expressions/StaticExpressionTemplates.h   |  14 +-
 .../Expressions/StaticVerticalOperations.h    | 165 +++++
 .../Expressions/VerticalOperations.h          | 145 +----
 src/TNL/Containers/VectorExpressions.h        |  58 +-
 src/TNL/Containers/VectorViewExpressions.h    |  38 +-
 12 files changed, 1523 insertions(+), 440 deletions(-)
 create mode 100644 src/TNL/Containers/Expressions/StaticComparison.h
 create mode 100644 src/TNL/Containers/Expressions/StaticVerticalOperations.h

diff --git a/src/TNL/Containers/DistributedArray.h b/src/TNL/Containers/DistributedArray.h
index bfaef291ce..b02d521137 100644
--- a/src/TNL/Containers/DistributedArray.h
+++ b/src/TNL/Containers/DistributedArray.h
@@ -24,7 +24,6 @@ template< typename Value,
           typename Communicator = Communicators::MpiCommunicator >
 class DistributedArray
 {
-   using CommunicationGroup = typename Communicator::CommunicationGroup;
 public:
    using ValueType = Value;
    using DeviceType = Device;
@@ -38,6 +37,7 @@ public:
    using CudaType = DistributedArray< Value, Devices::Cuda, Index, Communicator >;
    using ViewType = DistributedArrayView< Value, Device, Index, Communicator >;
    using ConstViewType = DistributedArrayView< std::add_const_t< Value >, Device, Index, Communicator >;
+   using CommunicationGroup = typename Communicator::CommunicationGroup;
 
    DistributedArray() = default;
 
diff --git a/src/TNL/Containers/DistributedVectorExpressions.h b/src/TNL/Containers/DistributedVectorExpressions.h
index 98f425b484..e18bbd7949 100644
--- a/src/TNL/Containers/DistributedVectorExpressions.h
+++ b/src/TNL/Containers/DistributedVectorExpressions.h
@@ -678,7 +678,7 @@ lpNorm( const Containers::DistributedVector< Real, Device, Index, Communicator >
    using CommunicatorType = typename Containers::DistributedVector< Real, Device, Index, Communicator >::CommunicatorType;
    Real result = ( Real ) 0.0;
    if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
-      const Real localResult = TNL::pow( Containers::Expressions::ExpressionLpNorm( a.getLocalVectorView(), p ), p );
+      const Real localResult = Containers::Expressions::ExpressionLpNorm( a.getLocalVectorView(), p );
       CommunicatorType::template Allreduce< Real >( &localResult, &result, 1, MPI_SUM, a.getCommunicationGroup() );
    }
    return TNL::pow( result, 1.0 / p );
@@ -777,7 +777,6 @@ operator,( const Containers::DistributedVector< Real, Device, Index, Communicato
 }
 
 template< typename ET, typename Real, typename Device, typename Index, typename Communicator >
-//Real
 auto
 operator,( const ET& a, const Containers::DistributedVector< Real, Device, Index, Communicator >& b ) 
 -> decltype( TNL::sum( a * b.getLocalVectorView() ) )
diff --git a/src/TNL/Containers/DistributedVectorViewExpressions.h b/src/TNL/Containers/DistributedVectorViewExpressions.h
index 3bf694f67d..258c7d75f8 100644
--- a/src/TNL/Containers/DistributedVectorViewExpressions.h
+++ b/src/TNL/Containers/DistributedVectorViewExpressions.h
@@ -679,7 +679,7 @@ lpNorm( const Containers::DistributedVectorView< Real, Device, Index, Communicat
    using CommunicatorType = typename Containers::DistributedVectorView< Real, Device, Index, Communicator >::CommunicatorType;
    Real result = ( Real ) 0.0;
    if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
-      const Real localResult = TNL::pow( Containers::Expressions::ExpressionLpNorm( a.getLocalVectorView(), p ), p );
+      const Real localResult = Containers::Expressions::ExpressionLpNorm( a.getLocalVectorView(), p );
       CommunicatorType::template Allreduce< Real >( &localResult, &result, 1, MPI_SUM, a.getCommunicationGroup() );
    }
    return TNL::pow( result, 1.0 / p );
diff --git a/src/TNL/Containers/Expressions/Comparison.h b/src/TNL/Containers/Expressions/Comparison.h
index 6719867ebd..d96819a19d 100644
--- a/src/TNL/Containers/Expressions/Comparison.h
+++ b/src/TNL/Containers/Expressions/Comparison.h
@@ -17,254 +17,197 @@ namespace TNL {
    namespace Containers {
       namespace Expressions {
 
+////
+// Non-static comparison
 template< typename T1,
           typename T2,
           ExpressionVariableType T1Type = ExpressionVariableTypeGetter< T1 >::value,
           ExpressionVariableType T2Type = ExpressionVariableTypeGetter< T2 >::value >
-struct StaticComparison
+struct Comparison
 {
 };
 
 /////
-// Static comparison of vector expressions
+// Comparison of two vector expressions
 template< typename T1,
           typename T2 >
-struct StaticComparison< T1, T2, VectorVariable, VectorVariable >
+struct Comparison< T1, T2, VectorVariable, VectorVariable >
 {
 
-   __cuda_callable__
-   static bool EQ( const T1& a, const T2& b )
+
+   bool EQ( const T1& a, const T2& b )
    {
-      static_assert( a.getSize() == b.getSize(), "Sizes of expressions to be compared do not fit." );
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a[ i ] != b[ i ] )
-            return false;
-      return true;
+      if( a.getSize() != b.getSize() )
+         return false;
+      if( a.getSize() == 0 )
+         return true;
+
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] == b[ i ] ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool NE( const T1& a, const T2& b )
+   bool NE( const T1& a, const T2& b )
    {
       return ! EQ( a, b );
    }
 
-   __cuda_callable__
-   static bool GT( const T1& a, const T2& b )
+   bool GT( const T1& a, const T2& b )
    {
-      static_assert( a.getSize() == b.getSize(), "Sizes of expressions to be compared do not fit." );
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a[ i ] <= b[ i ] )
-            return false;
-      return true;
+      TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
+
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] > b[ i ] ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool LE( const T1& a, const T2& b )
+   bool LE( const T1& a, const T2& b )
    {
       return ! GT( a, b );
    }
 
-   __cuda_callable__
-   static bool LT( const T1& a, const T2& b )
+   bool LT( const T1& a, const T2& b )
    {
-      static_assert( a.getSize() == b.getSize(), "Sizes of expressions to be compared do not fit." );
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a[ i ] >= b[ i ] )
-            return false;
-      return true;
+      TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
+
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] < b[ i ] ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool GE( const T1& a, const T2& b )
+   bool GE( const T1& a, const T2& b )
    {
       return ! LT( a, b );
    }
 };
 
 /////
-// Static comparison of number and vector expressions
+// Comparison of number and vector expression
 template< typename T1,
           typename T2 >
-struct StaticComparison< T1, T2, ArithmeticVariable, VectorVariable >
+struct Comparison< T1, T2, ArithmeticVariable, VectorVariable >
 {
 
-   __cuda_callable__
-   static bool EQ( const T1& a, const T2& b )
+
+   bool EQ( const T1& a, const T2& b )
    {
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a != b[ i ] )
-            return false;
-      return true;
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a == b[ i ] ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool NE( const T1& a, const T2& b )
+   bool NE( const T1& a, const T2& b )
    {
       return ! EQ( a, b );
    }
 
-   __cuda_callable__
-   static bool GT( const T1& a, const T2& b )
+   bool GT( const T1& a, const T2& b )
    {
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a <= b[ i ] )
-            return false;
-      return true;
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a > b[ i ] ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool LE( const T1& a, const T2& b )
+   bool LE( const T1& a, const T2& b )
    {
       return ! GT( a, b );
    }
 
-   __cuda_callable__
-   static bool LT( const T1& a, const T2& b )
+   bool LT( const T1& a, const T2& b )
    {
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a >= b[ i ] )
-            return false;
-      return true;
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a < b[ i ] ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool GE( const T1& a, const T2& b )
+   bool GE( const T1& a, const T2& b )
    {
       return ! LT( a, b );
    }
 };
 
 /////
-// Static comparison of vector expressions and number
+// Comparison of vector expressions and number
 template< typename T1,
           typename T2 >
-struct StaticComparison< T1, T2, VectorVariable, ArithmeticVariable >
+struct Comparison< T1, T2, VectorVariable, ArithmeticVariable >
 {
 
-   __cuda_callable__
-   static bool EQ( const T1& a, const T2& b )
+
+   bool EQ( const T1& a, const T2& b )
    {
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a[ i ] != b )
-            return false;
-      return true;
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] == b ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool NE( const T1& a, const T2& b )
+   bool NE( const T1& a, const T2& b )
    {
       return ! EQ( a, b );
    }
 
-   __cuda_callable__
-   static bool GT( const T1& a, const T2& b )
+   bool GT( const T1& a, const T2& b )
    {
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a[ i ] <= b )
-            return false;
-      return true;
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] > b ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool LE( const T1& a, const T2& b )
+   bool LE( const T1& a, const T2& b )
    {
       return ! GT( a, b );
    }
 
-   __cuda_callable__
-   static bool LT( const T1& a, const T2& b )
+   bool LT( const T1& a, const T2& b )
    {
-      for( int i = 0; i < a.getSize(); i++ )
-         if( a[ i ] >= b )
-            return false;
-      return true;
+      using DeviceType = typename T1::DeviceType;
+      using IndexType = typename T1::IndexType;
+
+      auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] < b ); };
+      auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
+      auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
+      return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
    }
 
-   __cuda_callable__
-   static bool GE( const T1& a, const T2& b )
+   bool GE( const T1& a, const T2& b )
    {
       return ! LT( a, b );
    }
 };
 
 
-////
-// Non-static comparison
-template< typename T1,
-          typename T2 >
-__cuda_callable__
-bool ComparisonEQ( const T1& a, const T2& b )
-{
-   if( a.getSize() != b.getSize() )
-      return false;
-   if( a.getSize() == 0 )
-      return true;
-
-   using DeviceType = typename T1::DeviceType;
-   using IndexType = typename T1::IndexType;
-
-   auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] == b[ i ] ); };
-   auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
-   auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
-   return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
-}
-
-template< typename T1,
-          typename T2 >
-__cuda_callable__
-bool ComparisonNE( const T1& a, const T2& b )
-{
-   return ! ComparisonEQ( a, b );
-}
-
-template< typename T1,
-          typename T2 >
-__cuda_callable__
-bool ComparisonGT( const T1& a, const T2& b )
-{
-   TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
-
-   using DeviceType = typename T1::DeviceType;
-   using IndexType = typename T1::IndexType;
-
-   auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] > b[ i ] ); };
-   auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
-   auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
-   return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
-}
-
-template< typename T1,
-          typename T2 >
-__cuda_callable__
-bool ComparisonLE( const T1& a, const T2& b )
-{
-   return ! ComparisonGT( a, b );
-}
-
-template< typename T1,
-          typename T2 >
-__cuda_callable__
-bool ComparisonLT( const T1& a, const T2& b )
-{
-   TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
-
-   using DeviceType = typename T1::DeviceType;
-   using IndexType = typename T1::IndexType;
-
-   auto fetch = [=] __cuda_callable__ ( IndexType i ) -> bool { return  ( a[ i ] < b[ i ] ); };
-   auto reduction = [=] __cuda_callable__ ( bool& a, const bool& b ) { a &= b; };
-   auto volatileReduction = [=] __cuda_callable__ ( volatile bool& a, volatile bool& b ) { a &= b; };
-   return Algorithms::Reduction< DeviceType >::reduce( a.getSize(), reduction, volatileReduction, fetch, true );
-}
-
-template< typename T1,
-          typename T2 >
-__cuda_callable__
-bool ComparisonGE( const T1& a, const T2& b )
-{
-   return ! ComparisonLT( a, b );
-}
-
       } //namespace Expressions
    } // namespace Containers
 } // namespace TNL
diff --git a/src/TNL/Containers/Expressions/DistributedExpressionTemplates.h b/src/TNL/Containers/Expressions/DistributedExpressionTemplates.h
index db4df9068e..f5c9b8169d 100644
--- a/src/TNL/Containers/Expressions/DistributedExpressionTemplates.h
+++ b/src/TNL/Containers/Expressions/DistributedExpressionTemplates.h
@@ -13,7 +13,7 @@
 #include <iostream>
 #include <TNL/Containers/Expressions/ExpressionTemplatesOperations.h>
 #include <TNL/Containers/Expressions/ExpressionVariableType.h>
-#include <TNL/Containers/Expressions/Comparison.h>
+#include <TNL/Containers/Expressions/DistributedComparison.h>
 #include <TNL/Containers/Expressions/IsStatic.h>
 
 namespace TNL {
@@ -50,8 +50,11 @@ struct DistributedBinaryExpressionTemplate< T1, T2, Operation, VectorVariable, V
    using DeviceType = typename T1::DeviceType;
    using IndexType = typename T1::IndexType;
    using IsExpressionTemplate = bool;
+   using CommunicatorType = typename T1::CommunicatorType;
+   using CommunicationGroup = typename T1::CommunicationGroup;
 
    static_assert( std::is_same< typename T1::DeviceType, typename T2::DeviceType >::value, "Attempt to mix operands allocated on different device types." );
+   static_assert( std::is_same< typename T1::CommunicatorType, typename T2::CommunicatorType >::value, "Attempt to mix operands using different communicators." );
    static_assert( IsStaticType< T1 >::value == IsStaticType< T2 >::value, "Attempt to mix static and non-static operands in binary expression templates." );
    static constexpr bool isStatic() { return false; }
 
@@ -79,6 +82,12 @@ struct DistributedBinaryExpressionTemplate< T1, T2, Operation, VectorVariable, V
        return op1.getSize();
    }
 
+   CommunicationGroup getCommunicationGroup() const
+   {
+      TNL_ASSERT_EQ( op1.getCommunicationGroup(), op2.getCommunicationGroup(), "Cannot create expression from operands using different communication groups." );
+      return op1.getCommunicationGroup();
+   }
+
    protected:
       const T1 op1;
       const T2 op2;
@@ -94,6 +103,8 @@ struct DistributedBinaryExpressionTemplate< T1, T2, Operation, VectorVariable, A
    using RealType = typename T1::RealType;
    using DeviceType = typename T1::DeviceType;
    using IndexType = typename T1::IndexType;
+   using CommunicatorType = typename T1::CommunicatorType;
+   using CommunicationGroup = typename T1::CommunicationGroup;
 
    using IsExpressionTemplate = bool;
    static constexpr bool isStatic() { return false; }
@@ -122,6 +133,11 @@ struct DistributedBinaryExpressionTemplate< T1, T2, Operation, VectorVariable, A
        return op1.getSize();
    }
 
+   CommunicationGroup getCommunicationGroup() const
+   {
+      return op1.getCommunicationGroup();
+   }
+
    protected:
       const T1 op1;
       const T2 op2;
@@ -137,6 +153,8 @@ struct DistributedBinaryExpressionTemplate< T1, T2, Operation, ArithmeticVariabl
    using RealType = typename T2::RealType;
    using DeviceType = typename T2::DeviceType;
    using IndexType = typename T2::IndexType;
+   using CommunicatorType = typename T2::CommunicatorType;
+   using CommunicationGroup = typename T2::CommunicationGroup;
 
    using IsExpressionTemplate = bool;
    static constexpr bool isStatic() { return false; }
@@ -165,6 +183,11 @@ struct DistributedBinaryExpressionTemplate< T1, T2, Operation, ArithmeticVariabl
        return op2.getSize();
    }
 
+   CommunicationGroup getCommunicationGroup() const
+   {
+      return op2.getCommunicationGroup();
+   }
+
    protected:
       const T1 op1;
       const T2 op2;
@@ -186,6 +209,8 @@ struct DistributedUnaryExpressionTemplate< T1, Operation, Parameter, VectorVaria
    using DeviceType = typename T1::DeviceType;
    using IndexType = typename T1::IndexType;
    using IsExpressionTemplate = bool;
+   using CommunicatorType = typename T1::CommunicatorType;
+   using CommunicationGroup = typename T1::CommunicationGroup;
    static constexpr bool isStatic() { return false; }
 
    DistributedUnaryExpressionTemplate( const T1& a, const Parameter& p )
@@ -217,6 +242,11 @@ struct DistributedUnaryExpressionTemplate< T1, Operation, Parameter, VectorVaria
 
    const Parameter& get() { return parameter; }
 
+   CommunicationGroup getCommunicationGroup() const
+   {
+      return operand.getCommunicationGroup();
+   }
+
    protected:
       const T1 operand;
       //typename OperandType< T1, DeviceType >::type operand;
@@ -233,6 +263,8 @@ struct DistributedUnaryExpressionTemplate< T1, Operation, void, VectorVariable >
    using DeviceType = typename T1::DeviceType;
    using IndexType = typename T1::IndexType;
    using IsExpressionTemplate = bool;
+   using CommunicatorType = typename T1::CommunicatorType;
+   using CommunicationGroup = typename T1::CommunicationGroup;
    static constexpr bool isStatic() { return false; }
 
    DistributedUnaryExpressionTemplate( const T1& a ): operand( a ){}
@@ -259,6 +291,11 @@ struct DistributedUnaryExpressionTemplate< T1, Operation, void, VectorVariable >
        return operand.getSize();
    }
 
+   CommunicationGroup getCommunicationGroup() const
+   {
+      return operand.getCommunicationGroup();
+   }
+
    protected:
       const T1 operand; // TODO: fix
       //typename std::add_const< typename OperandType< T1, DeviceType >::type >::type operand;
@@ -1096,21 +1133,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator == ( const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator == ( const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator == ( const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
 }
 
 template< typename L1,
@@ -1118,11 +1198,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
 }
 
 template< typename L1,
@@ -1130,11 +1213,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::EQ( a, b );
 }
 
 ////
@@ -1145,21 +1231,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator != ( const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator != ( const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator != ( const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
 }
 
 template< typename L1,
@@ -1167,11 +1296,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
 }
 
 template< typename L1,
@@ -1179,11 +1311,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::NE( a, b );
 }
 
 ////
@@ -1194,21 +1329,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& a,
              const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator < ( const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& a,
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator < ( const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+             const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator < ( const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& a,
+             const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
 }
 
 template< typename L1,
@@ -1216,11 +1394,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >& a,
-              const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
+             const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
 }
 
 template< typename L1,
@@ -1228,11 +1409,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LT( a, b );
 }
 
 ////
@@ -1243,21 +1427,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
-             const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
+              const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& a,
-             const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
+              const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator <= ( const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator <= ( const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator <= ( const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
 }
 
 template< typename L1,
@@ -1265,11 +1492,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
 }
 
 template< typename L1,
@@ -1277,11 +1507,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::LE( a, b );
 }
 
 ////
@@ -1292,21 +1525,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& a,
              const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator > ( const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& a,
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator > ( const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+             const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator > ( const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& a,
+             const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
 }
 
 template< typename L1,
@@ -1314,11 +1590,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >& a,
              const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
 }
 
 template< typename L1,
@@ -1326,11 +1605,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
-             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GT( a, b );
 }
 
 ////
@@ -1341,21 +1623,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
               const Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator >= ( const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator >= ( const typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator >= ( const typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
 }
 
 template< typename L1,
@@ -1363,11 +1688,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::DistributedBinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
 }
 
 template< typename L1,
@@ -1375,11 +1703,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::DistributedUnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::DistributedUnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::DistributedComparison< Left, Right >::GE( a, b );
 }
 
 ////
@@ -1940,7 +2271,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 min( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
 {
-   return ExpressionMin( a );
+   using CommunicatorType = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = std::numeric_limits< Real >::max();
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionMin( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_MIN, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -1949,7 +2286,35 @@ template< typename L1,
 typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 min( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
 {
-   return ExpressionMin( a );
+   using CommunicatorType = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = std::numeric_limits< Real >::max();
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionMin( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_MIN, a.getCommunicationGroup() );
+   }
+   return result;
+}
+
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation,
+          typename Index >
+auto
+argMin( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a, Index& arg ) -> decltype( ExpressionArgMin( a, arg ) )
+{
+   throw Exceptions::NotImplementedError( "agrMin for distributed vector view is not implemented yet." );
+   return ExpressionArgMin( a, arg );
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter,
+          typename Index >
+auto
+argMin( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a, Index& arg ) -> decltype( ExpressionMin( a, arg ) )
+{
+   throw Exceptions::NotImplementedError( "agrMin for distributed vector view is not implemented yet." );
+   return ExpressionArgMin( a );
 }
 
 template< typename L1,
@@ -1958,7 +2323,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 max( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
 {
-   return ExpressionMax( a );
+   using CommunicatorType = Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = std::numeric_limits< Real >::min();
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionMax( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_MAX, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -1967,7 +2338,35 @@ template< typename L1,
 typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 max( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
 {
-   return ExpressionMax( a );
+   using CommunicatorType = Containers::Expressions::DistributedUnaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = std::numeric_limits< Real >::min();
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionMax( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_MAX, a.getCommunicationGroup() );
+   }
+   return result;
+}
+
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation,
+          typename Index >
+auto
+argMax( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a, Index& arg ) -> decltype( ExpressionArgMax( a, arg ) )
+{
+   throw Exceptions::NotImplementedError( "agrMax for distributed vector view is not implemented yet." );
+   return ExpressionArgMax( a, arg );
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter,
+          typename Index >
+auto
+argMax( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a, Index& arg ) -> decltype( ExpressionMax( a, arg ) )
+{
+   throw Exceptions::NotImplementedError( "agrMax for distributed vector view is not implemented yet." );
+   return ExpressionArgMax( a );
 }
 
 template< typename L1,
@@ -1976,7 +2375,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 sum( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
 {
-   return ExpressionSum( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = 0.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionSum( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_SUM, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -1985,7 +2390,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 sum( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
 {
-   return ExpressionSum( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = 0.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionSum( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_SUM, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -1995,7 +2406,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 lpNorm( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a, const Real& p )
 {
-   return ExpressionLpNorm( a, p );
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = ( Real ) 0.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = TNL::pow( Containers::Expressions::ExpressionLpNorm( a, p ), p );
+      CommunicatorType::template Allreduce< Real >( &localResult, &result, 1, MPI_SUM, a.getCommunicationGroup() );
+   }
+   return TNL::pow( result, 1.0 / p );
 }
 
 template< typename L1,
@@ -2005,7 +2422,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 lpNorm( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a, const Real& p )
 {
-   return ExpressionLpNorm( a, p );
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = ( Real ) 0.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = TNL::pow( Containers::Expressions::ExpressionLpNorm( a, p ), p );
+      CommunicatorType::template Allreduce< Real >( &localResult, &result, 1, MPI_SUM, a.getCommunicationGroup() );
+   }
+   return TNL::pow( result, 1.0 / p );
 }
 
 template< typename L1,
@@ -2014,7 +2437,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 product( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
 {
-   return ExpressionProduct( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = ( Real ) 1.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionProduct( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_PROD, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -2023,7 +2452,13 @@ template< typename L1,
 typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 product( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
 {
-   return ExpressionProduct( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = ( Real ) 1.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionProduct( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_PROD, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -2032,7 +2467,44 @@ template< typename L1,
 bool
 logicalOr( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
 {
-   return ExpressionLogicalOr( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   bool result = false;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionLogicalOr( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_LOR, a.getCommunicationGroup() );
+   }
+   return result;
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter >
+bool
+logicalOr( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
+{
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::CommunicatorType;
+   bool result = false;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionLogicalOr( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_LOR, a.getCommunicationGroup() );
+   }
+   return result;
+}
+
+
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation >
+bool
+logicalAnd( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
+{
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = std::numeric_limits< Real >::max();
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionLogicalAnd( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_LAND, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -2041,7 +2513,13 @@ template< typename L1,
 bool
 logicalAnd( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
 {
-   return ExpressionLogicalAnd( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::CommunicatorType;
+   Real result = std::numeric_limits< Real >::max();
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionLogicalAnd( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_LAND, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
@@ -2050,18 +2528,59 @@ template< typename L1,
 typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 binaryOr( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
 {
-   return ExpressionBinaryOr( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   Real result = ( Real ) 0.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionBinaryOr( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_BOR, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter >
 typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
-binaryAnd( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
+binaryOr( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
 {
-   return ExpressionBinaryAnd( a );
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::CommunicatorType;
+   Real result = ( Real ) 0.0;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionBinaryOr( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_BOR, a.getCommunicationGroup() );
+   }
+   return result;
 }
 
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation >
+typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::RealType
+binaryAnd( const Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >& a )
+{
+   using CommunicatorType = typename Containers::Expressions::DistributedBinaryExpressionTemplate< L1, L2, LOperation >::CommunicatorType;
+   bool result = true;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionBinaryAnd( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_BAND, a.getCommunicationGroup() );
+   }
+   return result;
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter >
+typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
+binaryAnd( const Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >& a )
+{
+   using CommunicatorType = typename Containers::Expressions::DistributedUnaryExpressionTemplate< L1, LOperation, Parameter >::CommunicatorType;
+   bool result = true;
+   if( a.getCommunicationGroup() != CommunicatorType::NullGroup ) {
+      const Real localResult = Containers::Expressions::ExpressionBinaryAnd( a );
+      CommunicatorType::Allreduce( &localResult, &result, 1, MPI_BAND, a.getCommunicationGroup() );
+   }
+   return result;
+}
 
 ////
 // Scalar product
diff --git a/src/TNL/Containers/Expressions/ExpressionTemplates.h b/src/TNL/Containers/Expressions/ExpressionTemplates.h
index 659e4c69c1..f80da43aa8 100644
--- a/src/TNL/Containers/Expressions/ExpressionTemplates.h
+++ b/src/TNL/Containers/Expressions/ExpressionTemplates.h
@@ -53,7 +53,7 @@ struct BinaryExpressionTemplate< T1, T2, Operation, VectorVariable, VectorVariab
 
    static_assert( std::is_same< typename T1::DeviceType, typename T2::DeviceType >::value, "Attempt to mix operands allocated on different device types." );
    static_assert( IsStaticType< T1 >::value == IsStaticType< T2 >::value, "Attempt to mix static and non-static operands in binary expression templates." );
-   static constexpr bool isStatic() { return false; }
+   static constexpr bool is() { return false; }
 
    BinaryExpressionTemplate( const T1& a, const T2& b ): op1( a ), op2( b ){}
 
@@ -96,7 +96,7 @@ struct BinaryExpressionTemplate< T1, T2, Operation, VectorVariable, ArithmeticVa
    using IndexType = typename T1::IndexType;
 
    using IsExpressionTemplate = bool;
-   static constexpr bool isStatic() { return false; }
+   static constexpr bool is() { return false; }
 
    BinaryExpressionTemplate( const T1& a, const T2& b ): op1( a ), op2( b ){}
 
@@ -139,7 +139,7 @@ struct BinaryExpressionTemplate< T1, T2, Operation, ArithmeticVariable, VectorVa
    using IndexType = typename T2::IndexType;
 
    using IsExpressionTemplate = bool;
-   static constexpr bool isStatic() { return false; }
+   static constexpr bool is() { return false; }
 
    BinaryExpressionTemplate( const T1& a, const T2& b ): op1( a ), op2( b ){}
 
@@ -186,7 +186,7 @@ struct UnaryExpressionTemplate< T1, Operation, Parameter, VectorVariable >
    using DeviceType = typename T1::DeviceType;
    using IndexType = typename T1::IndexType;
    using IsExpressionTemplate = bool;
-   static constexpr bool isStatic() { return false; }
+   static constexpr bool is() { return false; }
 
    UnaryExpressionTemplate( const T1& a, const Parameter& p )
    : operand( a ), parameter( p ) {}
@@ -233,7 +233,7 @@ struct UnaryExpressionTemplate< T1, Operation, void, VectorVariable >
    using DeviceType = typename T1::DeviceType;
    using IndexType = typename T1::IndexType;
    using IsExpressionTemplate = bool;
-   static constexpr bool isStatic() { return false; }
+   static constexpr bool is() { return false; }
 
    UnaryExpressionTemplate( const T1& a ): operand( a ){}
 
@@ -1096,21 +1096,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator == ( const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator == ( const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator == ( const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
 }
 
 template< typename L1,
@@ -1118,11 +1161,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
 }
 
 template< typename L1,
@@ -1130,11 +1176,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator == ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
              const typename Containers::Expressions::UnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::EQ( a, b );
 }
 
 ////
@@ -1145,21 +1194,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator != ( const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator != ( const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator != ( const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
 }
 
 template< typename L1,
@@ -1167,11 +1259,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
 }
 
 template< typename L1,
@@ -1179,11 +1274,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator != ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::UnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::NE( a, b );
 }
 
 ////
@@ -1194,21 +1292,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& a,
              const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator < ( const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& a,
+             const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator < ( const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+             const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator < ( const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& a,
+             const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
 }
 
 template< typename L1,
@@ -1216,11 +1357,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a,
-              const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
+             const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
 }
 
 template< typename L1,
@@ -1228,11 +1372,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator < ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::UnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::LT( a, b );
 }
 
 ////
@@ -1243,21 +1390,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
-             const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
+              const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& a,
-             const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
+              const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator <= ( const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator <= ( const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator <= ( const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
 }
 
 template< typename L1,
@@ -1265,11 +1455,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
 }
 
 template< typename L1,
@@ -1277,11 +1470,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator <= ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::UnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::LE( a, b );
 }
 
 ////
@@ -1292,21 +1488,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
              const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& a,
              const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator > ( const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& a,
+             const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator > ( const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+             const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator > ( const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& a,
+             const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
 }
 
 template< typename L1,
@@ -1314,11 +1553,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a,
              const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
 }
 
 template< typename L1,
@@ -1326,11 +1568,14 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator > ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
-             const typename Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::UnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::GT( a, b );
 }
 
 ////
@@ -1341,21 +1586,64 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
               const Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
 }
 
 template< typename T1,
           typename T2,
           template< typename, typename > class Operation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   using Right = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator >= ( const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& a,
+              const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& b )
+{
+   using Left = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   using Right = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
+}
+
+template< typename T1,
+          typename T2,
+          template< typename, typename > class Operation >
+__cuda_callable__
+bool
+operator >= ( const typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType& a,
+              const Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >& b )
+{
+   using Left = typename Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >::RealType;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< T1, T2, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
+}
+
+template< typename T1,
+          template< typename > class Operation >
+__cuda_callable__
+bool
+operator >= ( const typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType& a,
+              const Containers::Expressions::UnaryExpressionTemplate< T1, Operation >& b )
+{
+   using Left = typename Containers::Expressions::UnaryExpressionTemplate< T1, Operation >::RealType;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< T1, Operation >;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
 }
 
 template< typename L1,
@@ -1363,11 +1651,14 @@ template< typename L1,
           typename R1,
           typename R2,
           template< typename, typename > class ROperation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a,
               const typename Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >;
+   using Right = Containers::Expressions::BinaryExpressionTemplate< R1, R2, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
 }
 
 template< typename L1,
@@ -1375,17 +1666,20 @@ template< typename L1,
           template< typename, typename > class LOperation,
           typename R1,
           template< typename > class ROperation >
+__cuda_callable__
 bool
 operator >= ( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a,
-              const typename Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >& b )
+             const typename Containers::Expressions::UnaryExpressionTemplate< R1,ROperation >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   using Left = Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >;
+   using Right = Containers::Expressions::UnaryExpressionTemplate< R1, ROperation >;
+   return Containers::Expressions::Comparison< Left, Right >::GE( a, b );
 }
 
+
 ////
 // Unary operations
 
-
 ////
 // Minus
 template< typename L1,
@@ -1937,7 +2231,6 @@ exp( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation >& a
 template< typename L1,
           typename L2,
           template< typename, typename > class LOperation >
-//typename Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >::RealType
 auto
 min( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a ) -> decltype( ExpressionMin( a ) )
 {
@@ -1947,13 +2240,32 @@ min( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter >
-//typename Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 auto
 min( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a ) -> decltype( ExpressionMin( a ) )
 {
    return ExpressionMin( a );
 }
 
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation,
+          typename Index >
+auto
+argMin( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a, Index& arg ) -> decltype( ExpressionArgMin( a, arg ) )
+{
+   return ExpressionArgMin( a, arg );
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter,
+          typename Index >
+auto
+argMin( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a, Index& arg ) -> decltype( ExpressionMin( a, arg ) )
+{
+   return ExpressionArgMin( a );
+}
+
 template< typename L1,
           typename L2,
           template< typename, typename > class LOperation >
@@ -1974,10 +2286,29 @@ max( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Par
    return ExpressionMax( a );
 }
 
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation,
+          typename Index >
+auto
+argMax( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a, Index& arg ) -> decltype( ExpressionArgMax( a, arg ) )
+{
+   return ExpressionArgMax( a, arg );
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter,
+          typename Index >
+auto
+argMax( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a, Index& arg ) -> decltype( ExpressionMax( a, arg ) )
+{
+   return ExpressionArgMax( a );
+}
+
 template< typename L1,
           typename L2,
           template< typename, typename > class LOperation >
-//typename Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >::RealType
 auto
 sum( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a ) -> decltype( ExpressionSum( a ) )
 {
@@ -1987,7 +2318,6 @@ sum( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter >
-//typename Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 auto
 sum( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a ) -> decltype( ExpressionSum( a ) )
 {
@@ -1998,28 +2328,33 @@ template< typename L1,
           typename L2,
           template< typename, typename > class LOperation,
           typename Real >
-//typename Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >::RealType
 auto
 lpNorm( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a, const Real& p ) -> decltype( ExpressionLpNorm( a, p ) )
 {
-   return ExpressionLpNorm( a, p );
+   if( p == 1.0 )
+      return ExpressionLpNorm( a, p );
+   if( p == 2.0 )
+      return TNL::sqrt( ExpressionLpNorm( a, p ) );
+   return TNL::pow( ExpressionLpNorm( a, p ), 1.0 / p );
 }
 
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter,
           typename Real >
-//typename Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 auto
 lpNorm( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a, const Real& p ) -> decltype( ExpressionLpNorm( a, p ) )
 {
-   return ExpressionLpNorm( a, p );
+   if( p == 1.0 )
+      return ExpressionLpNorm( a, p );
+   if( p == 2.0 )
+      return TNL::sqrt( ExpressionLpNorm( a, p ) );
+   return TNL::pow( ExpressionLpNorm( a, p ), 1.0 / p );
 }
 
 template< typename L1,
           typename L2,
           template< typename, typename > class LOperation >
-//typename Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >::RealType
 auto
 product( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a ) -> decltype( ExpressionProduct( a ) )
 {
@@ -2029,7 +2364,6 @@ product( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOpera
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter >
-//typename Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 auto
 product( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a ) -> decltype( ExpressionProduct( a ) )
 {
@@ -2045,6 +2379,24 @@ logicalOr( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOpe
    return ExpressionLogicalOr( a );
 }
 
+template< typename L1,
+          template< typename, typename > class LOperation,
+          typename Parameter >
+bool
+logicalOr( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a )
+{
+   return ExpressionLogicalOr( a );
+}
+
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation >
+bool
+logicalAnd( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a )
+{
+   return ExpressionLogicalAnd( a );
+}
+
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter >
@@ -2057,7 +2409,6 @@ logicalAnd( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperati
 template< typename L1,
           typename L2,
           template< typename, typename > class LOperation >
-//typename Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >::RealType
 auto
 binaryOr( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOperation >& a ) -> decltype( ExpressionBinaryOr( a ) )
 {
@@ -2067,7 +2418,24 @@ binaryOr( const Containers::Expressions::BinaryExpressionTemplate< L1, L2, LOper
 template< typename L1,
           template< typename > class LOperation,
           typename Parameter >
-//typename Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
+auto
+binaryOr( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a ) -> decltype( ExpressionBinaryOr( a ) )
+{
+   return ExpressionBinaryOr( a );
+}
+
+template< typename L1,
+          typename L2,
+          template< typename, typename > class LOperation >
+auto
+binaryAnd( const Containers::Expressions::BiaryExpressionTemplate< L1, L2, LOperation >& a ) -> decltype( ExpressionBinaryAnd( a ) )
+{
+   return ExpressionBinaryAnd( a );
+}
+
+template< typename L1,
+          template< typename > class LOperation,
+          typename Parameter >
 auto
 binaryAnd( const Containers::Expressions::UnaryExpressionTemplate< L1, LOperation, Parameter >& a ) -> decltype( ExpressionBinaryAnd( a ) )
 {
diff --git a/src/TNL/Containers/Expressions/StaticComparison.h b/src/TNL/Containers/Expressions/StaticComparison.h
new file mode 100644
index 0000000000..388e554333
--- /dev/null
+++ b/src/TNL/Containers/Expressions/StaticComparison.h
@@ -0,0 +1,192 @@
+/***************************************************************************
+                          StaticComparison.h  -  description
+                             -------------------
+    begin                : Jul 3, 2019
+    copyright            : (C) 2019 by Tomas Oberhuber
+    email                : tomas.oberhuber@fjfi.cvut.cz
+ ***************************************************************************/
+
+/* See Copyright Notice in tnl/Copyright */
+
+#pragma once
+
+#include <TNL/Assert.h>
+#include <TNL/Containers/Algorithms/Reduction.h>
+
+namespace TNL {
+   namespace Containers {
+      namespace Expressions {
+
+template< typename T1,
+          typename T2,
+          ExpressionVariableType T1Type = ExpressionVariableTypeGetter< T1 >::value,
+          ExpressionVariableType T2Type = ExpressionVariableTypeGetter< T2 >::value >
+struct StaticComparison
+{
+};
+
+/////
+// Static comparison of vector expressions
+template< typename T1,
+          typename T2 >
+struct StaticComparison< T1, T2, VectorVariable, VectorVariable >
+{
+
+   __cuda_callable__
+   static bool EQ( const T1& a, const T2& b )
+   {
+      TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a[ i ] != b[ i ] )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool NE( const T1& a, const T2& b )
+   {
+      return ! EQ( a, b );
+   }
+
+   __cuda_callable__
+   static bool GT( const T1& a, const T2& b )
+   {
+      TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a[ i ] <= b[ i ] )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool LE( const T1& a, const T2& b )
+   {
+      return ! GT( a, b );
+   }
+
+   __cuda_callable__
+   static bool LT( const T1& a, const T2& b )
+   {
+      TNL_ASSERT_EQ( a.getSize(), b.getSize(), "Sizes of expressions to be compared do not fit." );
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a[ i ] >= b[ i ] )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool GE( const T1& a, const T2& b )
+   {
+      return ! LT( a, b );
+   }
+};
+
+/////
+// Static comparison of number and vector expressions
+template< typename T1,
+          typename T2 >
+struct StaticComparison< T1, T2, ArithmeticVariable, VectorVariable >
+{
+
+   __cuda_callable__
+   static bool EQ( const T1& a, const T2& b )
+   {
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a != b[ i ] )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool NE( const T1& a, const T2& b )
+   {
+      return ! EQ( a, b );
+   }
+
+   __cuda_callable__
+   static bool GT( const T1& a, const T2& b )
+   {
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a <= b[ i ] )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool LE( const T1& a, const T2& b )
+   {
+      return ! GT( a, b );
+   }
+
+   __cuda_callable__
+   static bool LT( const T1& a, const T2& b )
+   {
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a >= b[ i ] )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool GE( const T1& a, const T2& b )
+   {
+      return ! LT( a, b );
+   }
+};
+
+/////
+// Static comparison of vector expressions and number
+template< typename T1,
+          typename T2 >
+struct StaticComparison< T1, T2, VectorVariable, ArithmeticVariable >
+{
+
+   __cuda_callable__
+   static bool EQ( const T1& a, const T2& b )
+   {
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a[ i ] != b )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool NE( const T1& a, const T2& b )
+   {
+      return ! EQ( a, b );
+   }
+
+   __cuda_callable__
+   static bool GT( const T1& a, const T2& b )
+   {
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a[ i ] <= b )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool LE( const T1& a, const T2& b )
+   {
+      return ! GT( a, b );
+   }
+
+   __cuda_callable__
+   static bool LT( const T1& a, const T2& b )
+   {
+      for( int i = 0; i < a.getSize(); i++ )
+         if( a[ i ] >= b )
+            return false;
+      return true;
+   }
+
+   __cuda_callable__
+   static bool GE( const T1& a, const T2& b )
+   {
+      return ! LT( a, b );
+   }
+};
+
+      } //namespace Expressions
+   } // namespace Containers
+} // namespace TNL
diff --git a/src/TNL/Containers/Expressions/StaticExpressionTemplates.h b/src/TNL/Containers/Expressions/StaticExpressionTemplates.h
index d9e5c655c4..4427bc7a58 100644
--- a/src/TNL/Containers/Expressions/StaticExpressionTemplates.h
+++ b/src/TNL/Containers/Expressions/StaticExpressionTemplates.h
@@ -13,7 +13,7 @@
 #include <iostream>
 #include <TNL/Containers/Expressions/ExpressionTemplatesOperations.h>
 #include <TNL/Containers/Expressions/ExpressionVariableType.h>
-#include <TNL/Containers/Expressions/Comparison.h>
+#include <TNL/Containers/Expressions/StaticComparison.h>
 #include <TNL/Containers/Expressions/IsStatic.h>
 #include <TNL/Containers/Expressions/VerticalOperations.h>
 
@@ -2462,7 +2462,11 @@ __cuda_callable__
 typename Containers::Expressions::StaticBinaryExpressionTemplate< L1, L2, LOperation >::RealType
 lpNorm( const Containers::Expressions::StaticBinaryExpressionTemplate< L1, L2, LOperation >& a, const Real& p )
 {
-   return StaticExpressionLpNorm( a, p );
+   if( p == 1.0 )
+      return StaticExpressionLpNorm( a, p );
+   if( p == 2.0 )
+      return TNL::sqrt( StaticExpressionLpNorm( a, p ) );
+   return TNL::pow( StaticExpressionLpNorm( a, p ), 1.0 / p );
 }
 
 template< typename L1,
@@ -2473,7 +2477,11 @@ __cuda_callable__
 typename Containers::Expressions::StaticUnaryExpressionTemplate< L1, LOperation, Parameter >::RealType
 lpNorm( const Containers::Expressions::StaticUnaryExpressionTemplate< L1, LOperation, Parameter >& a, const Real& p )
 {
-   return StaticExpressionLpNorm( a, p );
+   if( p == 1.0 )
+      return StaticExpressionLpNorm( a, p );
+   if( p == 2.0 )
+      return TNL::sqrt( StaticExpressionLpNorm( a, p ) );
+   return TNL::pow( StaticExpressionLpNorm( a, p ), 1.0 / p );
 }
 
 template< typename L1,
diff --git a/src/TNL/Containers/Expressions/StaticVerticalOperations.h b/src/TNL/Containers/Expressions/StaticVerticalOperations.h
new file mode 100644
index 0000000000..4f33a1a7cd
--- /dev/null
+++ b/src/TNL/Containers/Expressions/StaticVerticalOperations.h
@@ -0,0 +1,165 @@
+/***************************************************************************
+                          StaticVerticalOperations.h  -  description
+                             -------------------
+    begin                : Jul 3, 2019
+    copyright            : (C) 2019 by Tomas Oberhuber
+    email                : tomas.oberhuber@fjfi.cvut.cz
+ ***************************************************************************/
+
+/* See Copyright Notice in tnl/Copyright */
+
+#pragma once
+
+#include <TNL/Assert.h>
+#include <TNL/Containers/Algorithms/Reduction.h>
+
+////
+// By vertical operations we mean those applied across vector elements or
+// vector expression elements. It means for example minim/maximum of all
+// vector elements etc.
+namespace TNL {
+   namespace Containers {
+      namespace Expressions {
+
+template< typename Expression >
+__cuda_callable__
+auto StaticExpressionMin( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux = TNL::min( aux, expression[ i ] );
+   return aux;
+}
+
+template< typename Expression, typename Real >
+__cuda_callable__
+auto StaticExpressionArgMin( const Expression& expression, int& arg ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto value = expression[ 0 ];
+   arg = 0;
+   for( int i = 1; i < expression.getSize(); i++ )
+   {
+      if( expression[ i ] < value )
+      {
+         value = expression[ i ];
+         arg = i;
+      }
+   }
+   return value;
+}
+
+template< typename Expression >
+__cuda_callable__
+auto StaticExpressionMax( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux = TNL::max( aux, expression[ i ] );
+   return aux;
+}
+
+template< typename Expression, typename Real >
+__cuda_callable__
+auto StaticExpressionArgMax( const Expression& expression, int& arg ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto value = expression[ 0 ];
+   arg = 0;
+   for( int i = 1; i < expression.getSize(); i++ )
+   {
+      if( expression[ i ] > value )
+      {
+         value = expression[ i ];
+         arg = i;
+      }
+   }
+   return value;
+}
+
+template< typename Expression >
+__cuda_callable__
+auto StaticExpressionSum( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux += expression[ i ];
+   return aux;
+}
+
+template< typename Expression, typename Real >
+__cuda_callable__
+auto StaticExpressionLpNorm( const Expression& expression, const Real& p ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   if( p == ( Real ) 1.0 )
+   {
+      auto aux = TNL::abs( expression[ 0 ] );
+      for( int i = 1; i < expression.getSize(); i++ )
+         aux += TNL::abs( expression[ i ] );
+      return aux;
+   }
+   if( p == ( Real ) 2.0 )
+   {
+      auto aux = expression[ 0 ] * expression[ 0 ];
+      for( int i = 1; i < expression.getSize(); i++ )
+         aux += expression[ i ] * expression[ i ];
+      return aux;
+   }
+   auto aux = TNL::pow( expression[ 0 ], p );
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux += TNL::pow( expression[ i ], p );
+   return aux;
+}
+
+
+template< typename Expression >
+__cuda_callable__
+auto StaticExpressionProduct( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux *= expression[ i ];
+   return aux;
+}
+
+template< typename Expression >
+__cuda_callable__
+bool StaticExpressionLogicalAnd( const Expression& expression )
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux = aux && expression[ i ];
+   return aux;
+}
+
+template< typename Expression >
+__cuda_callable__
+bool StaticExpressionLogicalOr( const Expression& expression )
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux = aux || expression[ i ];
+   return aux;
+}
+
+template< typename Expression >
+__cuda_callable__
+auto StaticExpressionBinaryAnd( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux = aux & expression[ i ];
+   return aux;
+}
+
+template< typename Expression >
+__cuda_callable__
+auto StaticExpressionBinaryOr( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
+{
+   auto aux = expression[ 0 ];
+   for( int i = 1; i < expression.getSize(); i++ )
+      aux = aux | expression[ i ];
+   return aux;
+}
+
+      } //namespace Expressions
+   } // namespace Containers
+} // namespace TNL
diff --git a/src/TNL/Containers/Expressions/VerticalOperations.h b/src/TNL/Containers/Expressions/VerticalOperations.h
index e57b2c1cf0..408b9925fd 100644
--- a/src/TNL/Containers/Expressions/VerticalOperations.h
+++ b/src/TNL/Containers/Expressions/VerticalOperations.h
@@ -21,147 +21,8 @@ namespace TNL {
    namespace Containers {
       namespace Expressions {
 
-template< typename Expression >
-__cuda_callable__
-auto StaticExpressionMin( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux = TNL::min( aux, expression[ i ] );
-   return aux;
-}
-
-template< typename Expression, typename Real >
-__cuda_callable__
-auto StaticExpressionArgMin( const Expression& expression, int& arg ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto value = expression[ 0 ];
-   arg = 0;
-   for( int i = 1; i < expression.getSize(); i++ )
-   {
-      if( expression[ i ] < value )
-      {
-         value = expression[ i ];
-         arg = i;
-      }
-   }
-   return value;
-}
-
-template< typename Expression >
-__cuda_callable__
-auto StaticExpressionMax( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux = TNL::max( aux, expression[ i ] );
-   return aux;
-}
-
-template< typename Expression, typename Real >
-__cuda_callable__
-auto StaticExpressionArgMax( const Expression& expression, int& arg ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto value = expression[ 0 ];
-   arg = 0;
-   for( int i = 1; i < expression.getSize(); i++ )
-   {
-      if( expression[ i ] > value )
-      {
-         value = expression[ i ];
-         arg = i;
-      }
-   }
-   return value;
-}
-
-template< typename Expression >
-__cuda_callable__
-auto StaticExpressionSum( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux += expression[ i ];
-   return aux;
-}
-
-template< typename Expression, typename Real >
-__cuda_callable__
-auto StaticExpressionLpNorm( const Expression& expression, const Real& p ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   if( p == ( Real ) 1.0 )
-   {
-      auto aux = TNL::abs( expression[ 0 ] );
-      for( int i = 1; i < expression.getSize(); i++ )
-         aux += TNL::abs( expression[ i ] );
-      return aux;
-   }
-   if( p == ( Real ) 2.0 )
-   {
-      auto aux = expression[ 0 ] * expression[ 0 ];
-      for( int i = 1; i < expression.getSize(); i++ )
-         aux += expression[ i ] * expression[ i ];
-      return TNL::sqrt( aux );
-   }
-   auto aux = TNL::pow( expression[ 0 ], p );
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux += TNL::pow( expression[ i ], p );
-   return TNL::pow( aux, 1.0 / p );
-}
-
-
-template< typename Expression >
-__cuda_callable__
-auto StaticExpressionProduct( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux *= expression[ i ];
-   return aux;
-}
-
-template< typename Expression >
-__cuda_callable__
-bool StaticExpressionLogicalAnd( const Expression& expression )
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux = aux && expression[ i ];
-   return aux;
-}
-
-template< typename Expression >
-__cuda_callable__
-bool StaticExpressionLogicalOr( const Expression& expression )
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux = aux || expression[ i ];
-   return aux;
-}
-
-template< typename Expression >
-__cuda_callable__
-auto StaticExpressionBinaryAnd( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux = aux & expression[ i ];
-   return aux;
-}
-
-template< typename Expression >
-__cuda_callable__
-auto StaticExpressionBinaryOr( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
-{
-   auto aux = expression[ 0 ];
-   for( int i = 1; i < expression.getSize(); i++ )
-      aux = aux | expression[ i ];
-   return aux;
-}
-
 ////
-// Non-static operations
+// Vertical operations
 template< typename Expression >
 auto ExpressionMin( const Expression& expression ) -> typename std::remove_reference< decltype( expression[ 0 ] ) >::type
 {
@@ -282,12 +143,12 @@ auto ExpressionLpNorm( const Expression& expression, const Real& p ) ->
       auto fetch = [=] __cuda_callable__ ( IndexType i ) { return  expression[ i ] * expression[ i ]; };
       auto reduction = [=] __cuda_callable__ ( ResultType& a, const ResultType& b ) { a += b; };
       auto volatileReduction = [=] __cuda_callable__ ( volatile ResultType& a, volatile ResultType& b ) { a += b; };
-      return TNL::sqrt( Algorithms::Reduction< typename Expression::DeviceType >::reduce( expression.getSize(), reduction, volatileReduction, fetch, ( ResultType ) 0.0 ) );
+      return Algorithms::Reduction< typename Expression::DeviceType >::reduce( expression.getSize(), reduction, volatileReduction, fetch, ( ResultType ) 0.0 );
    }
    auto fetch = [=] __cuda_callable__ ( IndexType i ) { return  TNL::pow( expression[ i ], p ); };
    auto reduction = [=] __cuda_callable__ ( ResultType& a, const ResultType& b ) { a += b; };
    auto volatileReduction = [=] __cuda_callable__ ( volatile ResultType& a, volatile ResultType& b ) { a += b; };
-   return TNL::pow( Algorithms::Reduction< typename Expression::DeviceType >::reduce( expression.getSize(), reduction, volatileReduction, fetch, ( ResultType ) 0.0 ), ( Real ) 1.0 / p );
+   return Algorithms::Reduction< typename Expression::DeviceType >::reduce( expression.getSize(), reduction, volatileReduction, fetch, ( ResultType ) 0.0 );
 }
 
 template< typename Expression >
diff --git a/src/TNL/Containers/VectorExpressions.h b/src/TNL/Containers/VectorExpressions.h
index aa11d3b5ba..f6fdd51d04 100644
--- a/src/TNL/Containers/VectorExpressions.h
+++ b/src/TNL/Containers/VectorExpressions.h
@@ -203,13 +203,15 @@ max( const Containers::Vector< Real1, Device, Index >& a, const Containers::Vect
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator==( const Containers::Vector< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonEQ( a.getView(), b );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView, ET >::EQ( a.getView(), b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator==( const ET& a, const Containers::Vector< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b.getView() );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ET, ConstView >::EQ( a, b.getView() );
 }
 
 template< typename Real1, typename Real2, typename Device1, typename Device2, typename Index >
@@ -230,13 +232,15 @@ bool operator==( const Containers::Vector< Real1, Device1, Index >& a, const Con
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator!=( const Containers::Vector< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonNE( a.getView(), b );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView, ET >::NE( a.getView(), b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator!=( const ET& a, const Containers::Vector< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b.getView() );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ET, ConstView >::NE( a, b.getView() );
 }
 
 template< typename Real1, typename Real2, typename Device1, typename Device2, typename Index >
@@ -257,19 +261,23 @@ bool operator!=( const Containers::Vector< Real1, Device1, Index >& a, const Con
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator<( const Containers::Vector< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonLT( a.getView(), b );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView, ET >::LT( a.getView(), b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator<( const ET& a, const Containers::Vector< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b.getView() );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ET, ConstView >::LT( a, b.getView() );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator<( const Containers::Vector< Real1, Device, Index >& a, const Containers::Vector< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLT( a.getView(), b.getView() );
+   using ConstView1 = typename Containers::VectorView< Real1, Device, Index >::ConstViewType;
+   using ConstView2 = typename Containers::VectorView< Real2, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView1, ConstView2 >::LT( a.getView(), b.getView() );
 }
 
 ////
@@ -277,19 +285,23 @@ bool operator<( const Containers::Vector< Real1, Device, Index >& a, const Conta
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator<=( const Containers::Vector< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonLE( a.getView(), b );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView, ET >::LE( a.getView(), b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator<=( const ET& a, const Containers::Vector< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b.getView() );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ET, ConstView >::LE( a, b.getView() );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator<=( const Containers::Vector< Real1, Device, Index >& a, const Containers::Vector< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLE( a.getView(), b.getView() );
+   using ConstView1 = typename Containers::VectorView< Real1, Device, Index >::ConstViewType;
+   using ConstView2 = typename Containers::VectorView< Real2, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView1, ConstView2 >::LE( a.getView(), b.getView() );
 }
 
 ////
@@ -297,19 +309,23 @@ bool operator<=( const Containers::Vector< Real1, Device, Index >& a, const Cont
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator>( const Containers::Vector< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonGT( a.getView(), b );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView, ET >::GT( a.getView(), b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator>( const ET& a, const Containers::Vector< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b.getView() );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ET, ConstView >::GT( a, b.getView() );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator>( const Containers::Vector< Real1, Device, Index >& a, const Containers::Vector< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGT( a.getView(), b.getView() );
+   using ConstView1 = typename Containers::VectorView< Real1, Device, Index >::ConstViewType;
+   using ConstView2 = typename Containers::VectorView< Real2, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView1, ConstView2 >::GT( a.getView(), b.getView() );
 }
 
 ////
@@ -317,19 +333,23 @@ bool operator>( const Containers::Vector< Real1, Device, Index >& a, const Conta
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator>=( const Containers::Vector< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonGE( a.getView(), b );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView, ET >::GE( a.getView(), b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator>=( const ET& a, const Containers::Vector< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b.getView() );
+   using ConstView = typename Containers::VectorView< Real, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ET, ConstView >::GE( a, b.getView() );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator>=( const Containers::Vector< Real1, Device, Index >& a, const Containers::Vector< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGE( a.getView(), b.getView() );
+   using ConstView1 = typename Containers::VectorView< Real1, Device, Index >::ConstViewType;
+   using ConstView2 = typename Containers::VectorView< Real2, Device, Index >::ConstViewType;
+   return Containers::Expressions::Comparison< ConstView1, ConstView2 >::GE( a.getView(), b.getView() );
 }
 
 ////
@@ -586,7 +606,11 @@ template< typename Real,
 auto
 lpNorm( const Containers::Vector< Real, Device, Index >& a, const Real2& p ) -> decltype( Containers::Expressions::ExpressionLpNorm( a.getView(), p ) )
 {
-   return Containers::Expressions::ExpressionLpNorm( a.getView(), p );
+   if( p == 1.0 )
+      return Containers::Expressions::ExpressionLpNorm( a.getView(), p );
+   if( p == 2.0 )
+      return TNL::sqrt( Containers::Expressions::ExpressionLpNorm( a.getView(), p ) );
+   return TNL::pow( Containers::Expressions::ExpressionLpNorm( a.getView(), p ), 1.0 / p );
 }
 
 template< typename Real,
diff --git a/src/TNL/Containers/VectorViewExpressions.h b/src/TNL/Containers/VectorViewExpressions.h
index 44b7c4cca2..4584fd4f8b 100644
--- a/src/TNL/Containers/VectorViewExpressions.h
+++ b/src/TNL/Containers/VectorViewExpressions.h
@@ -165,13 +165,13 @@ max( const Containers::VectorView< Real1, Device, Index >& a, const Containers::
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator==( const Containers::VectorView< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real, Device, Index >, ET >::EQ( a, b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator==( const ET& a, const Containers::VectorView< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonEQ( a, b );
+   return Containers::Expressions::Comparison< ET, Containers::VectorView< Real, Device, Index > >::EQ( a, b );
 }
 
 template< typename Real1, typename Real2, typename Device1, typename Device2, typename Index >
@@ -192,13 +192,13 @@ bool operator==( const Containers::VectorView< Real1, Device1, Index >& a, const
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator!=( const Containers::VectorView< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real, Device, Index >, ET >::NE( a, b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator!=( const ET& a, const Containers::VectorView< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonNE( a, b );
+   return Containers::Expressions::Comparison< ET, Containers::VectorView< Real, Device, Index > >::NE( a, b );
 }
 
 template< typename Real1, typename Real2, typename Device1, typename Device2, typename Index >
@@ -219,19 +219,19 @@ bool operator!=( const Containers::VectorView< Real1, Device1, Index >& a, const
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator<( const Containers::VectorView< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real, Device, Index >, ET >::LT( a, b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator<( const ET& a, const Containers::VectorView< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   return Containers::Expressions::Comparison< ET, Containers::VectorView< Real, Device, Index > >::LT( a, b );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator<( const Containers::VectorView< Real1, Device, Index >& a, const Containers::VectorView< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLT( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real1, Device, Index >, Containers::VectorView< Real2, Device, Index > >::LT( a, b );
 }
 
 ////
@@ -239,19 +239,19 @@ bool operator<( const Containers::VectorView< Real1, Device, Index >& a, const C
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator<=( const Containers::VectorView< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real, Device, Index >, ET >::LE( a, b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator<=( const ET& a, const Containers::VectorView< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   return Containers::Expressions::Comparison< ET, Containers::VectorView< Real, Device, Index > >::LE( a, b );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator<=( const Containers::VectorView< Real1, Device, Index >& a, const Containers::VectorView< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonLE( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real1, Device, Index >, Containers::VectorView< Real2, Device, Index > >::LE( a, b );
 }
 
 ////
@@ -259,19 +259,19 @@ bool operator<=( const Containers::VectorView< Real1, Device, Index >& a, const
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator>( const Containers::VectorView< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real, Device, Index >, ET >::GT( a, b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator>( const ET& a, const Containers::VectorView< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   return Containers::Expressions::Comparison< ET, Containers::VectorView< Real, Device, Index > >::GT( a, b );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator>( const Containers::VectorView< Real1, Device, Index >& a, const Containers::VectorView< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGT( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real1, Device, Index >, Containers::VectorView< Real2, Device, Index > >::GT( a, b );
 }
 
 ////
@@ -279,19 +279,19 @@ bool operator>( const Containers::VectorView< Real1, Device, Index >& a, const C
 template< typename Real, typename Device, typename Index, typename ET >
 bool operator>=( const Containers::VectorView< Real, Device, Index >& a, const ET& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real, Device, Index >, ET >::GE( a, b );
 }
 
 template< typename ET, typename Real, typename Device, typename Index >
 bool operator>=( const ET& a, const Containers::VectorView< Real, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   return Containers::Expressions::Comparison< ET, Containers::VectorView< Real, Device, Index > >::GE( a, b );
 }
 
 template< typename Real1, typename Real2, typename Device, typename Index >
 bool operator>=( const Containers::VectorView< Real1, Device, Index >& a, const Containers::VectorView< Real2, Device, Index >& b )
 {
-   return Containers::Expressions::ComparisonGE( a, b );
+   return Containers::Expressions::Comparison< Containers::VectorView< Real1, Device, Index >, Containers::VectorView< Real2, Device, Index > >::GE( a, b );
 }
 
 ////
@@ -530,7 +530,11 @@ auto
 lpNorm( const Containers::VectorView< Real, Device, Index >& a, const Real2& p )
 -> decltype( Containers::Expressions::ExpressionLpNorm( a, p ) )
 {
-   return Containers::Expressions::ExpressionLpNorm( a, p );
+   if( p == 1.0 )
+      return Containers::Expressions::ExpressionLpNorm( a, p );
+   if( p == 2.0 )
+      return TNL::sqrt( Containers::Expressions::ExpressionLpNorm( a, p ) );
+   return TNL::pow( Containers::Expressions::ExpressionLpNorm( a, p ), 1.0 / p );
 }
 
 template< typename Real,
-- 
GitLab