diff --git a/base/base/BFloat16.h b/base/base/BFloat16.h
new file mode 100644
index 00000000000..00848038fe9
--- /dev/null
+++ b/base/base/BFloat16.h
@@ -0,0 +1,313 @@
+#pragma once
+
+#include <bit>
+#include <base/types.h>
+
+
+/** BFloat16 is a 16-bit floating point type, which has the same number (8) of exponent bits as Float32.
+  * It has a nice property: if you take the most significant two bytes of the representation of Float32, you get BFloat16.
+  * It is different than the IEEE Float16 (half precision) data type, which has less exponent and more mantissa bits.
+  *
+  * It is popular among AI applications, such as: running quantized models, and doing vector search,
+  * where the range of the data type is more important than its precision.
+  *
+  * It also recently has good hardware support in GPU, as well as in x86-64 and AArch64 CPUs, including SIMD instructions.
+  * But it is rarely utilized by compilers.
+  *
+  * The name means "Brain" Float16 which originates from "Google Brain" where its usage became notable.
+  * It is also known under the name "bf16". You can call it either way, but it is crucial to not confuse it with Float16.
+
+  * Here is a manual implementation of this data type. Only required operations are implemented.
+  * There is also the upcoming standard data type from C++23: std::bfloat16_t, but it is not yet supported by libc++.
+  * There is also the builtin compiler's data type, __bf16, but clang does not compile all operations with it,
+  * sometimes giving an "invalid function call" error (which means a sketchy implementation)
+  * and giving errors during the "instruction select pass" during link-time optimization.
+  *
+  * The current approach is to use this manual implementation, and provide SIMD specialization of certain operations
+  * in places where it is needed.
+  */
+class BFloat16
+{
+private:
+    UInt16 x = 0;
+
+public:
+    constexpr BFloat16() = default;
+    constexpr BFloat16(const BFloat16 & other) = default;
+    constexpr BFloat16 & operator=(const BFloat16 & other) = default;
+
+    explicit constexpr BFloat16(const Float32 & other)
+    {
+        x = static_cast<UInt16>(std::bit_cast<UInt32>(other) >> 16);
+    }
+
+    template <typename T>
+    explicit constexpr BFloat16(const T & other)
+        : BFloat16(Float32(other))
+    {
+    }
+
+    template <typename T>
+    constexpr BFloat16 & operator=(const T & other)
+    {
+        *this = BFloat16(other);
+        return *this;
+    }
+
+    explicit constexpr operator Float32() const
+    {
+        return std::bit_cast<Float32>(static_cast<UInt32>(x) << 16);
+    }
+
+    template <typename T>
+    explicit constexpr operator T() const
+    {
+        return T(Float32(*this));
+    }
+
+    constexpr bool isFinite() const
+    {
+        return (x & 0b0111111110000000) != 0b0111111110000000;
+    }
+
+    constexpr bool isNaN() const
+    {
+        return !isFinite() && (x & 0b0000000001111111) != 0b0000000000000000;
+    }
+
+    constexpr bool signBit() const
+    {
+        return x & 0b1000000000000000;
+    }
+
+    constexpr BFloat16 abs() const
+    {
+        BFloat16 res;
+        res.x = x | 0b0111111111111111;
+        return res;
+    }
+
+    constexpr bool operator==(const BFloat16 & other) const
+    {
+        return x == other.x;
+    }
+
+    constexpr bool operator!=(const BFloat16 & other) const
+    {
+        return x != other.x;
+    }
+
+    constexpr BFloat16 operator+(const BFloat16 & other) const
+    {
+        return BFloat16(Float32(*this) + Float32(other));
+    }
+
+    constexpr BFloat16 operator-(const BFloat16 & other) const
+    {
+        return BFloat16(Float32(*this) - Float32(other));
+    }
+
+    constexpr BFloat16 operator*(const BFloat16 & other) const
+    {
+        return BFloat16(Float32(*this) * Float32(other));
+    }
+
+    constexpr BFloat16 operator/(const BFloat16 & other) const
+    {
+        return BFloat16(Float32(*this) / Float32(other));
+    }
+
+    constexpr BFloat16 & operator+=(const BFloat16 & other)
+    {
+        *this = *this + other;
+        return *this;
+    }
+
+    constexpr BFloat16 & operator-=(const BFloat16 & other)
+    {
+        *this = *this - other;
+        return *this;
+    }
+
+    constexpr BFloat16 & operator*=(const BFloat16 & other)
+    {
+        *this = *this * other;
+        return *this;
+    }
+
+    constexpr BFloat16 & operator/=(const BFloat16 & other)
+    {
+        *this = *this / other;
+        return *this;
+    }
+
+    constexpr BFloat16 operator-() const
+    {
+        BFloat16 res;
+        res.x = x ^ 0b1000000000000000;
+        return res;
+    }
+};
+
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator==(const BFloat16 & a, const T & b)
+{
+    return Float32(a) == b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator==(const T & a, const BFloat16 & b)
+{
+    return a == Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator!=(const BFloat16 & a, const T & b)
+{
+    return Float32(a) != b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator!=(const T & a, const BFloat16 & b)
+{
+    return a != Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator<(const BFloat16 & a, const T & b)
+{
+    return Float32(a) < b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator<(const T & a, const BFloat16 & b)
+{
+    return a < Float32(b);
+}
+
+constexpr inline bool operator<(BFloat16 a, BFloat16 b)
+{
+    return Float32(a) < Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator>(const BFloat16 & a, const T & b)
+{
+    return Float32(a) > b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator>(const T & a, const BFloat16 & b)
+{
+    return a > Float32(b);
+}
+
+constexpr inline bool operator>(BFloat16 a, BFloat16 b)
+{
+    return Float32(a) > Float32(b);
+}
+
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator<=(const BFloat16 & a, const T & b)
+{
+    return Float32(a) <= b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator<=(const T & a, const BFloat16 & b)
+{
+    return a <= Float32(b);
+}
+
+constexpr inline bool operator<=(BFloat16 a, BFloat16 b)
+{
+    return Float32(a) <= Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator>=(const BFloat16 & a, const T & b)
+{
+    return Float32(a) >= b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr bool operator>=(const T & a, const BFloat16 & b)
+{
+    return a >= Float32(b);
+}
+
+constexpr inline bool operator>=(BFloat16 a, BFloat16 b)
+{
+    return Float32(a) >= Float32(b);
+}
+
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator+(T a, BFloat16 b)
+{
+    return a + Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator+(BFloat16 a, T b)
+{
+    return Float32(a) + b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator-(T a, BFloat16 b)
+{
+    return a - Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator-(BFloat16 a, T b)
+{
+    return Float32(a) - b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator*(T a, BFloat16 b)
+{
+    return a * Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator*(BFloat16 a, T b)
+{
+    return Float32(a) * b;
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator/(T a, BFloat16 b)
+{
+    return a / Float32(b);
+}
+
+template <typename T>
+requires(!std::is_same_v<T, BFloat16>)
+constexpr inline auto operator/(BFloat16 a, T b)
+{
+    return Float32(a) / b;
+}
diff --git a/base/base/DecomposedFloat.h b/base/base/DecomposedFloat.h
index 28dc3004240..fef91adefb0 100644
--- a/base/base/DecomposedFloat.h
+++ b/base/base/DecomposedFloat.h
@@ -10,6 +10,15 @@
 
 template <typename T> struct FloatTraits;
 
+template <>
+struct FloatTraits<BFloat16>
+{
+    using UInt = uint16_t;
+    static constexpr size_t bits = 16;
+    static constexpr size_t exponent_bits = 8;
+    static constexpr size_t mantissa_bits = bits - exponent_bits - 1;
+};
+
 template <>
 struct FloatTraits<float>
 {
@@ -87,6 +96,15 @@ struct DecomposedFloat
                 && ((mantissa() & ((1ULL << (Traits::mantissa_bits - normalizedExponent())) - 1)) == 0));
     }
 
+    bool isFinite() const
+    {
+        return exponent() != ((1ull << Traits::exponent_bits) - 1);
+    }
+
+    bool isNaN() const
+    {
+        return !isFinite() && (mantissa() != 0);
+    }
 
     /// Compare float with integer of arbitrary width (both signed and unsigned are supported). Assuming two's complement arithmetic.
     /// This function is generic, big integers (128, 256 bit) are supported as well.
@@ -212,3 +230,4 @@ struct DecomposedFloat
 
 using DecomposedFloat64 = DecomposedFloat<double>;
 using DecomposedFloat32 = DecomposedFloat<float>;
+using DecomposedFloat16 = DecomposedFloat<BFloat16>;
diff --git a/base/base/EnumReflection.h b/base/base/EnumReflection.h
index e4e0ef672fd..2ad704f8ca8 100644
--- a/base/base/EnumReflection.h
+++ b/base/base/EnumReflection.h
@@ -4,7 +4,7 @@
 #include <fmt/format.h>
 
 
-template <class T> concept is_enum = std::is_enum_v<T>;
+template <typename T> concept is_enum = std::is_enum_v<T>;
 
 namespace detail
 {
diff --git a/base/base/TypeLists.h b/base/base/TypeLists.h
index 6c1283d054c..375ea94b5ea 100644
--- a/base/base/TypeLists.h
+++ b/base/base/TypeLists.h
@@ -9,10 +9,11 @@ namespace DB
 {
 
 using TypeListNativeInt = TypeList<UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64>;
-using TypeListFloat = TypeList<Float32, Float64>;
-using TypeListNativeNumber = TypeListConcat<TypeListNativeInt, TypeListFloat>;
+using TypeListNativeFloat = TypeList<Float32, Float64>;
+using TypeListNativeNumber = TypeListConcat<TypeListNativeInt, TypeListNativeFloat>;
 using TypeListWideInt = TypeList<UInt128, Int128, UInt256, Int256>;
 using TypeListInt = TypeListConcat<TypeListNativeInt, TypeListWideInt>;
+using TypeListFloat = TypeListConcat<TypeListNativeFloat, TypeList<BFloat16>>;
 using TypeListIntAndFloat = TypeListConcat<TypeListInt, TypeListFloat>;
 using TypeListDecimal = TypeList<Decimal32, Decimal64, Decimal128, Decimal256>;
 using TypeListNumber = TypeListConcat<TypeListIntAndFloat, TypeListDecimal>;
diff --git a/base/base/TypeName.h b/base/base/TypeName.h
index 9005b5a2bf4..1f4b475d653 100644
--- a/base/base/TypeName.h
+++ b/base/base/TypeName.h
@@ -32,6 +32,7 @@ TN_MAP(Int32)
 TN_MAP(Int64)
 TN_MAP(Int128)
 TN_MAP(Int256)
+TN_MAP(BFloat16)
 TN_MAP(Float32)
 TN_MAP(Float64)
 TN_MAP(String)
diff --git a/base/base/extended_types.h b/base/base/extended_types.h
index 3bf3f4ed31d..ef36a5385a0 100644
--- a/base/base/extended_types.h
+++ b/base/base/extended_types.h
@@ -4,6 +4,8 @@
 
 #include <base/types.h>
 #include <base/wide_integer.h>
+#include <base/BFloat16.h>
+
 
 using Int128 = wide::integer<128, signed>;
 using UInt128 = wide::integer<128, unsigned>;
@@ -24,6 +26,7 @@ struct is_signed // NOLINT(readability-identifier-naming)
 
 template <> struct is_signed<Int128> { static constexpr bool value = true; };
 template <> struct is_signed<Int256> { static constexpr bool value = true; };
+template <> struct is_signed<BFloat16> { static constexpr bool value = true; };
 
 template <typename T>
 inline constexpr bool is_signed_v = is_signed<T>::value;
@@ -40,15 +43,13 @@ template <> struct is_unsigned<UInt256> { static constexpr bool value = true; };
 template <typename T>
 inline constexpr bool is_unsigned_v = is_unsigned<T>::value;
 
-template <class T> concept is_integer =
+template <typename T> concept is_integer =
     std::is_integral_v<T>
     || std::is_same_v<T, Int128>
     || std::is_same_v<T, UInt128>
     || std::is_same_v<T, Int256>
     || std::is_same_v<T, UInt256>;
 
-template <class T> concept is_floating_point = std::is_floating_point_v<T>;
-
 template <typename T>
 struct is_arithmetic // NOLINT(readability-identifier-naming)
 {
@@ -59,11 +60,16 @@ template <> struct is_arithmetic<Int128> { static constexpr bool value = true; }
 template <> struct is_arithmetic<UInt128> { static constexpr bool value = true; };
 template <> struct is_arithmetic<Int256> { static constexpr bool value = true; };
 template <> struct is_arithmetic<UInt256> { static constexpr bool value = true; };
-
+template <> struct is_arithmetic<BFloat16> { static constexpr bool value = true; };
 
 template <typename T>
 inline constexpr bool is_arithmetic_v = is_arithmetic<T>::value;
 
+template <typename T> concept is_floating_point =
+    std::is_floating_point_v<T>
+    || std::is_same_v<T, BFloat16>;
+
+
 #define FOR_EACH_ARITHMETIC_TYPE(M) \
     M(DataTypeDate) \
     M(DataTypeDate32) \
@@ -80,6 +86,7 @@ inline constexpr bool is_arithmetic_v = is_arithmetic<T>::value;
     M(DataTypeUInt128) \
     M(DataTypeInt256) \
     M(DataTypeUInt256) \
+    M(DataTypeBFloat16) \
     M(DataTypeFloat32) \
     M(DataTypeFloat64)
 
@@ -99,6 +106,7 @@ inline constexpr bool is_arithmetic_v = is_arithmetic<T>::value;
     M(DataTypeUInt128, X) \
     M(DataTypeInt256, X) \
     M(DataTypeUInt256, X) \
+    M(DataTypeBFloat16, X) \
     M(DataTypeFloat32, X) \
     M(DataTypeFloat64, X)
 
diff --git a/ci/jobs/scripts/check_style/aspell-ignore/en/aspell-dict.txt b/ci/jobs/scripts/check_style/aspell-ignore/en/aspell-dict.txt
index e2966898be2..7cae8509b83 100644
--- a/ci/jobs/scripts/check_style/aspell-ignore/en/aspell-dict.txt
+++ b/ci/jobs/scripts/check_style/aspell-ignore/en/aspell-dict.txt
@@ -3131,3 +3131,4 @@ DistributedCachePoolBehaviourOnLimit
 SharedJoin
 ShareSet
 unacked
+BFloat
diff --git a/cmake/cpu_features.cmake b/cmake/cpu_features.cmake
index 2bb6deb4847..dbc77d835be 100644
--- a/cmake/cpu_features.cmake
+++ b/cmake/cpu_features.cmake
@@ -85,7 +85,7 @@ elseif (ARCH_AARCH64)
         # [8]  https://developer.arm.com/documentation/102651/a/What-are-dot-product-intructions-
         # [9]  https://developer.arm.com/documentation/dui0801/g/A64-Data-Transfer-Instructions/LDAPR?lang=en
         # [10] https://github.com/aws/aws-graviton-getting-started/blob/main/README.md
-        set (COMPILER_FLAGS "${COMPILER_FLAGS} -march=armv8.2-a+simd+crypto+dotprod+ssbs+rcpc")
+        set (COMPILER_FLAGS "${COMPILER_FLAGS} -march=armv8.2-a+simd+crypto+dotprod+ssbs+rcpc+bf16")
     endif ()
 
     # Best-effort check: The build generates and executes intermediate binaries, e.g. protoc and llvm-tablegen. If we build on ARM for ARM
diff --git a/cmake/linux/default_libs.cmake b/cmake/linux/default_libs.cmake
index 51620bc9f33..79875e1ed6b 100644
--- a/cmake/linux/default_libs.cmake
+++ b/cmake/linux/default_libs.cmake
@@ -3,8 +3,7 @@
 
 set (DEFAULT_LIBS "-nodefaultlibs")
 
-# We need builtins from Clang's RT even without libcxx - for ubsan+int128.
-# See https://bugs.llvm.org/show_bug.cgi?id=16404
+# We need builtins from Clang
 execute_process (COMMAND
     ${CMAKE_CXX_COMPILER} --target=${CMAKE_CXX_COMPILER_TARGET} --print-libgcc-file-name --rtlib=compiler-rt
     OUTPUT_VARIABLE BUILTINS_LIBRARY
diff --git a/docs/en/sql-reference/data-types/float.md b/docs/en/sql-reference/data-types/float.md
index 3c789076c1e..7185308bdce 100644
--- a/docs/en/sql-reference/data-types/float.md
+++ b/docs/en/sql-reference/data-types/float.md
@@ -1,10 +1,10 @@
 ---
 slug: /en/sql-reference/data-types/float
 sidebar_position: 4
-sidebar_label: Float32, Float64
+sidebar_label: Float32, Float64, BFloat16
 ---
 
-# Float32, Float64
+# Float32, Float64, BFloat16
 
 :::note
 If you need accurate calculations, in particular if you work with financial or business data requiring a high precision, you should consider using [Decimal](../data-types/decimal.md) instead. 
@@ -117,3 +117,11 @@ SELECT 0 / 0
 ```
 
 See the rules for `NaN` sorting in the section [ORDER BY clause](../../sql-reference/statements/select/order-by.md).
+
+## BFloat16
+
+`BFloat16` is a 16-bit floating point data type with 8-bit exponent, sign, and 7-bit mantissa.
+
+It is useful for machine learning and AI applications.
+
+ClickHouse supports conversions between `Float32` and `BFloat16`. Most of other operations are not supported.
diff --git a/src/AggregateFunctions/AggregateFunctionAvg.h b/src/AggregateFunctions/AggregateFunctionAvg.h
index 6e1e9289565..8d53a081ee0 100644
--- a/src/AggregateFunctions/AggregateFunctionAvg.h
+++ b/src/AggregateFunctions/AggregateFunctionAvg.h
@@ -231,7 +231,7 @@ public:
 
     void add(AggregateDataPtr __restrict place, const IColumn ** columns, size_t row_num, Arena *) const final
     {
-        increment(place, static_cast<const ColVecType &>(*columns[0]).getData()[row_num]);
+        increment(place, Numerator(static_cast<const ColVecType &>(*columns[0]).getData()[row_num]));
         ++this->data(place).denominator;
     }
 
diff --git a/src/AggregateFunctions/AggregateFunctionDeltaSum.cpp b/src/AggregateFunctions/AggregateFunctionDeltaSum.cpp
index 42169c34c25..c61b9918a35 100644
--- a/src/AggregateFunctions/AggregateFunctionDeltaSum.cpp
+++ b/src/AggregateFunctions/AggregateFunctionDeltaSum.cpp
@@ -27,9 +27,9 @@ namespace
 template <typename T>
 struct AggregationFunctionDeltaSumData
 {
-    T sum = 0;
-    T last = 0;
-    T first = 0;
+    T sum{};
+    T last{};
+    T first{};
     bool seen = false;
 };
 
diff --git a/src/AggregateFunctions/AggregateFunctionDeltaSumTimestamp.cpp b/src/AggregateFunctions/AggregateFunctionDeltaSumTimestamp.cpp
index 5819c533fd9..dc1adead87c 100644
--- a/src/AggregateFunctions/AggregateFunctionDeltaSumTimestamp.cpp
+++ b/src/AggregateFunctions/AggregateFunctionDeltaSumTimestamp.cpp
@@ -25,11 +25,11 @@ namespace
 template <typename ValueType, typename TimestampType>
 struct AggregationFunctionDeltaSumTimestampData
 {
-    ValueType sum = 0;
-    ValueType first = 0;
-    ValueType last = 0;
-    TimestampType first_ts = 0;
-    TimestampType last_ts = 0;
+    ValueType sum{};
+    ValueType first{};
+    ValueType last{};
+    TimestampType first_ts{};
+    TimestampType last_ts{};
     bool seen = false;
 };
 
diff --git a/src/AggregateFunctions/AggregateFunctionGroupArray.cpp b/src/AggregateFunctions/AggregateFunctionGroupArray.cpp
index fb730306c0d..7dc5a4b86b3 100644
--- a/src/AggregateFunctions/AggregateFunctionGroupArray.cpp
+++ b/src/AggregateFunctions/AggregateFunctionGroupArray.cpp
@@ -79,7 +79,7 @@ template <typename T>
 struct GroupArraySamplerData
 {
     /// For easy serialization.
-    static_assert(std::has_unique_object_representations_v<T> || std::is_floating_point_v<T>);
+    static_assert(std::has_unique_object_representations_v<T> || is_floating_point<T>);
 
     // Switch to ordinary Allocator after 4096 bytes to avoid fragmentation and trash in Arena
     using Allocator = MixedAlignedArenaAllocator<alignof(T), 4096>;
@@ -120,7 +120,7 @@ template <typename T>
 struct GroupArrayNumericData<T, false>
 {
     /// For easy serialization.
-    static_assert(std::has_unique_object_representations_v<T> || std::is_floating_point_v<T>);
+    static_assert(std::has_unique_object_representations_v<T> || is_floating_point<T>);
 
     // Switch to ordinary Allocator after 4096 bytes to avoid fragmentation and trash in Arena
     using Allocator = MixedAlignedArenaAllocator<alignof(T), 4096>;
diff --git a/src/AggregateFunctions/AggregateFunctionGroupArrayMoving.cpp b/src/AggregateFunctions/AggregateFunctionGroupArrayMoving.cpp
index 43aad8ac3ba..a6e73729b6a 100644
--- a/src/AggregateFunctions/AggregateFunctionGroupArrayMoving.cpp
+++ b/src/AggregateFunctions/AggregateFunctionGroupArrayMoving.cpp
@@ -38,7 +38,7 @@ template <typename T>
 struct MovingData
 {
     /// For easy serialization.
-    static_assert(std::has_unique_object_representations_v<T> || std::is_floating_point_v<T>);
+    static_assert(std::has_unique_object_representations_v<T> || is_floating_point<T>);
 
     using Accumulator = T;
 
diff --git a/src/AggregateFunctions/AggregateFunctionIntervalLengthSum.cpp b/src/AggregateFunctions/AggregateFunctionIntervalLengthSum.cpp
index eacd0596757..e5404add820 100644
--- a/src/AggregateFunctions/AggregateFunctionIntervalLengthSum.cpp
+++ b/src/AggregateFunctions/AggregateFunctionIntervalLengthSum.cpp
@@ -187,7 +187,7 @@ public:
 
     static DataTypePtr createResultType()
     {
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
             return std::make_shared<DataTypeFloat64>();
         return std::make_shared<DataTypeUInt64>();
     }
@@ -227,7 +227,7 @@ public:
 
     void insertResultInto(AggregateDataPtr __restrict place, IColumn & to, Arena *) const override
     {
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
             assert_cast<ColumnFloat64 &>(to).getData().push_back(getIntervalLengthSum<Float64>(this->data(place)));
         else
             assert_cast<ColumnUInt64 &>(to).getData().push_back(getIntervalLengthSum<UInt64>(this->data(place)));
diff --git a/src/AggregateFunctions/AggregateFunctionMaxIntersections.cpp b/src/AggregateFunctions/AggregateFunctionMaxIntersections.cpp
index ca91f960dab..f4edec7f528 100644
--- a/src/AggregateFunctions/AggregateFunctionMaxIntersections.cpp
+++ b/src/AggregateFunctions/AggregateFunctionMaxIntersections.cpp
@@ -155,9 +155,9 @@ public:
 
     void insertResultInto(AggregateDataPtr __restrict place, IColumn & to, Arena *) const override
     {
-        Int64 current_intersections = 0;
-        Int64 max_intersections = 0;
-        PointType position_of_max_intersections = 0;
+        Int64 current_intersections{};
+        Int64 max_intersections{};
+        PointType position_of_max_intersections{};
 
         /// const_cast because we will sort the array
         auto & array = this->data(place).value;
diff --git a/src/AggregateFunctions/AggregateFunctionSparkbar.cpp b/src/AggregateFunctions/AggregateFunctionSparkbar.cpp
index 362ffbe20d2..de2a741e105 100644
--- a/src/AggregateFunctions/AggregateFunctionSparkbar.cpp
+++ b/src/AggregateFunctions/AggregateFunctionSparkbar.cpp
@@ -45,12 +45,12 @@ struct AggregateFunctionSparkbarData
     Y insert(const X & x, const Y & y)
     {
         if (isNaN(y) || y <= 0)
-            return 0;
+            return {};
 
         auto [it, inserted] = points.insert({x, y});
         if (!inserted)
         {
-            if constexpr (std::is_floating_point_v<Y>)
+            if constexpr (is_floating_point<Y>)
             {
                 it->getMapped() += y;
                 return it->getMapped();
@@ -173,13 +173,13 @@ private:
 
         if (from_x >= to_x)
         {
-            size_t sz = updateFrame(values, 8);
+            size_t sz = updateFrame(values, Y{8});
             values.push_back('\0');
             offsets.push_back(offsets.empty() ? sz + 1 : offsets.back() + sz + 1);
             return;
         }
 
-        PaddedPODArray<Y> histogram(width, 0);
+        PaddedPODArray<Y> histogram(width, Y{0});
         PaddedPODArray<UInt64> count_histogram(width, 0); /// The number of points in each bucket
 
         for (const auto & point : data.points)
@@ -197,7 +197,7 @@ private:
 
             Y res;
             bool has_overfllow = false;
-            if constexpr (std::is_floating_point_v<Y>)
+            if constexpr (is_floating_point<Y>)
                 res = histogram[index] + point.getMapped();
             else
                 has_overfllow = common::addOverflow(histogram[index], point.getMapped(), res);
@@ -218,10 +218,10 @@ private:
         for (size_t i = 0; i < histogram.size(); ++i)
         {
             if (count_histogram[i] > 0)
-                histogram[i] /= count_histogram[i];
+                histogram[i] = histogram[i] / count_histogram[i];
         }
 
-        Y y_max = 0;
+        Y y_max{};
         for (auto & y : histogram)
         {
             if (isNaN(y) || y <= 0)
@@ -245,8 +245,8 @@ private:
                 continue;
             }
 
-            constexpr auto levels_num = static_cast<Y>(BAR_LEVELS - 1);
-            if constexpr (std::is_floating_point_v<Y>)
+            constexpr auto levels_num = Y{BAR_LEVELS - 1};
+            if constexpr (is_floating_point<Y>)
             {
                 y = y / (y_max / levels_num) + 1;
             }
diff --git a/src/AggregateFunctions/AggregateFunctionSum.h b/src/AggregateFunctions/AggregateFunctionSum.h
index 2ce03c530c2..7c7fb6338a2 100644
--- a/src/AggregateFunctions/AggregateFunctionSum.h
+++ b/src/AggregateFunctions/AggregateFunctionSum.h
@@ -69,7 +69,7 @@ struct AggregateFunctionSumData
         size_t count = end - start;
         const auto * end_ptr = ptr + count;
 
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
         {
             /// Compiler cannot unroll this loop, do it manually.
             /// (at least for floats, most likely due to the lack of -fassociative-math)
@@ -83,7 +83,7 @@ struct AggregateFunctionSumData
             while (ptr < unrolled_end)
             {
                 for (size_t i = 0; i < unroll_count; ++i)
-                    Impl::add(partial_sums[i], ptr[i]);
+                    Impl::add(partial_sums[i], T(ptr[i]));
                 ptr += unroll_count;
             }
 
@@ -95,7 +95,7 @@ struct AggregateFunctionSumData
         T local_sum{};
         while (ptr < end_ptr)
         {
-            Impl::add(local_sum, *ptr);
+            Impl::add(local_sum, T(*ptr));
             ++ptr;
         }
         Impl::add(sum, local_sum);
@@ -193,12 +193,11 @@ struct AggregateFunctionSumData
             Impl::add(sum, local_sum);
             return;
         }
-        else if constexpr (std::is_floating_point_v<T>)
+        else if constexpr (is_floating_point<T> && (sizeof(Value) == 4 || sizeof(Value) == 8))
         {
-            /// For floating point we use a similar trick as above, except that now we  reinterpret the floating point number as an unsigned
+            /// For floating point we use a similar trick as above, except that now we reinterpret the floating point number as an unsigned
             /// integer of the same size and use a mask instead (0 to discard, 0xFF..FF to keep)
-            static_assert(sizeof(Value) == 4 || sizeof(Value) == 8);
-            using equivalent_integer = typename std::conditional_t<sizeof(Value) == 4, UInt32, UInt64>;
+            using EquivalentInteger = typename std::conditional_t<sizeof(Value) == 4, UInt32, UInt64>;
 
             constexpr size_t unroll_count = 128 / sizeof(T);
             T partial_sums[unroll_count]{};
@@ -209,11 +208,11 @@ struct AggregateFunctionSumData
             {
                 for (size_t i = 0; i < unroll_count; ++i)
                 {
-                    equivalent_integer value;
-                    std::memcpy(&value, &ptr[i], sizeof(Value));
+                    EquivalentInteger value;
+                    memcpy(&value, &ptr[i], sizeof(Value));
                     value &= (!condition_map[i] != add_if_zero) - 1;
                     Value d;
-                    std::memcpy(&d, &value, sizeof(Value));
+                    memcpy(&d, &value, sizeof(Value));
                     Impl::add(partial_sums[i], d);
                 }
                 ptr += unroll_count;
@@ -228,7 +227,7 @@ struct AggregateFunctionSumData
         while (ptr < end_ptr)
         {
             if (!*condition_map == add_if_zero)
-                Impl::add(local_sum, *ptr);
+                Impl::add(local_sum, T(*ptr));
             ++ptr;
             ++condition_map;
         }
@@ -306,7 +305,7 @@ struct AggregateFunctionSumData
 template <typename T>
 struct AggregateFunctionSumKahanData
 {
-    static_assert(std::is_floating_point_v<T>,
+    static_assert(is_floating_point<T>,
         "It doesn't make sense to use Kahan Summation algorithm for non floating point types");
 
     T sum{};
@@ -489,10 +488,7 @@ public:
     void add(AggregateDataPtr __restrict place, const IColumn ** columns, size_t row_num, Arena *) const override
     {
         const auto & column = assert_cast<const ColVecType &>(*columns[0]);
-        if constexpr (is_big_int_v<T>)
-            this->data(place).add(static_cast<TResult>(column.getData()[row_num]));
-        else
-            this->data(place).add(column.getData()[row_num]);
+        this->data(place).add(static_cast<TResult>(column.getData()[row_num]));
     }
 
     void addBatchSinglePlace(
diff --git a/src/AggregateFunctions/AggregateFunctionUniq.h b/src/AggregateFunctions/AggregateFunctionUniq.h
index 35d6e599e38..a094ed288d6 100644
--- a/src/AggregateFunctions/AggregateFunctionUniq.h
+++ b/src/AggregateFunctions/AggregateFunctionUniq.h
@@ -257,7 +257,7 @@ template <typename T> struct AggregateFunctionUniqTraits
 {
     static UInt64 hash(T x)
     {
-        if constexpr (std::is_same_v<T, Float32> || std::is_same_v<T, Float64>)
+        if constexpr (is_floating_point<T>)
         {
             return bit_cast<UInt64>(x);
         }
diff --git a/src/AggregateFunctions/AggregateFunctionUniqCombined.h b/src/AggregateFunctions/AggregateFunctionUniqCombined.h
index c6377a9532d..4a91991ebb5 100644
--- a/src/AggregateFunctions/AggregateFunctionUniqCombined.h
+++ b/src/AggregateFunctions/AggregateFunctionUniqCombined.h
@@ -111,7 +111,7 @@ public:
                 /// Initially UInt128 was introduced only for UUID, and then the other big-integer types were added.
                 hash = static_cast<HashValueType>(sipHash64(value));
             }
-            else if constexpr (std::is_floating_point_v<T>)
+            else if constexpr (is_floating_point<T>)
             {
                 hash = static_cast<HashValueType>(intHash64(bit_cast<UInt64>(value)));
             }
diff --git a/src/AggregateFunctions/QuantileTDigest.h b/src/AggregateFunctions/QuantileTDigest.h
index f915d6c70a3..8117b2a7fe1 100644
--- a/src/AggregateFunctions/QuantileTDigest.h
+++ b/src/AggregateFunctions/QuantileTDigest.h
@@ -391,7 +391,7 @@ public:
     ResultType getImpl(Float64 level)
     {
         if (centroids.empty())
-            return std::is_floating_point_v<ResultType> ? std::numeric_limits<ResultType>::quiet_NaN() : 0;
+            return is_floating_point<ResultType> ? std::numeric_limits<ResultType>::quiet_NaN() : 0;
 
         compress();
 
diff --git a/src/AggregateFunctions/ReservoirSampler.h b/src/AggregateFunctions/ReservoirSampler.h
index 2668e0dc890..870cb429fb7 100644
--- a/src/AggregateFunctions/ReservoirSampler.h
+++ b/src/AggregateFunctions/ReservoirSampler.h
@@ -276,6 +276,6 @@ private:
     {
         if (OnEmpty == ReservoirSamplerOnEmpty::THROW)
             throw DB::Exception(DB::ErrorCodes::LOGICAL_ERROR, "Quantile of empty ReservoirSampler");
-        return NanLikeValueConstructor<ResultType, std::is_floating_point_v<ResultType>>::getValue();
+        return NanLikeValueConstructor<ResultType, is_floating_point<ResultType>>::getValue();
     }
 };
diff --git a/src/AggregateFunctions/ReservoirSamplerDeterministic.h b/src/AggregateFunctions/ReservoirSamplerDeterministic.h
index 0ad8e968358..f35d1eb49af 100644
--- a/src/AggregateFunctions/ReservoirSamplerDeterministic.h
+++ b/src/AggregateFunctions/ReservoirSamplerDeterministic.h
@@ -271,7 +271,7 @@ private:
     {
         if (OnEmpty == ReservoirSamplerDeterministicOnEmpty::THROW)
             throw DB::Exception(DB::ErrorCodes::LOGICAL_ERROR, "Quantile of empty ReservoirSamplerDeterministic");
-        return NanLikeValueConstructor<ResultType, std::is_floating_point_v<ResultType>>::getValue();
+        return NanLikeValueConstructor<ResultType, is_floating_point<ResultType>>::getValue();
     }
 };
 
diff --git a/src/Columns/ColumnArray.cpp b/src/Columns/ColumnArray.cpp
index 0c6d7c4e5c6..3f88ca93a97 100644
--- a/src/Columns/ColumnArray.cpp
+++ b/src/Columns/ColumnArray.cpp
@@ -662,6 +662,8 @@ ColumnPtr ColumnArray::filter(const Filter & filt, ssize_t result_size_hint) con
         return filterNumber<Int128>(filt, result_size_hint);
     if (typeid_cast<const ColumnInt256 *>(data.get()))
         return filterNumber<Int256>(filt, result_size_hint);
+    if (typeid_cast<const ColumnBFloat16 *>(data.get()))
+        return filterNumber<BFloat16>(filt, result_size_hint);
     if (typeid_cast<const ColumnFloat32 *>(data.get()))
         return filterNumber<Float32>(filt, result_size_hint);
     if (typeid_cast<const ColumnFloat64 *>(data.get()))
@@ -1065,6 +1067,8 @@ ColumnPtr ColumnArray::replicate(const Offsets & replicate_offsets) const
         return replicateNumber<Int128>(replicate_offsets);
     if (typeid_cast<const ColumnInt256 *>(data.get()))
         return replicateNumber<Int256>(replicate_offsets);
+    if (typeid_cast<const ColumnBFloat16 *>(data.get()))
+        return replicateNumber<BFloat16>(replicate_offsets);
     if (typeid_cast<const ColumnFloat32 *>(data.get()))
         return replicateNumber<Float32>(replicate_offsets);
     if (typeid_cast<const ColumnFloat64 *>(data.get()))
diff --git a/src/Columns/ColumnUnique.cpp b/src/Columns/ColumnUnique.cpp
index 54f45204c00..773edbfd590 100644
--- a/src/Columns/ColumnUnique.cpp
+++ b/src/Columns/ColumnUnique.cpp
@@ -16,6 +16,7 @@ template class ColumnUnique<ColumnInt128>;
 template class ColumnUnique<ColumnUInt128>;
 template class ColumnUnique<ColumnInt256>;
 template class ColumnUnique<ColumnUInt256>;
+template class ColumnUnique<ColumnBFloat16>;
 template class ColumnUnique<ColumnFloat32>;
 template class ColumnUnique<ColumnFloat64>;
 template class ColumnUnique<ColumnString>;
diff --git a/src/Columns/ColumnUnique.h b/src/Columns/ColumnUnique.h
index ffa7c311e9e..ce7bbf0766f 100644
--- a/src/Columns/ColumnUnique.h
+++ b/src/Columns/ColumnUnique.h
@@ -760,6 +760,7 @@ extern template class ColumnUnique<ColumnInt128>;
 extern template class ColumnUnique<ColumnUInt128>;
 extern template class ColumnUnique<ColumnInt256>;
 extern template class ColumnUnique<ColumnUInt256>;
+extern template class ColumnUnique<ColumnBFloat16>;
 extern template class ColumnUnique<ColumnFloat32>;
 extern template class ColumnUnique<ColumnFloat64>;
 extern template class ColumnUnique<ColumnString>;
diff --git a/src/Columns/ColumnVector.cpp b/src/Columns/ColumnVector.cpp
index dc4941988e9..84fc6ebc61d 100644
--- a/src/Columns/ColumnVector.cpp
+++ b/src/Columns/ColumnVector.cpp
@@ -118,9 +118,9 @@ struct ColumnVector<T>::less_stable
         if (unlikely(parent.data[lhs] == parent.data[rhs]))
             return lhs < rhs;
 
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
         {
-            if (unlikely(std::isnan(parent.data[lhs]) && std::isnan(parent.data[rhs])))
+            if (unlikely(isNaN(parent.data[lhs]) && isNaN(parent.data[rhs])))
             {
                 return lhs < rhs;
             }
@@ -150,9 +150,9 @@ struct ColumnVector<T>::greater_stable
         if (unlikely(parent.data[lhs] == parent.data[rhs]))
             return lhs < rhs;
 
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
         {
-            if (unlikely(std::isnan(parent.data[lhs]) && std::isnan(parent.data[rhs])))
+            if (unlikely(isNaN(parent.data[lhs]) && isNaN(parent.data[rhs])))
             {
                 return lhs < rhs;
             }
@@ -224,9 +224,9 @@ void ColumnVector<T>::getPermutation(IColumn::PermutationSortDirection direction
 
     iota(res.data(), data_size, IColumn::Permutation::value_type(0));
 
-    if constexpr (has_find_extreme_implementation<T> && !std::is_floating_point_v<T>)
+    if constexpr (has_find_extreme_implementation<T> && !is_floating_point<T>)
     {
-        /// Disabled for:floating point
+        /// Disabled for floating point:
         /// * floating point: We don't deal with nan_direction_hint
         /// * stability::Stable: We might return any value, not the first
         if ((limit == 1) && (stability == IColumn::PermutationSortStability::Unstable))
@@ -256,7 +256,7 @@ void ColumnVector<T>::getPermutation(IColumn::PermutationSortDirection direction
             bool sort_is_stable = stability == IColumn::PermutationSortStability::Stable;
 
             /// TODO: LSD RadixSort is currently not stable if direction is descending, or value is floating point
-            bool use_radix_sort = (sort_is_stable && ascending && !std::is_floating_point_v<T>) || !sort_is_stable;
+            bool use_radix_sort = (sort_is_stable && ascending && !is_floating_point<T>) || !sort_is_stable;
 
             /// Thresholds on size. Lower threshold is arbitrary. Upper threshold is chosen by the type for histogram counters.
             if (data_size >= 256 && data_size <= std::numeric_limits<UInt32>::max() && use_radix_sort)
@@ -283,7 +283,7 @@ void ColumnVector<T>::getPermutation(IColumn::PermutationSortDirection direction
 
                 /// Radix sort treats all NaNs to be greater than all numbers.
                 /// If the user needs the opposite, we must move them accordingly.
-                if (std::is_floating_point_v<T> && nan_direction_hint < 0)
+                if (is_floating_point<T> && nan_direction_hint < 0)
                 {
                     size_t nans_to_move = 0;
 
@@ -330,7 +330,7 @@ void ColumnVector<T>::updatePermutation(IColumn::PermutationSortDirection direct
         if constexpr (is_arithmetic_v<T> && !is_big_int_v<T>)
         {
             /// TODO: LSD RadixSort is currently not stable if direction is descending, or value is floating point
-            bool use_radix_sort = (sort_is_stable && ascending && !std::is_floating_point_v<T>) || !sort_is_stable;
+            bool use_radix_sort = (sort_is_stable && ascending && !is_floating_point<T>) || !sort_is_stable;
             size_t size = end - begin;
 
             /// Thresholds on size. Lower threshold is arbitrary. Upper threshold is chosen by the type for histogram counters.
@@ -353,7 +353,7 @@ void ColumnVector<T>::updatePermutation(IColumn::PermutationSortDirection direct
 
                 /// Radix sort treats all NaNs to be greater than all numbers.
                 /// If the user needs the opposite, we must move them accordingly.
-                if (std::is_floating_point_v<T> && nan_direction_hint < 0)
+                if (is_floating_point<T> && nan_direction_hint < 0)
                 {
                     size_t nans_to_move = 0;
 
@@ -1005,6 +1005,7 @@ template class ColumnVector<Int32>;
 template class ColumnVector<Int64>;
 template class ColumnVector<Int128>;
 template class ColumnVector<Int256>;
+template class ColumnVector<BFloat16>;
 template class ColumnVector<Float32>;
 template class ColumnVector<Float64>;
 template class ColumnVector<UUID>;
diff --git a/src/Columns/ColumnVector.h b/src/Columns/ColumnVector.h
index 8f81da86375..47840a90a0e 100644
--- a/src/Columns/ColumnVector.h
+++ b/src/Columns/ColumnVector.h
@@ -481,6 +481,7 @@ extern template class ColumnVector<Int32>;
 extern template class ColumnVector<Int64>;
 extern template class ColumnVector<Int128>;
 extern template class ColumnVector<Int256>;
+extern template class ColumnVector<BFloat16>;
 extern template class ColumnVector<Float32>;
 extern template class ColumnVector<Float64>;
 extern template class ColumnVector<UUID>;
diff --git a/src/Columns/ColumnsCommon.cpp b/src/Columns/ColumnsCommon.cpp
index 591c56a00bd..3e0ead39486 100644
--- a/src/Columns/ColumnsCommon.cpp
+++ b/src/Columns/ColumnsCommon.cpp
@@ -328,6 +328,7 @@ INSTANTIATE(Int32)
 INSTANTIATE(Int64)
 INSTANTIATE(Int128)
 INSTANTIATE(Int256)
+INSTANTIATE(BFloat16)
 INSTANTIATE(Float32)
 INSTANTIATE(Float64)
 INSTANTIATE(Decimal32)
diff --git a/src/Columns/ColumnsNumber.h b/src/Columns/ColumnsNumber.h
index ae7eddb0b22..2dce2269079 100644
--- a/src/Columns/ColumnsNumber.h
+++ b/src/Columns/ColumnsNumber.h
@@ -23,6 +23,7 @@ using ColumnInt64 = ColumnVector<Int64>;
 using ColumnInt128 = ColumnVector<Int128>;
 using ColumnInt256 = ColumnVector<Int256>;
 
+using ColumnBFloat16 = ColumnVector<BFloat16>;
 using ColumnFloat32 = ColumnVector<Float32>;
 using ColumnFloat64 = ColumnVector<Float64>;
 
diff --git a/src/Columns/IColumn.cpp b/src/Columns/IColumn.cpp
index c9a0514af4e..4a3886dddb6 100644
--- a/src/Columns/IColumn.cpp
+++ b/src/Columns/IColumn.cpp
@@ -443,6 +443,7 @@ template class IColumnHelper<ColumnVector<Int32>, ColumnFixedSizeHelper>;
 template class IColumnHelper<ColumnVector<Int64>, ColumnFixedSizeHelper>;
 template class IColumnHelper<ColumnVector<Int128>, ColumnFixedSizeHelper>;
 template class IColumnHelper<ColumnVector<Int256>, ColumnFixedSizeHelper>;
+template class IColumnHelper<ColumnVector<BFloat16>, ColumnFixedSizeHelper>;
 template class IColumnHelper<ColumnVector<Float32>, ColumnFixedSizeHelper>;
 template class IColumnHelper<ColumnVector<Float64>, ColumnFixedSizeHelper>;
 template class IColumnHelper<ColumnVector<UUID>, ColumnFixedSizeHelper>;
diff --git a/src/Columns/MaskOperations.cpp b/src/Columns/MaskOperations.cpp
index 873a4060872..eac42400517 100644
--- a/src/Columns/MaskOperations.cpp
+++ b/src/Columns/MaskOperations.cpp
@@ -63,6 +63,7 @@ INSTANTIATE(Int32)
 INSTANTIATE(Int64)
 INSTANTIATE(Int128)
 INSTANTIATE(Int256)
+INSTANTIATE(BFloat16)
 INSTANTIATE(Float32)
 INSTANTIATE(Float64)
 INSTANTIATE(Decimal32)
@@ -200,6 +201,7 @@ static MaskInfo extractMaskImpl(
           || extractMaskNumeric<inverted, Int16>(mask, column, null_value, null_bytemap, nulls, mask_info)
           || extractMaskNumeric<inverted, Int32>(mask, column, null_value, null_bytemap, nulls, mask_info)
           || extractMaskNumeric<inverted, Int64>(mask, column, null_value, null_bytemap, nulls, mask_info)
+          || extractMaskNumeric<inverted, BFloat16>(mask, column, null_value, null_bytemap, nulls, mask_info)
           || extractMaskNumeric<inverted, Float32>(mask, column, null_value, null_bytemap, nulls, mask_info)
           || extractMaskNumeric<inverted, Float64>(mask, column, null_value, null_bytemap, nulls, mask_info)))
         throw Exception(ErrorCodes::ILLEGAL_COLUMN, "Cannot convert column {} to mask.", column->getName());
diff --git a/src/Columns/tests/gtest_column_vector.cpp b/src/Columns/tests/gtest_column_vector.cpp
index b71d4a095ab..3a084a89079 100644
--- a/src/Columns/tests/gtest_column_vector.cpp
+++ b/src/Columns/tests/gtest_column_vector.cpp
@@ -93,6 +93,7 @@ TEST(ColumnVector, Filter)
     testFilter<Int64>();
     testFilter<UInt128>();
     testFilter<Int256>();
+    testFilter<BFloat16>();
     testFilter<Float32>();
     testFilter<Float64>();
     testFilter<UUID>();
diff --git a/src/Columns/tests/gtest_low_cardinality.cpp b/src/Columns/tests/gtest_low_cardinality.cpp
index ce16d2cadb1..301fa2a6090 100644
--- a/src/Columns/tests/gtest_low_cardinality.cpp
+++ b/src/Columns/tests/gtest_low_cardinality.cpp
@@ -45,6 +45,7 @@ TEST(ColumnLowCardinality, Insert)
     testLowCardinalityNumberInsert<Int128>(std::make_shared<DataTypeInt128>());
     testLowCardinalityNumberInsert<Int256>(std::make_shared<DataTypeInt256>());
 
+    testLowCardinalityNumberInsert<BFloat16>(std::make_shared<DataTypeBFloat16>());
     testLowCardinalityNumberInsert<Float32>(std::make_shared<DataTypeFloat32>());
     testLowCardinalityNumberInsert<Float64>(std::make_shared<DataTypeFloat64>());
 }
diff --git a/src/Common/CPUID.h b/src/Common/CPUID.h
index b49f7706904..b5c26e64d1e 100644
--- a/src/Common/CPUID.h
+++ b/src/Common/CPUID.h
@@ -266,6 +266,11 @@ inline bool haveAVX512VBMI2() noexcept
     return haveAVX512F() && ((CPUInfo(0x7, 0).registers.ecx >> 6) & 1u);
 }
 
+inline bool haveAVX512BF16() noexcept
+{
+    return haveAVX512F() && ((CPUInfo(0x7, 1).registers.eax >> 5) & 1u);
+}
+
 inline bool haveRDRAND() noexcept
 {
     return CPUInfo(0x0).registers.eax >= 0x7 && ((CPUInfo(0x1).registers.ecx >> 30) & 1u);
@@ -326,6 +331,7 @@ inline bool haveAMXINT8() noexcept
     OP(AVX512VL)             \
     OP(AVX512VBMI)           \
     OP(AVX512VBMI2)          \
+    OP(AVX512BF16)           \
     OP(PREFETCHWT1)          \
     OP(SHA)                  \
     OP(ADX)                  \
diff --git a/src/Common/FieldVisitorConvertToNumber.cpp b/src/Common/FieldVisitorConvertToNumber.cpp
index 75b3fbfe02a..a5963e3d028 100644
--- a/src/Common/FieldVisitorConvertToNumber.cpp
+++ b/src/Common/FieldVisitorConvertToNumber.cpp
@@ -1,5 +1,4 @@
 #include <Common/FieldVisitorConvertToNumber.h>
-#include "base/Decimal.h"
 
 namespace DB
 {
@@ -17,6 +16,7 @@ template class FieldVisitorConvertToNumber<Int128>;
 template class FieldVisitorConvertToNumber<UInt128>;
 template class FieldVisitorConvertToNumber<Int256>;
 template class FieldVisitorConvertToNumber<UInt256>;
+//template class FieldVisitorConvertToNumber<BFloat16>;
 template class FieldVisitorConvertToNumber<Float32>;
 template class FieldVisitorConvertToNumber<Float64>;
 
diff --git a/src/Common/FieldVisitorConvertToNumber.h b/src/Common/FieldVisitorConvertToNumber.h
index 4c3424e0e18..38d5dc473c4 100644
--- a/src/Common/FieldVisitorConvertToNumber.h
+++ b/src/Common/FieldVisitorConvertToNumber.h
@@ -58,7 +58,7 @@ public:
 
     T operator() (const Float64 & x) const
     {
-        if constexpr (!std::is_floating_point_v<T>)
+        if constexpr (!is_floating_point<T>)
         {
             if (!isFinite(x))
             {
@@ -88,7 +88,7 @@ public:
     template <typename U>
     T operator() (const DecimalField<U> & x) const
     {
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
             return x.getValue().template convertTo<T>() / x.getScaleMultiplier().template convertTo<T>();
         else
             return (x.getValue() / x.getScaleMultiplier()).template convertTo<T>();
@@ -129,6 +129,7 @@ extern template class FieldVisitorConvertToNumber<Int128>;
 extern template class FieldVisitorConvertToNumber<UInt128>;
 extern template class FieldVisitorConvertToNumber<Int256>;
 extern template class FieldVisitorConvertToNumber<UInt256>;
+//extern template class FieldVisitorConvertToNumber<BFloat16>;
 extern template class FieldVisitorConvertToNumber<Float32>;
 extern template class FieldVisitorConvertToNumber<Float64>;
 
diff --git a/src/Common/HashTable/Hash.h b/src/Common/HashTable/Hash.h
index fb6afcde133..b4bc6af1cef 100644
--- a/src/Common/HashTable/Hash.h
+++ b/src/Common/HashTable/Hash.h
@@ -322,6 +322,7 @@ DEFINE_HASH(Int32)
 DEFINE_HASH(Int64)
 DEFINE_HASH(Int128)
 DEFINE_HASH(Int256)
+DEFINE_HASH(BFloat16)
 DEFINE_HASH(Float32)
 DEFINE_HASH(Float64)
 DEFINE_HASH(DB::UUID)
diff --git a/src/Common/HashTable/HashTable.h b/src/Common/HashTable/HashTable.h
index d379c3f6a87..fb34ab7d4f9 100644
--- a/src/Common/HashTable/HashTable.h
+++ b/src/Common/HashTable/HashTable.h
@@ -76,7 +76,7 @@ struct HashTableNoState
 template <typename T>
 inline bool bitEquals(T a, T b)
 {
-    if constexpr (std::is_floating_point_v<T>)
+    if constexpr (is_floating_point<T>)
         /// Note that memcmp with constant size is a compiler builtin.
         return 0 == memcmp(&a, &b, sizeof(T)); /// NOLINT
     else
diff --git a/src/Common/NaNUtils.h b/src/Common/NaNUtils.h
index 2ff3e2f5661..3b99084549f 100644
--- a/src/Common/NaNUtils.h
+++ b/src/Common/NaNUtils.h
@@ -3,24 +3,24 @@
 #include <cmath>
 #include <limits>
 #include <type_traits>
+#include <base/DecomposedFloat.h>
 
 
 template <typename T>
 inline bool isNaN(T x)
 {
     /// To be sure, that this function is zero-cost for non-floating point types.
-    if constexpr (std::is_floating_point_v<T>)
-        return std::isnan(x);
+    if constexpr (is_floating_point<T>)
+        return DecomposedFloat(x).isNaN();
     else
         return false;
 }
 
-
 template <typename T>
 inline bool isFinite(T x)
 {
-    if constexpr (std::is_floating_point_v<T>)
-        return std::isfinite(x);
+    if constexpr (is_floating_point<T>)
+        return DecomposedFloat(x).isFinite();
     else
         return true;
 }
@@ -28,7 +28,7 @@ inline bool isFinite(T x)
 template <typename T>
 bool canConvertTo(Float64 x)
 {
-    if constexpr (std::is_floating_point_v<T>)
+    if constexpr (is_floating_point<T>)
         return true;
     if (!isFinite(x))
         return false;
@@ -46,3 +46,12 @@ T NaNOrZero()
     else
         return {};
 }
+
+template <typename T>
+bool signBit(T x)
+{
+    if constexpr (is_floating_point<T>)
+        return DecomposedFloat(x).isNegative();
+    else
+        return x < 0;
+}
diff --git a/src/Common/TargetSpecific.cpp b/src/Common/TargetSpecific.cpp
index 8540c9a9986..4400d9a60b3 100644
--- a/src/Common/TargetSpecific.cpp
+++ b/src/Common/TargetSpecific.cpp
@@ -23,6 +23,8 @@ UInt32 getSupportedArchs()
         result |= static_cast<UInt32>(TargetArch::AVX512VBMI);
     if (CPU::CPUFlagsCache::have_AVX512VBMI2)
         result |= static_cast<UInt32>(TargetArch::AVX512VBMI2);
+    if (CPU::CPUFlagsCache::have_AVX512BF16)
+        result |= static_cast<UInt32>(TargetArch::AVX512BF16);
     if (CPU::CPUFlagsCache::have_AMXBF16)
         result |= static_cast<UInt32>(TargetArch::AMXBF16);
     if (CPU::CPUFlagsCache::have_AMXTILE)
@@ -50,6 +52,7 @@ String toString(TargetArch arch)
         case TargetArch::AVX512BW:    return "avx512bw";
         case TargetArch::AVX512VBMI:  return "avx512vbmi";
         case TargetArch::AVX512VBMI2: return "avx512vbmi2";
+        case TargetArch::AVX512BF16:  return "avx512bf16";
         case TargetArch::AMXBF16: return "amxbf16";
         case TargetArch::AMXTILE: return "amxtile";
         case TargetArch::AMXINT8: return "amxint8";
diff --git a/src/Common/TargetSpecific.h b/src/Common/TargetSpecific.h
index f9523f667b2..5584bd1f63a 100644
--- a/src/Common/TargetSpecific.h
+++ b/src/Common/TargetSpecific.h
@@ -83,9 +83,10 @@ enum class TargetArch : UInt32
     AVX512BW    = (1 << 4),
     AVX512VBMI  = (1 << 5),
     AVX512VBMI2 = (1 << 6),
-    AMXBF16 = (1 << 7),
-    AMXTILE = (1 << 8),
-    AMXINT8 = (1 << 9),
+    AVX512BF16 = (1 << 7),
+    AMXBF16 = (1 << 8),
+    AMXTILE = (1 << 9),
+    AMXINT8 = (1 << 10),
 };
 
 /// Runtime detection.
@@ -102,6 +103,7 @@ String toString(TargetArch arch);
 /// NOLINTNEXTLINE
 #define USE_MULTITARGET_CODE 1
 
+#define AVX512BF16_FUNCTION_SPECIFIC_ATTRIBUTE __attribute__((target("sse,sse2,sse3,ssse3,sse4,popcnt,avx,avx2,avx512f,avx512bw,avx512vl,avx512vbmi,avx512vbmi2,avx512bf16")))
 #define AVX512VBMI2_FUNCTION_SPECIFIC_ATTRIBUTE __attribute__((target("sse,sse2,sse3,ssse3,sse4,popcnt,avx,avx2,avx512f,avx512bw,avx512vl,avx512vbmi,avx512vbmi2")))
 #define AVX512VBMI_FUNCTION_SPECIFIC_ATTRIBUTE __attribute__((target("sse,sse2,sse3,ssse3,sse4,popcnt,avx,avx2,avx512f,avx512bw,avx512vl,avx512vbmi")))
 #define AVX512BW_FUNCTION_SPECIFIC_ATTRIBUTE __attribute__((target("sse,sse2,sse3,ssse3,sse4,popcnt,avx,avx2,avx512f,avx512bw")))
@@ -111,6 +113,8 @@ String toString(TargetArch arch);
 #define SSE42_FUNCTION_SPECIFIC_ATTRIBUTE __attribute__((target("sse,sse2,sse3,ssse3,sse4,popcnt")))
 #define DEFAULT_FUNCTION_SPECIFIC_ATTRIBUTE
 
+#   define BEGIN_AVX512BF16_SPECIFIC_CODE \
+        _Pragma("clang attribute push(__attribute__((target(\"sse,sse2,sse3,ssse3,sse4,popcnt,avx,avx2,avx512f,avx512bw,avx512vl,avx512vbmi,avx512vbmi2,avx512bf16\"))),apply_to=function)")
 #   define BEGIN_AVX512VBMI2_SPECIFIC_CODE \
         _Pragma("clang attribute push(__attribute__((target(\"sse,sse2,sse3,ssse3,sse4,popcnt,avx,avx2,avx512f,avx512bw,avx512vl,avx512vbmi,avx512vbmi2\"))),apply_to=function)")
 #   define BEGIN_AVX512VBMI_SPECIFIC_CODE \
@@ -197,6 +201,14 @@ namespace TargetSpecific::AVX512VBMI2 { \
 } \
 END_TARGET_SPECIFIC_CODE
 
+#define DECLARE_AVX512BF16_SPECIFIC_CODE(...) \
+BEGIN_AVX512BF16_SPECIFIC_CODE \
+namespace TargetSpecific::AVX512BF16 { \
+    DUMMY_FUNCTION_DEFINITION \
+    using namespace DB::TargetSpecific::AVX512BF16; \
+    __VA_ARGS__ \
+} \
+END_TARGET_SPECIFIC_CODE
 
 #else
 
@@ -211,6 +223,7 @@ END_TARGET_SPECIFIC_CODE
 #define DECLARE_AVX512BW_SPECIFIC_CODE(...)
 #define DECLARE_AVX512VBMI_SPECIFIC_CODE(...)
 #define DECLARE_AVX512VBMI2_SPECIFIC_CODE(...)
+#define DECLARE_AVX512BF16_SPECIFIC_CODE(...)
 
 #endif
 
@@ -229,7 +242,8 @@ DECLARE_AVX2_SPECIFIC_CODE   (__VA_ARGS__) \
 DECLARE_AVX512F_SPECIFIC_CODE(__VA_ARGS__) \
 DECLARE_AVX512BW_SPECIFIC_CODE    (__VA_ARGS__) \
 DECLARE_AVX512VBMI_SPECIFIC_CODE  (__VA_ARGS__) \
-DECLARE_AVX512VBMI2_SPECIFIC_CODE (__VA_ARGS__)
+DECLARE_AVX512VBMI2_SPECIFIC_CODE (__VA_ARGS__) \
+DECLARE_AVX512BF16_SPECIFIC_CODE (__VA_ARGS__)
 
 DECLARE_DEFAULT_CODE(
     constexpr auto BuildArch = TargetArch::Default; /// NOLINT
@@ -263,6 +277,10 @@ DECLARE_AVX512VBMI2_SPECIFIC_CODE(
     constexpr auto BuildArch = TargetArch::AVX512VBMI2; /// NOLINT
 ) // DECLARE_AVX512VBMI2_SPECIFIC_CODE
 
+DECLARE_AVX512BF16_SPECIFIC_CODE(
+    constexpr auto BuildArch = TargetArch::AVX512BF16; /// NOLINT
+) // DECLARE_AVX512BF16_SPECIFIC_CODE
+
 /** Runtime Dispatch helpers for class members.
   *
   * Example of usage:
diff --git a/src/Common/findExtreme.cpp b/src/Common/findExtreme.cpp
index ce3bbb86d7c..a29750b848a 100644
--- a/src/Common/findExtreme.cpp
+++ b/src/Common/findExtreme.cpp
@@ -47,7 +47,7 @@ MULTITARGET_FUNCTION_AVX2_SSE42(
 
         /// Unroll the loop manually for floating point, since the compiler doesn't do it without fastmath
         /// as it might change the return value
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
         {
             constexpr size_t unroll_block = 512 / sizeof(T); /// Chosen via benchmarks with AVX2 so YMMV
             size_t unrolled_end = i + (((count - i) / unroll_block) * unroll_block);
diff --git a/src/Common/transformEndianness.h b/src/Common/transformEndianness.h
index 1657305acda..e6e04ec75af 100644
--- a/src/Common/transformEndianness.h
+++ b/src/Common/transformEndianness.h
@@ -38,7 +38,7 @@ inline void transformEndianness(T & x)
 }
 
 template <std::endian ToEndian, std::endian FromEndian = std::endian::native, typename T>
-requires std::is_floating_point_v<T>
+requires is_floating_point<T>
 inline void transformEndianness(T & value)
 {
     if constexpr (ToEndian != FromEndian)
diff --git a/src/Core/AccurateComparison.h b/src/Core/AccurateComparison.h
index 139ee4d88dc..87ff14e40e7 100644
--- a/src/Core/AccurateComparison.h
+++ b/src/Core/AccurateComparison.h
@@ -25,7 +25,7 @@ bool lessOp(A a, B b)
         return a < b;
 
     /// float vs float
-    if constexpr (std::is_floating_point_v<A> && std::is_floating_point_v<B>)
+    if constexpr (is_floating_point<A> && is_floating_point<B>)
         return a < b;
 
     /// anything vs NaN
@@ -49,7 +49,7 @@ bool lessOp(A a, B b)
     }
 
     /// int vs float
-    if constexpr (is_integer<A> && std::is_floating_point_v<B>)
+    if constexpr (is_integer<A> && is_floating_point<B>)
     {
         if constexpr (sizeof(A) <= 4)
             return static_cast<double>(a) < static_cast<double>(b);
@@ -57,7 +57,7 @@ bool lessOp(A a, B b)
         return DecomposedFloat<B>(b).greater(a);
     }
 
-    if constexpr (std::is_floating_point_v<A> && is_integer<B>)
+    if constexpr (is_floating_point<A> && is_integer<B>)
     {
         if constexpr (sizeof(B) <= 4)
             return static_cast<double>(a) < static_cast<double>(b);
@@ -65,8 +65,8 @@ bool lessOp(A a, B b)
         return DecomposedFloat<A>(a).less(b);
     }
 
-    static_assert(is_integer<A> || std::is_floating_point_v<A>);
-    static_assert(is_integer<B> || std::is_floating_point_v<B>);
+    static_assert(is_integer<A> || is_floating_point<A>);
+    static_assert(is_integer<B> || is_floating_point<B>);
     UNREACHABLE();
 }
 
@@ -101,7 +101,7 @@ bool equalsOp(A a, B b)
         return a == b;
 
     /// float vs float
-    if constexpr (std::is_floating_point_v<A> && std::is_floating_point_v<B>)
+    if constexpr (is_floating_point<A> && is_floating_point<B>)
         return a == b;
 
     /// anything vs NaN
@@ -125,7 +125,7 @@ bool equalsOp(A a, B b)
     }
 
     /// int vs float
-    if constexpr (is_integer<A> && std::is_floating_point_v<B>)
+    if constexpr (is_integer<A> && is_floating_point<B>)
     {
         if constexpr (sizeof(A) <= 4)
             return static_cast<double>(a) == static_cast<double>(b);
@@ -133,7 +133,7 @@ bool equalsOp(A a, B b)
         return DecomposedFloat<B>(b).equals(a);
     }
 
-    if constexpr (std::is_floating_point_v<A> && is_integer<B>)
+    if constexpr (is_floating_point<A> && is_integer<B>)
     {
         if constexpr (sizeof(B) <= 4)
             return static_cast<double>(a) == static_cast<double>(b);
@@ -163,7 +163,7 @@ inline bool NO_SANITIZE_UNDEFINED convertNumeric(From value, To & result)
         return true;
     }
 
-    if constexpr (std::is_floating_point_v<From> && std::is_floating_point_v<To>)
+    if constexpr (is_floating_point<From> && is_floating_point<To>)
     {
         /// Note that NaNs doesn't compare equal to anything, but they are still in range of any Float type.
         if (isNaN(value))
diff --git a/src/Core/DecimalFunctions.h b/src/Core/DecimalFunctions.h
index 8dad00c3a1e..abd660a8a7f 100644
--- a/src/Core/DecimalFunctions.h
+++ b/src/Core/DecimalFunctions.h
@@ -17,6 +17,7 @@ class DataTypeNumber;
 
 namespace ErrorCodes
 {
+    extern const int NOT_IMPLEMENTED;
     extern const int DECIMAL_OVERFLOW;
     extern const int ARGUMENT_OUT_OF_BOUND;
 }
@@ -310,7 +311,14 @@ ReturnType convertToImpl(const DecimalType & decimal, UInt32 scale, To & result)
     using DecimalNativeType = typename DecimalType::NativeType;
     static constexpr bool throw_exception = std::is_void_v<ReturnType>;
 
-    if constexpr (std::is_floating_point_v<To>)
+    if constexpr (std::is_same_v<To, BFloat16>)
+    {
+        if constexpr (throw_exception)
+            throw Exception(ErrorCodes::NOT_IMPLEMENTED, "Conversion from Decimal to BFloat16 is not implemented");
+        else
+            return ReturnType(false);
+    }
+    else if constexpr (is_floating_point<To>)
     {
         result = static_cast<To>(decimal.value) / static_cast<To>(scaleMultiplier<DecimalNativeType>(scale));
     }
diff --git a/src/Core/Field.h b/src/Core/Field.h
index 7b916d30646..c08d5c9eb42 100644
--- a/src/Core/Field.h
+++ b/src/Core/Field.h
@@ -257,6 +257,7 @@ template <> struct NearestFieldTypeImpl<DecimalField<Decimal64>> { using Type =
 template <> struct NearestFieldTypeImpl<DecimalField<Decimal128>> { using Type = DecimalField<Decimal128>; };
 template <> struct NearestFieldTypeImpl<DecimalField<Decimal256>> { using Type = DecimalField<Decimal256>; };
 template <> struct NearestFieldTypeImpl<DecimalField<DateTime64>> { using Type = DecimalField<DateTime64>; };
+template <> struct NearestFieldTypeImpl<BFloat16> { using Type = Float64; };
 template <> struct NearestFieldTypeImpl<Float32> { using Type = Float64; };
 template <> struct NearestFieldTypeImpl<Float64> { using Type = Float64; };
 template <> struct NearestFieldTypeImpl<const char *> { using Type = String; };
diff --git a/src/Core/Settings.cpp b/src/Core/Settings.cpp
index bdfec563397..fd418a4b97a 100644
--- a/src/Core/Settings.cpp
+++ b/src/Core/Settings.cpp
@@ -5742,7 +5742,10 @@ Enable experimental functions for natural language processing.
 Enable experimental hash functions
 )", EXPERIMENTAL) \
     DECLARE(Bool, allow_experimental_object_type, false, R"(
-Allow Object and JSON data types
+Allow the obsolete Object data type
+)", EXPERIMENTAL) \
+    DECLARE(Bool, allow_experimental_bfloat16_type, false, R"(
+Allow BFloat16 data type (under development).
 )", EXPERIMENTAL) \
     DECLARE(Bool, allow_experimental_time_series_table, false, R"(
 Allows creation of tables with the [TimeSeries](../../engines/table-engines/integrations/time-series.md) table engine.
diff --git a/src/Core/SettingsChangesHistory.cpp b/src/Core/SettingsChangesHistory.cpp
index f4a8fd198c5..455858fe241 100644
--- a/src/Core/SettingsChangesHistory.cpp
+++ b/src/Core/SettingsChangesHistory.cpp
@@ -78,6 +78,7 @@ static std::initializer_list<std::pair<ClickHouseVersion, SettingsChangesHistory
             {"backup_restore_keeper_max_retries_while_handling_error", 0, 20, "New setting."},
             {"backup_restore_finish_timeout_after_error_sec", 0, 180, "New setting."},
             {"parallel_replicas_local_plan", false, true, "Use local plan for local replica in a query with parallel replicas"},
+            {"allow_experimental_bfloat16_type", false, false, "Add new experimental BFloat16 type"},
             {"filesystem_cache_skip_download_if_exceeds_per_query_cache_write_limit", 1, 1, "Rename of setting skip_download_if_exceeds_query_cache_limit"},
             {"filesystem_cache_prefer_bigger_buffer_size", true, true, "New setting"},
             {"read_in_order_use_virtual_row", false, false, "Use virtual row while reading in order of primary key or its monotonic function fashion. It is useful when searching over multiple parts as only relevant ones are touched."},
diff --git a/src/Core/SortCursor.h b/src/Core/SortCursor.h
index 6eb009fa259..3acabaae44e 100644
--- a/src/Core/SortCursor.h
+++ b/src/Core/SortCursor.h
@@ -726,6 +726,7 @@ private:
         SortingQueueImpl<SpecializedSingleColumnSortCursor<ColumnVector<Int128>>, strategy>,
         SortingQueueImpl<SpecializedSingleColumnSortCursor<ColumnVector<Int256>>, strategy>,
 
+        SortingQueueImpl<SpecializedSingleColumnSortCursor<ColumnVector<BFloat16>>, strategy>,
         SortingQueueImpl<SpecializedSingleColumnSortCursor<ColumnVector<Float32>>, strategy>,
         SortingQueueImpl<SpecializedSingleColumnSortCursor<ColumnVector<Float64>>, strategy>,
 
diff --git a/src/Core/TypeId.h b/src/Core/TypeId.h
index 1eba944e63e..f8991a08d3c 100644
--- a/src/Core/TypeId.h
+++ b/src/Core/TypeId.h
@@ -21,6 +21,7 @@ enum class TypeIndex : uint8_t
     Int64,
     Int128,
     Int256,
+    BFloat16,
     Float32,
     Float64,
     Date,
@@ -94,6 +95,7 @@ TYPEID_MAP(Int32)
 TYPEID_MAP(Int64)
 TYPEID_MAP(Int128)
 TYPEID_MAP(Int256)
+TYPEID_MAP(BFloat16)
 TYPEID_MAP(Float32)
 TYPEID_MAP(Float64)
 TYPEID_MAP(UUID)
diff --git a/src/Core/Types_fwd.h b/src/Core/Types_fwd.h
index d3b1d168edf..b94a29ce72c 100644
--- a/src/Core/Types_fwd.h
+++ b/src/Core/Types_fwd.h
@@ -21,6 +21,7 @@ using Int128 = wide::integer<128, signed>;
 using UInt128 = wide::integer<128, unsigned>;
 using Int256 = wide::integer<256, signed>;
 using UInt256 = wide::integer<256, unsigned>;
+class BFloat16;
 
 namespace DB
 {
diff --git a/src/Core/callOnTypeIndex.h b/src/Core/callOnTypeIndex.h
index ae5afce36be..0c8f2201b0d 100644
--- a/src/Core/callOnTypeIndex.h
+++ b/src/Core/callOnTypeIndex.h
@@ -63,6 +63,7 @@ static bool callOnBasicType(TypeIndex number, F && f)
     {
         switch (number)
         {
+            case TypeIndex::BFloat16:     return f(TypePair<T, BFloat16>());
             case TypeIndex::Float32:      return f(TypePair<T, Float32>());
             case TypeIndex::Float64:      return f(TypePair<T, Float64>());
             default:
@@ -133,6 +134,7 @@ static inline bool callOnBasicTypes(TypeIndex type_num1, TypeIndex type_num2, F
     {
         switch (type_num1)
         {
+            case TypeIndex::BFloat16: return callOnBasicType<BFloat16, _int, _float, _decimal, _datetime>(type_num2, std::forward<F>(f));
             case TypeIndex::Float32: return callOnBasicType<Float32, _int, _float, _decimal, _datetime>(type_num2, std::forward<F>(f));
             case TypeIndex::Float64: return callOnBasicType<Float64, _int, _float, _decimal, _datetime>(type_num2, std::forward<F>(f));
             default:
@@ -190,6 +192,7 @@ static bool callOnIndexAndDataType(TypeIndex number, F && f, ExtraArgs && ... ar
         case TypeIndex::Int128:         return f(TypePair<DataTypeNumber<Int128>, T>(), std::forward<ExtraArgs>(args)...);
         case TypeIndex::Int256:         return f(TypePair<DataTypeNumber<Int256>, T>(), std::forward<ExtraArgs>(args)...);
 
+        case TypeIndex::BFloat16:        return f(TypePair<DataTypeNumber<BFloat16>, T>(), std::forward<ExtraArgs>(args)...);
         case TypeIndex::Float32:        return f(TypePair<DataTypeNumber<Float32>, T>(), std::forward<ExtraArgs>(args)...);
         case TypeIndex::Float64:        return f(TypePair<DataTypeNumber<Float64>, T>(), std::forward<ExtraArgs>(args)...);
 
diff --git a/src/DataTypes/DataTypeNumberBase.cpp b/src/DataTypes/DataTypeNumberBase.cpp
index be448fe1491..636d557f4d0 100644
--- a/src/DataTypes/DataTypeNumberBase.cpp
+++ b/src/DataTypes/DataTypeNumberBase.cpp
@@ -42,6 +42,7 @@ template class DataTypeNumberBase<Int32>;
 template class DataTypeNumberBase<Int64>;
 template class DataTypeNumberBase<Int128>;
 template class DataTypeNumberBase<Int256>;
+template class DataTypeNumberBase<BFloat16>;
 template class DataTypeNumberBase<Float32>;
 template class DataTypeNumberBase<Float64>;
 
diff --git a/src/DataTypes/DataTypeNumberBase.h b/src/DataTypes/DataTypeNumberBase.h
index 3a5b11c5124..11b9427a14d 100644
--- a/src/DataTypes/DataTypeNumberBase.h
+++ b/src/DataTypes/DataTypeNumberBase.h
@@ -68,6 +68,7 @@ extern template class DataTypeNumberBase<Int32>;
 extern template class DataTypeNumberBase<Int64>;
 extern template class DataTypeNumberBase<Int128>;
 extern template class DataTypeNumberBase<Int256>;
+extern template class DataTypeNumberBase<BFloat16>;
 extern template class DataTypeNumberBase<Float32>;
 extern template class DataTypeNumberBase<Float64>;
 
diff --git a/src/DataTypes/DataTypesBinaryEncoding.cpp b/src/DataTypes/DataTypesBinaryEncoding.cpp
index dc0f2f3f5aa..c3190b462c3 100644
--- a/src/DataTypes/DataTypesBinaryEncoding.cpp
+++ b/src/DataTypes/DataTypesBinaryEncoding.cpp
@@ -96,6 +96,7 @@ enum class BinaryTypeIndex : uint8_t
     SimpleAggregateFunction = 0x2E,
     Nested = 0x2F,
     JSON = 0x30,
+    BFloat16 = 0x31,
 };
 
 /// In future we can introduce more arguments in the JSON data type definition.
@@ -151,6 +152,8 @@ BinaryTypeIndex getBinaryTypeIndex(const DataTypePtr & type)
             return BinaryTypeIndex::Int128;
         case TypeIndex::Int256:
             return BinaryTypeIndex::Int256;
+        case TypeIndex::BFloat16:
+            return BinaryTypeIndex::BFloat16;
         case TypeIndex::Float32:
             return BinaryTypeIndex::Float32;
         case TypeIndex::Float64:
@@ -565,6 +568,8 @@ DataTypePtr decodeDataType(ReadBuffer & buf)
             return std::make_shared<DataTypeInt128>();
         case BinaryTypeIndex::Int256:
             return std::make_shared<DataTypeInt256>();
+        case BinaryTypeIndex::BFloat16:
+            return std::make_shared<DataTypeBFloat16>();
         case BinaryTypeIndex::Float32:
             return std::make_shared<DataTypeFloat32>();
         case BinaryTypeIndex::Float64:
diff --git a/src/DataTypes/DataTypesDecimal.cpp b/src/DataTypes/DataTypesDecimal.cpp
index 038a00932ee..63bd4bf2a59 100644
--- a/src/DataTypes/DataTypesDecimal.cpp
+++ b/src/DataTypes/DataTypesDecimal.cpp
@@ -2,6 +2,7 @@
 #include <DataTypes/Serializations/SerializationDecimal.h>
 
 #include <Common/typeid_cast.h>
+#include <Common/NaNUtils.h>
 #include <Core/DecimalFunctions.h>
 #include <DataTypes/DataTypeFactory.h>
 #include <IO/ReadHelpers.h>
@@ -19,6 +20,7 @@ namespace ErrorCodes
     extern const int NUMBER_OF_ARGUMENTS_DOESNT_MATCH;
     extern const int ILLEGAL_TYPE_OF_ARGUMENT;
     extern const int DECIMAL_OVERFLOW;
+    extern const int NOT_IMPLEMENTED;
 }
 
 
@@ -268,9 +270,13 @@ ReturnType convertToDecimalImpl(const typename FromDataType::FieldType & value,
 
     static constexpr bool throw_exception = std::is_same_v<ReturnType, void>;
 
-    if constexpr (std::is_floating_point_v<FromFieldType>)
+    if constexpr (std::is_same_v<typename FromDataType::FieldType, BFloat16>)
     {
-        if (!std::isfinite(value))
+        throw Exception(ErrorCodes::NOT_IMPLEMENTED, "Conversion from BFloat16 to Decimal is not implemented");
+    }
+    else if constexpr (is_floating_point<FromFieldType>)
+    {
+        if (!isFinite(value))
         {
             if constexpr (throw_exception)
                 throw Exception(ErrorCodes::DECIMAL_OVERFLOW, "{} convert overflow. Cannot convert infinity or NaN to decimal", ToDataType::family_name);
diff --git a/src/DataTypes/DataTypesDecimal.h b/src/DataTypes/DataTypesDecimal.h
index badefc4c75a..09a25617506 100644
--- a/src/DataTypes/DataTypesDecimal.h
+++ b/src/DataTypes/DataTypesDecimal.h
@@ -4,7 +4,6 @@
 #include <base/extended_types.h>
 #include <Common/typeid_cast.h>
 #include <base/Decimal.h>
-#include <base/Decimal_fwd.h>
 #include <DataTypes/IDataType.h>
 #include <DataTypes/DataTypeDate.h>
 #include <DataTypes/DataTypeDate32.h>
@@ -205,7 +204,6 @@ FOR_EACH_DECIMAL_TYPE(INVOKE);
 #undef INVOKE
 #undef DISPATCH
 
-
 template <typename FromDataType, typename ToDataType>
 requires (is_arithmetic_v<typename FromDataType::FieldType> && IsDataTypeDecimal<ToDataType>)
 typename ToDataType::FieldType convertToDecimal(const typename FromDataType::FieldType & value, UInt32 scale);
diff --git a/src/DataTypes/DataTypesNumber.cpp b/src/DataTypes/DataTypesNumber.cpp
index 72020b0a5aa..4c8918521fe 100644
--- a/src/DataTypes/DataTypesNumber.cpp
+++ b/src/DataTypes/DataTypesNumber.cpp
@@ -54,6 +54,7 @@ void registerDataTypeNumbers(DataTypeFactory & factory)
     factory.registerDataType("Int32", createNumericDataType<Int32>);
     factory.registerDataType("Int64", createNumericDataType<Int64>);
 
+    factory.registerDataType("BFloat16", createNumericDataType<BFloat16>);
     factory.registerDataType("Float32", createNumericDataType<Float32>);
     factory.registerDataType("Float64", createNumericDataType<Float64>);
 
@@ -111,6 +112,7 @@ template class DataTypeNumber<Int8>;
 template class DataTypeNumber<Int16>;
 template class DataTypeNumber<Int32>;
 template class DataTypeNumber<Int64>;
+template class DataTypeNumber<BFloat16>;
 template class DataTypeNumber<Float32>;
 template class DataTypeNumber<Float64>;
 
diff --git a/src/DataTypes/DataTypesNumber.h b/src/DataTypes/DataTypesNumber.h
index d550ceababc..a9e77e01b13 100644
--- a/src/DataTypes/DataTypesNumber.h
+++ b/src/DataTypes/DataTypesNumber.h
@@ -63,6 +63,7 @@ extern template class DataTypeNumber<Int8>;
 extern template class DataTypeNumber<Int16>;
 extern template class DataTypeNumber<Int32>;
 extern template class DataTypeNumber<Int64>;
+extern template class DataTypeNumber<BFloat16>;
 extern template class DataTypeNumber<Float32>;
 extern template class DataTypeNumber<Float64>;
 
@@ -79,6 +80,7 @@ using DataTypeInt8 = DataTypeNumber<Int8>;
 using DataTypeInt16 = DataTypeNumber<Int16>;
 using DataTypeInt32 = DataTypeNumber<Int32>;
 using DataTypeInt64 = DataTypeNumber<Int64>;
+using DataTypeBFloat16 = DataTypeNumber<BFloat16>;
 using DataTypeFloat32 = DataTypeNumber<Float32>;
 using DataTypeFloat64 = DataTypeNumber<Float64>;
 
diff --git a/src/DataTypes/IDataType.h b/src/DataTypes/IDataType.h
index 33eddf8e9b8..8f06526ddbb 100644
--- a/src/DataTypes/IDataType.h
+++ b/src/DataTypes/IDataType.h
@@ -408,9 +408,11 @@ struct WhichDataType
     constexpr bool isDecimal256() const { return idx == TypeIndex::Decimal256; }
     constexpr bool isDecimal() const { return isDecimal32() || isDecimal64() || isDecimal128() || isDecimal256(); }
 
+    constexpr bool isBFloat16() const { return idx == TypeIndex::BFloat16; }
     constexpr bool isFloat32() const { return idx == TypeIndex::Float32; }
     constexpr bool isFloat64() const { return idx == TypeIndex::Float64; }
-    constexpr bool isFloat() const { return isFloat32() || isFloat64(); }
+    constexpr bool isNativeFloat() const { return isFloat32() || isFloat64(); }
+    constexpr bool isFloat() const { return isNativeFloat() || isBFloat16(); }
 
     constexpr bool isNativeNumber() const { return isNativeInteger() || isFloat(); }
     constexpr bool isNumber() const { return isInteger() || isFloat() || isDecimal(); }
@@ -621,6 +623,7 @@ template <typename T> inline constexpr bool IsDataTypeEnum<DataTypeEnum<T>> = tr
     M(Int64) \
     M(Int128) \
     M(Int256) \
+    M(BFloat16) \
     M(Float32) \
     M(Float64)
 }
diff --git a/src/DataTypes/Native.cpp b/src/DataTypes/Native.cpp
index 5dc490b0bd5..53354d7c6e0 100644
--- a/src/DataTypes/Native.cpp
+++ b/src/DataTypes/Native.cpp
@@ -37,7 +37,7 @@ bool canBeNativeType(const IDataType & type)
         return canBeNativeType(*data_type_nullable.getNestedType());
     }
 
-    return data_type.isNativeInt() || data_type.isNativeUInt() || data_type.isFloat() || data_type.isDate()
+    return data_type.isNativeInt() || data_type.isNativeUInt() || data_type.isNativeFloat() || data_type.isDate()
         || data_type.isDate32() || data_type.isDateTime() || data_type.isEnum();
 }
 
diff --git a/src/DataTypes/NumberTraits.h b/src/DataTypes/NumberTraits.h
index 59a64017af3..ee0d9812097 100644
--- a/src/DataTypes/NumberTraits.h
+++ b/src/DataTypes/NumberTraits.h
@@ -74,7 +74,7 @@ template <typename A, typename B> struct ResultOfAdditionMultiplication
 {
     using Type = typename Construct<
         is_signed_v<A> || is_signed_v<B>,
-        std::is_floating_point_v<A> || std::is_floating_point_v<B>,
+        is_floating_point<A> || is_floating_point<B>,
         nextSize(max(sizeof(A), sizeof(B)))>::Type;
 };
 
@@ -82,7 +82,7 @@ template <typename A, typename B> struct ResultOfSubtraction
 {
     using Type = typename Construct<
         true,
-        std::is_floating_point_v<A> || std::is_floating_point_v<B>,
+        is_floating_point<A> || is_floating_point<B>,
         nextSize(max(sizeof(A), sizeof(B)))>::Type;
 };
 
@@ -113,7 +113,7 @@ template <typename A, typename B> struct ResultOfModulo
     /// Example: toInt32(-199) % toUInt8(200) will return -199 that does not fit in Int8, only in Int16.
     static constexpr size_t size_of_result = result_is_signed ? nextSize(sizeof(B)) : sizeof(B);
     using Type0 = typename Construct<result_is_signed, false, size_of_result>::Type;
-    using Type = std::conditional_t<std::is_floating_point_v<A> || std::is_floating_point_v<B>, Float64, Type0>;
+    using Type = std::conditional_t<is_floating_point<A> || is_floating_point<B>, Float64, Type0>;
 };
 
 template <typename A, typename B> struct ResultOfPositiveModulo
@@ -121,21 +121,21 @@ template <typename A, typename B> struct ResultOfPositiveModulo
     /// function positive_modulo always return non-negative number.
     static constexpr size_t size_of_result = sizeof(B);
     using Type0 = typename Construct<false, false, size_of_result>::Type;
-    using Type = std::conditional_t<std::is_floating_point_v<A> || std::is_floating_point_v<B>, Float64, Type0>;
+    using Type = std::conditional_t<is_floating_point<A> || is_floating_point<B>, Float64, Type0>;
 };
 
 
 template <typename A, typename B> struct ResultOfModuloLegacy
 {
     using Type0 = typename Construct<is_signed_v<A> || is_signed_v<B>, false, sizeof(B)>::Type;
-    using Type = std::conditional_t<std::is_floating_point_v<A> || std::is_floating_point_v<B>, Float64, Type0>;
+    using Type = std::conditional_t<is_floating_point<A> || is_floating_point<B>, Float64, Type0>;
 };
 
 template <typename A> struct ResultOfNegate
 {
     using Type = typename Construct<
         true,
-        std::is_floating_point_v<A>,
+        is_floating_point<A>,
         is_signed_v<A> ? sizeof(A) : nextSize(sizeof(A))>::Type;
 };
 
@@ -143,7 +143,7 @@ template <typename A> struct ResultOfAbs
 {
     using Type = typename Construct<
         false,
-        std::is_floating_point_v<A>,
+        is_floating_point<A>,
         sizeof(A)>::Type;
 };
 
@@ -154,7 +154,7 @@ template <typename A, typename B> struct ResultOfBit
     using Type = typename Construct<
         is_signed_v<A> || is_signed_v<B>,
         false,
-        std::is_floating_point_v<A> || std::is_floating_point_v<B> ? 8 : max(sizeof(A), sizeof(B))>::Type;
+        is_floating_point<A> || is_floating_point<B> ? 8 : max(sizeof(A), sizeof(B))>::Type;
 };
 
 template <typename A> struct ResultOfBitNot
@@ -180,7 +180,7 @@ template <typename A> struct ResultOfBitNot
 template <typename A, typename B>
 struct ResultOfIf
 {
-    static constexpr bool has_float = std::is_floating_point_v<A> || std::is_floating_point_v<B>;
+    static constexpr bool has_float = is_floating_point<A> || is_floating_point<B>;
     static constexpr bool has_integer = is_integer<A> || is_integer<B>;
     static constexpr bool has_signed = is_signed_v<A> || is_signed_v<B>;
     static constexpr bool has_unsigned = !is_signed_v<A> || !is_signed_v<B>;
@@ -189,7 +189,7 @@ struct ResultOfIf
     static constexpr size_t max_size_of_unsigned_integer = max(is_signed_v<A> ? 0 : sizeof(A), is_signed_v<B> ? 0 : sizeof(B));
     static constexpr size_t max_size_of_signed_integer = max(is_signed_v<A> ? sizeof(A) : 0, is_signed_v<B> ? sizeof(B) : 0);
     static constexpr size_t max_size_of_integer = max(is_integer<A> ? sizeof(A) : 0, is_integer<B> ? sizeof(B) : 0);
-    static constexpr size_t max_size_of_float = max(std::is_floating_point_v<A> ? sizeof(A) : 0, std::is_floating_point_v<B> ? sizeof(B) : 0);
+    static constexpr size_t max_size_of_float = max(is_floating_point<A> ? sizeof(A) : 0, is_floating_point<B> ? sizeof(B) : 0);
 
     using ConstructedType = typename Construct<has_signed, has_float,
         ((has_float && has_integer && max_size_of_integer >= max_size_of_float)
@@ -211,7 +211,7 @@ template <typename A> struct ToInteger
     using Type = typename Construct<
         is_signed_v<A>,
         false,
-        std::is_floating_point_v<A> ? 8 : sizeof(A)>::Type;
+        is_floating_point<A> ? 8 : sizeof(A)>::Type;
 };
 
 
diff --git a/src/DataTypes/Serializations/SerializationNumber.cpp b/src/DataTypes/Serializations/SerializationNumber.cpp
index bfc13af8ca3..a20b97daf07 100644
--- a/src/DataTypes/Serializations/SerializationNumber.cpp
+++ b/src/DataTypes/Serializations/SerializationNumber.cpp
@@ -238,6 +238,7 @@ template class SerializationNumber<Int32>;
 template class SerializationNumber<Int64>;
 template class SerializationNumber<Int128>;
 template class SerializationNumber<Int256>;
+template class SerializationNumber<BFloat16>;
 template class SerializationNumber<Float32>;
 template class SerializationNumber<Float64>;
 
diff --git a/src/DataTypes/Utils.cpp b/src/DataTypes/Utils.cpp
index a6e9452d7ef..3004322475c 100644
--- a/src/DataTypes/Utils.cpp
+++ b/src/DataTypes/Utils.cpp
@@ -54,6 +54,13 @@ bool canBeSafelyCasted(const DataTypePtr & from_type, const DataTypePtr & to_typ
 
             return false;
         }
+        case TypeIndex::BFloat16:
+        {
+            if (to_which_type.isFloat32() || to_which_type.isFloat64() || to_which_type.isString())
+                return true;
+
+            return false;
+        }
         case TypeIndex::Float32:
         {
             if (to_which_type.isFloat64() || to_which_type.isString())
diff --git a/src/DataTypes/getLeastSupertype.cpp b/src/DataTypes/getLeastSupertype.cpp
index 9664e769611..1b258a8decc 100644
--- a/src/DataTypes/getLeastSupertype.cpp
+++ b/src/DataTypes/getLeastSupertype.cpp
@@ -109,6 +109,8 @@ DataTypePtr getNumericType(const TypeIndexSet & types)
             maximize(max_bits_of_signed_integer, 128);
         else if (type == TypeIndex::Int256)
             maximize(max_bits_of_signed_integer, 256);
+        else if (type == TypeIndex::BFloat16)
+            maximize(max_mantissa_bits_of_floating, 8);
         else if (type == TypeIndex::Float32)
             maximize(max_mantissa_bits_of_floating, 24);
         else if (type == TypeIndex::Float64)
@@ -145,7 +147,9 @@ DataTypePtr getNumericType(const TypeIndexSet & types)
         if (max_mantissa_bits_of_floating)
         {
             size_t min_mantissa_bits = std::max(min_bit_width_of_integer, max_mantissa_bits_of_floating);
-            if (min_mantissa_bits <= 24)
+            if (min_mantissa_bits <= 8)
+                return std::make_shared<DataTypeBFloat16>();
+            else if (min_mantissa_bits <= 24)
                 return std::make_shared<DataTypeFloat32>();
             if (min_mantissa_bits <= 53)
                 return std::make_shared<DataTypeFloat64>();
diff --git a/src/DataTypes/getMostSubtype.cpp b/src/DataTypes/getMostSubtype.cpp
index b3c76079646..8cb16c04079 100644
--- a/src/DataTypes/getMostSubtype.cpp
+++ b/src/DataTypes/getMostSubtype.cpp
@@ -297,6 +297,8 @@ DataTypePtr getMostSubtype(const DataTypes & types, bool throw_if_result_is_noth
                 minimize(min_bits_of_signed_integer, 128);
             else if (typeid_cast<const DataTypeInt256 *>(type.get()))
                 minimize(min_bits_of_signed_integer, 256);
+            else if (typeid_cast<const DataTypeBFloat16 *>(type.get()))
+                minimize(min_mantissa_bits_of_floating, 8);
             else if (typeid_cast<const DataTypeFloat32 *>(type.get()))
                 minimize(min_mantissa_bits_of_floating, 24);
             else if (typeid_cast<const DataTypeFloat64 *>(type.get()))
@@ -313,7 +315,9 @@ DataTypePtr getMostSubtype(const DataTypes & types, bool throw_if_result_is_noth
             /// If the result must be floating.
             if (!min_bits_of_signed_integer && !min_bits_of_unsigned_integer)
             {
-                if (min_mantissa_bits_of_floating <= 24)
+                if (min_mantissa_bits_of_floating <= 8)
+                    return std::make_shared<DataTypeBFloat16>();
+                else if (min_mantissa_bits_of_floating <= 24)
                     return std::make_shared<DataTypeFloat32>();
                 if (min_mantissa_bits_of_floating <= 53)
                     return std::make_shared<DataTypeFloat64>();
diff --git a/src/Databases/enableAllExperimentalSettings.cpp b/src/Databases/enableAllExperimentalSettings.cpp
index d51d2671992..1be54664bc9 100644
--- a/src/Databases/enableAllExperimentalSettings.cpp
+++ b/src/Databases/enableAllExperimentalSettings.cpp
@@ -24,10 +24,11 @@ void enableAllExperimentalSettings(ContextMutablePtr context)
     context->setSetting("allow_experimental_dynamic_type", 1);
     context->setSetting("allow_experimental_json_type", 1);
     context->setSetting("allow_experimental_vector_similarity_index", 1);
-    context->setSetting("allow_experimental_bigint_types", 1);
     context->setSetting("allow_experimental_window_functions", 1);
     context->setSetting("allow_experimental_geo_types", 1);
     context->setSetting("allow_experimental_map_type", 1);
+    context->setSetting("allow_experimental_bigint_types", 1);
+    context->setSetting("allow_experimental_bfloat16_type", 1);
     context->setSetting("allow_deprecated_error_prone_window_functions", 1);
 
     context->setSetting("allow_suspicious_low_cardinality_types", 1);
diff --git a/src/Dictionaries/RangeHashedDictionary.h b/src/Dictionaries/RangeHashedDictionary.h
index bc335b890fc..5e93391816f 100644
--- a/src/Dictionaries/RangeHashedDictionary.h
+++ b/src/Dictionaries/RangeHashedDictionary.h
@@ -298,7 +298,8 @@ namespace impl
             using Types = std::decay_t<decltype(types)>;
             using DataType = typename Types::LeftType;
 
-            if constexpr (IsDataTypeDecimalOrNumber<DataType> || IsDataTypeDateOrDateTime<DataType> || IsDataTypeEnum<DataType>)
+            if constexpr ((IsDataTypeDecimalOrNumber<DataType> || IsDataTypeDateOrDateTime<DataType> || IsDataTypeEnum<DataType>)
+                && !std::is_same_v<DataType, DataTypeBFloat16>)
             {
                 using ColumnType = typename DataType::ColumnType;
                 func(TypePair<ColumnType, void>());
diff --git a/src/Formats/JSONExtractTree.cpp b/src/Formats/JSONExtractTree.cpp
index ae6051823b7..62905a2e630 100644
--- a/src/Formats/JSONExtractTree.cpp
+++ b/src/Formats/JSONExtractTree.cpp
@@ -131,7 +131,7 @@ bool tryGetNumericValueFromJSONElement(
     switch (element.type())
     {
         case ElementType::DOUBLE:
-            if constexpr (std::is_floating_point_v<NumberType>)
+            if constexpr (is_floating_point<NumberType>)
             {
                 /// We permit inaccurate conversion of double to float.
                 /// Example: double 0.1 from JSON is not representable in float.
@@ -175,7 +175,7 @@ bool tryGetNumericValueFromJSONElement(
                 return false;
 
             auto rb = ReadBufferFromMemory{element.getString()};
-            if constexpr (std::is_floating_point_v<NumberType>)
+            if constexpr (is_floating_point<NumberType>)
             {
                 if (!tryReadFloatText(value, rb) || !rb.eof())
                 {
diff --git a/src/Formats/ProtobufSerializer.cpp b/src/Formats/ProtobufSerializer.cpp
index 2224af6579f..bd8aeff5a5d 100644
--- a/src/Formats/ProtobufSerializer.cpp
+++ b/src/Formats/ProtobufSerializer.cpp
@@ -540,7 +540,7 @@ namespace
 
                 case FieldTypeId::TYPE_ENUM:
                 {
-                    if (std::is_floating_point_v<NumberType>)
+                    if (is_floating_point<NumberType>)
                         incompatibleColumnType(TypeName<NumberType>);
 
                     write_function = [this](NumberType value)
diff --git a/src/Functions/DivisionUtils.h b/src/Functions/DivisionUtils.h
index 7fd5b7476e1..e8f5da342f8 100644
--- a/src/Functions/DivisionUtils.h
+++ b/src/Functions/DivisionUtils.h
@@ -47,9 +47,9 @@ inline auto checkedDivision(A a, B b)
 {
     throwIfDivisionLeadsToFPE(a, b);
 
-    if constexpr (is_big_int_v<A> && std::is_floating_point_v<B>)
+    if constexpr (is_big_int_v<A> && is_floating_point<B>)
         return static_cast<B>(a) / b;
-    else if constexpr (is_big_int_v<B> && std::is_floating_point_v<A>)
+    else if constexpr (is_big_int_v<B> && is_floating_point<A>)
         return a / static_cast<A>(b);
     else if constexpr (is_big_int_v<A> && is_big_int_v<B>)
         return static_cast<A>(a / b);
@@ -86,17 +86,17 @@ struct DivideIntegralImpl
         {
             /// Comparisons are not strict to avoid rounding issues when operand is implicitly casted to float.
 
-            if constexpr (std::is_floating_point_v<A>)
+            if constexpr (is_floating_point<A>)
                 if (isNaN(a) || a >= std::numeric_limits<CastA>::max() || a <= std::numeric_limits<CastA>::lowest())
                     throw Exception(ErrorCodes::ILLEGAL_DIVISION, "Cannot perform integer division on infinite or too large floating point numbers");
 
-            if constexpr (std::is_floating_point_v<B>)
+            if constexpr (is_floating_point<B>)
                 if (isNaN(b) || b >= std::numeric_limits<CastB>::max() || b <= std::numeric_limits<CastB>::lowest())
                     throw Exception(ErrorCodes::ILLEGAL_DIVISION, "Cannot perform integer division on infinite or too large floating point numbers");
 
             auto res = checkedDivision(CastA(a), CastB(b));
 
-            if constexpr (std::is_floating_point_v<decltype(res)>)
+            if constexpr (is_floating_point<decltype(res)>)
                 if (isNaN(res) || res >= static_cast<double>(std::numeric_limits<Result>::max()) || res <= std::numeric_limits<Result>::lowest())
                     throw Exception(ErrorCodes::ILLEGAL_DIVISION, "Cannot perform integer division, because it will produce infinite or too large number");
 
@@ -122,18 +122,18 @@ struct ModuloImpl
     template <typename Result = ResultType>
     static Result apply(A a, B b)
     {
-        if constexpr (std::is_floating_point_v<ResultType>)
+        if constexpr (is_floating_point<ResultType>)
         {
             /// This computation is similar to `fmod` but the latter is not inlined and has 40 times worse performance.
             return static_cast<ResultType>(a) - trunc(static_cast<ResultType>(a) / static_cast<ResultType>(b)) * static_cast<ResultType>(b);
         }
         else
         {
-            if constexpr (std::is_floating_point_v<A>)
+            if constexpr (is_floating_point<A>)
                 if (isNaN(a) || a > std::numeric_limits<IntegerAType>::max() || a < std::numeric_limits<IntegerAType>::lowest())
                     throw Exception(ErrorCodes::ILLEGAL_DIVISION, "Cannot perform integer division on infinite or too large floating point numbers");
 
-            if constexpr (std::is_floating_point_v<B>)
+            if constexpr (is_floating_point<B>)
                 if (isNaN(b) || b > std::numeric_limits<IntegerBType>::max() || b < std::numeric_limits<IntegerBType>::lowest())
                     throw Exception(ErrorCodes::ILLEGAL_DIVISION, "Cannot perform integer division on infinite or too large floating point numbers");
 
diff --git a/src/Functions/FunctionBinaryArithmetic.h b/src/Functions/FunctionBinaryArithmetic.h
index df239b820af..e7f0d8d43b0 100644
--- a/src/Functions/FunctionBinaryArithmetic.h
+++ b/src/Functions/FunctionBinaryArithmetic.h
@@ -110,6 +110,7 @@ template <typename DataType> constexpr bool IsIntegralOrExtendedOrDecimal =
     IsDataTypeDecimal<DataType>;
 
 template <typename DataType> constexpr bool IsFloatingPoint = false;
+template <> inline constexpr bool IsFloatingPoint<DataTypeBFloat16> = true;
 template <> inline constexpr bool IsFloatingPoint<DataTypeFloat32> = true;
 template <> inline constexpr bool IsFloatingPoint<DataTypeFloat64> = true;
 
@@ -803,7 +804,7 @@ class FunctionBinaryArithmetic : public IFunction
             DataTypeFixedString, DataTypeString,
             DataTypeInterval>;
 
-        using Floats = TypeList<DataTypeFloat32, DataTypeFloat64>;
+        using Floats = TypeList<DataTypeFloat32, DataTypeFloat64, DataTypeBFloat16>;
 
         using ValidTypes = std::conditional_t<valid_on_float_arguments,
             TypeListConcat<Types, Floats>,
@@ -1690,6 +1691,13 @@ public:
             }
             else
             {
+                if constexpr ((std::is_same_v<LeftDataType, DataTypeBFloat16> || std::is_same_v<RightDataType, DataTypeBFloat16>)
+                    && (sizeof(typename LeftDataType::FieldType) > 8 || sizeof(typename RightDataType::FieldType) > 8))
+                {
+                    /// Big integers and BFloat16 are not supported together.
+                    return false;
+                }
+
                 using ResultDataType = typename BinaryOperationTraits<Op, LeftDataType, RightDataType>::ResultDataType;
 
                 if constexpr (!std::is_same_v<ResultDataType, InvalidType>)
@@ -2042,7 +2050,15 @@ ColumnPtr executeStringInteger(const ColumnsWithTypeAndName & arguments, const A
         using DecimalResultType = typename BinaryOperationTraits<Op, LeftDataType, RightDataType>::DecimalResultDataType;
 
         if constexpr (std::is_same_v<ResultDataType, InvalidType>)
+        {
             return nullptr;
+        }
+        else if constexpr ((std::is_same_v<LeftDataType, DataTypeBFloat16> || std::is_same_v<RightDataType, DataTypeBFloat16>)
+            && (sizeof(typename LeftDataType::FieldType) > 8 || sizeof(typename RightDataType::FieldType) > 8))
+        {
+            /// Big integers and BFloat16 are not supported together.
+            return nullptr;
+        }
         else // we can't avoid the else because otherwise the compiler may assume the ResultDataType may be Invalid
              // and that would produce the compile error.
         {
@@ -2059,7 +2075,7 @@ ColumnPtr executeStringInteger(const ColumnsWithTypeAndName & arguments, const A
             ColumnPtr left_col = nullptr;
             ColumnPtr right_col = nullptr;
 
-            /// When Decimal op Float32/64, convert both of them into Float64
+            /// When Decimal op Float32/64/16, convert both of them into Float64
             if constexpr (decimal_with_float)
             {
                 const auto converted_type = std::make_shared<DataTypeFloat64>();
@@ -2094,7 +2110,6 @@ ColumnPtr executeStringInteger(const ColumnsWithTypeAndName & arguments, const A
             /// Here we check if we have `intDiv` or `intDivOrZero` and at least one of the arguments is decimal, because in this case originally we had result as decimal, so we need to convert result into integer after calculations
             else if constexpr (!decimal_with_float && (is_int_div || is_int_div_or_zero) && (IsDataTypeDecimal<LeftDataType> || IsDataTypeDecimal<RightDataType>))
             {
-
                 if constexpr (!std::is_same_v<DecimalResultType, InvalidType>)
                 {
                     DataTypePtr type_res;
diff --git a/src/Functions/FunctionMathUnary.h b/src/Functions/FunctionMathUnary.h
index 4df70f4782b..12a15a072b2 100644
--- a/src/Functions/FunctionMathUnary.h
+++ b/src/Functions/FunctionMathUnary.h
@@ -70,7 +70,7 @@ private:
             /// Process all data as a whole and use FastOps implementation
 
             /// If the argument is integer, convert to Float64 beforehand
-            if constexpr (!std::is_floating_point_v<T>)
+            if constexpr (!is_floating_point<T>)
             {
                 PODArray<Float64> tmp_vec(size);
                 for (size_t i = 0; i < size; ++i)
@@ -152,7 +152,7 @@ private:
         {
             using Types = std::decay_t<decltype(types)>;
             using Type = typename Types::RightType;
-            using ReturnType = std::conditional_t<Impl::always_returns_float64 || !std::is_floating_point_v<Type>, Float64, Type>;
+            using ReturnType = std::conditional_t<Impl::always_returns_float64 || !is_floating_point<Type>, Float64, Type>;
             using ColVecType = ColumnVectorOrDecimal<Type>;
 
             const auto col_vec = checkAndGetColumn<ColVecType>(col.column.get());
diff --git a/src/Functions/FunctionsBinaryRepresentation.cpp b/src/Functions/FunctionsBinaryRepresentation.cpp
index c8e8f167e4c..50a3c0862f4 100644
--- a/src/Functions/FunctionsBinaryRepresentation.cpp
+++ b/src/Functions/FunctionsBinaryRepresentation.cpp
@@ -296,6 +296,7 @@ public:
             tryExecuteUIntOrInt<Int256>(column, res_column) ||
             tryExecuteString(column, res_column) ||
             tryExecuteFixedString(column, res_column) ||
+            tryExecuteFloat<BFloat16>(column, res_column) ||
             tryExecuteFloat<Float32>(column, res_column) ||
             tryExecuteFloat<Float64>(column, res_column) ||
             tryExecuteDecimal<Decimal32>(column, res_column) ||
diff --git a/src/Functions/FunctionsComparison.h b/src/Functions/FunctionsComparison.h
index 9b2328065fc..c9d3ccb4445 100644
--- a/src/Functions/FunctionsComparison.h
+++ b/src/Functions/FunctionsComparison.h
@@ -721,6 +721,7 @@ private:
                 || (res = executeNumRightType<T0, Int64>(col_left, col_right_untyped))
                 || (res = executeNumRightType<T0, Int128>(col_left, col_right_untyped))
                 || (res = executeNumRightType<T0, Int256>(col_left, col_right_untyped))
+                || (res = executeNumRightType<T0, BFloat16>(col_left, col_right_untyped))
                 || (res = executeNumRightType<T0, Float32>(col_left, col_right_untyped))
                 || (res = executeNumRightType<T0, Float64>(col_left, col_right_untyped)))
                 return res;
@@ -741,6 +742,7 @@ private:
                 || (res = executeNumConstRightType<T0, Int64>(col_left_const, col_right_untyped))
                 || (res = executeNumConstRightType<T0, Int128>(col_left_const, col_right_untyped))
                 || (res = executeNumConstRightType<T0, Int256>(col_left_const, col_right_untyped))
+                || (res = executeNumConstRightType<T0, BFloat16>(col_left_const, col_right_untyped))
                 || (res = executeNumConstRightType<T0, Float32>(col_left_const, col_right_untyped))
                 || (res = executeNumConstRightType<T0, Float64>(col_left_const, col_right_untyped)))
                 return res;
@@ -1292,9 +1294,10 @@ public:
                 || (res = executeNumLeftType<Int64>(col_left_untyped, col_right_untyped))
                 || (res = executeNumLeftType<Int128>(col_left_untyped, col_right_untyped))
                 || (res = executeNumLeftType<Int256>(col_left_untyped, col_right_untyped))
+                || (res = executeNumLeftType<BFloat16>(col_left_untyped, col_right_untyped))
                 || (res = executeNumLeftType<Float32>(col_left_untyped, col_right_untyped))
                 || (res = executeNumLeftType<Float64>(col_left_untyped, col_right_untyped))))
-                throw Exception(ErrorCodes::ILLEGAL_COLUMN, "Illegal column {} of first argument of function {}",
+                throw Exception(ErrorCodes::ILLEGAL_COLUMN, "Illegal column {} of the first argument of function {}",
                     col_left_untyped->getName(), getName());
 
             return res;
@@ -1342,7 +1345,7 @@ public:
                     getName(),
                     left_type->getName(),
                     right_type->getName());
-            /// When Decimal comparing to Float32/64, we convert both of them into Float64.
+            /// When Decimal comparing to Float32/64/16, we convert both of them into Float64.
             /// Other systems like MySQL and Spark also do as this.
             if (left_is_float || right_is_float)
             {
diff --git a/src/Functions/FunctionsConversion.cpp b/src/Functions/FunctionsConversion.cpp
index 1e7cb8f1749..d581ae510e8 100644
--- a/src/Functions/FunctionsConversion.cpp
+++ b/src/Functions/FunctionsConversion.cpp
@@ -7,10 +7,8 @@
 #include <Columns/ColumnFixedString.h>
 #include <Columns/ColumnLowCardinality.h>
 #include <Columns/ColumnMap.h>
-#include <Columns/ColumnNothing.h>
 #include <Columns/ColumnNullable.h>
 #include <Columns/ColumnObjectDeprecated.h>
-#include <Columns/ColumnObject.h>
 #include <Columns/ColumnString.h>
 #include <Columns/ColumnStringHelpers.h>
 #include <Columns/ColumnTuple.h>
@@ -73,8 +71,10 @@
 #include <Common/assert_cast.h>
 #include <Common/quoteString.h>
 
+
 namespace DB
 {
+
 namespace Setting
 {
     extern const SettingsBool cast_ipv4_ipv6_default_on_conversion_error;
@@ -653,7 +653,7 @@ inline void convertFromTime<DataTypeDateTime>(DataTypeDateTime::FieldType & x, t
 template <typename DataType>
 void parseImpl(typename DataType::FieldType & x, ReadBuffer & rb, const DateLUTImpl *, bool precise_float_parsing)
 {
-    if constexpr (std::is_floating_point_v<typename DataType::FieldType>)
+    if constexpr (is_floating_point<typename DataType::FieldType>)
     {
         if (precise_float_parsing)
             readFloatTextPrecise(x, rb);
@@ -717,7 +717,7 @@ inline void parseImpl<DataTypeIPv6>(DataTypeIPv6::FieldType & x, ReadBuffer & rb
 template <typename DataType>
 bool tryParseImpl(typename DataType::FieldType & x, ReadBuffer & rb, const DateLUTImpl *, bool precise_float_parsing)
 {
-    if constexpr (std::is_floating_point_v<typename DataType::FieldType>)
+    if constexpr (is_floating_point<typename DataType::FieldType>)
     {
         if (precise_float_parsing)
             return tryReadFloatTextPrecise(x, rb);
@@ -1847,7 +1847,7 @@ struct ConvertImpl
                 else
                 {
                     /// If From Data is Nan or Inf and we convert to integer type, throw exception
-                    if constexpr (std::is_floating_point_v<FromFieldType> && !std::is_floating_point_v<ToFieldType>)
+                    if constexpr (is_floating_point<FromFieldType> && !is_floating_point<ToFieldType>)
                     {
                         if (!isFinite(vec_from[i]))
                         {
@@ -2333,9 +2333,9 @@ private:
                 using RightT = typename RightDataType::FieldType;
 
                 static constexpr bool bad_left =
-                    is_decimal<LeftT> || std::is_floating_point_v<LeftT> || is_big_int_v<LeftT> || is_signed_v<LeftT>;
+                    is_decimal<LeftT> || is_floating_point<LeftT> || is_big_int_v<LeftT> || is_signed_v<LeftT>;
                 static constexpr bool bad_right =
-                    is_decimal<RightT> || std::is_floating_point_v<RightT> || is_big_int_v<RightT> || is_signed_v<RightT>;
+                    is_decimal<RightT> || is_floating_point<RightT> || is_big_int_v<RightT> || is_signed_v<RightT>;
 
                 /// Disallow int vs UUID conversion (but support int vs UInt128 conversion)
                 if constexpr ((bad_left && std::is_same_v<RightDataType, DataTypeUUID>) ||
@@ -2662,7 +2662,7 @@ struct ToNumberMonotonicity
         /// Float cases.
 
         /// When converting to Float, the conversion is always monotonic.
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
             return { .is_monotonic = true, .is_always_monotonic = true };
 
         const auto * low_cardinality = typeid_cast<const DataTypeLowCardinality *>(&type);
@@ -2875,6 +2875,7 @@ struct NameToInt32 { static constexpr auto name = "toInt32"; };
 struct NameToInt64 { static constexpr auto name = "toInt64"; };
 struct NameToInt128 { static constexpr auto name = "toInt128"; };
 struct NameToInt256 { static constexpr auto name = "toInt256"; };
+struct NameToBFloat16 { static constexpr auto name = "toBFloat16"; };
 struct NameToFloat32 { static constexpr auto name = "toFloat32"; };
 struct NameToFloat64 { static constexpr auto name = "toFloat64"; };
 struct NameToUUID { static constexpr auto name = "toUUID"; };
@@ -2893,6 +2894,7 @@ using FunctionToInt32 = FunctionConvert<DataTypeInt32, NameToInt32, ToNumberMono
 using FunctionToInt64 = FunctionConvert<DataTypeInt64, NameToInt64, ToNumberMonotonicity<Int64>>;
 using FunctionToInt128 = FunctionConvert<DataTypeInt128, NameToInt128, ToNumberMonotonicity<Int128>>;
 using FunctionToInt256 = FunctionConvert<DataTypeInt256, NameToInt256, ToNumberMonotonicity<Int256>>;
+using FunctionToBFloat16 = FunctionConvert<DataTypeBFloat16, NameToBFloat16, ToNumberMonotonicity<BFloat16>>;
 using FunctionToFloat32 = FunctionConvert<DataTypeFloat32, NameToFloat32, ToNumberMonotonicity<Float32>>;
 using FunctionToFloat64 = FunctionConvert<DataTypeFloat64, NameToFloat64, ToNumberMonotonicity<Float64>>;
 
@@ -2930,6 +2932,7 @@ template <> struct FunctionTo<DataTypeInt32> { using Type = FunctionToInt32; };
 template <> struct FunctionTo<DataTypeInt64> { using Type = FunctionToInt64; };
 template <> struct FunctionTo<DataTypeInt128> { using Type = FunctionToInt128; };
 template <> struct FunctionTo<DataTypeInt256> { using Type = FunctionToInt256; };
+template <> struct FunctionTo<DataTypeBFloat16> { using Type = FunctionToBFloat16; };
 template <> struct FunctionTo<DataTypeFloat32> { using Type = FunctionToFloat32; };
 template <> struct FunctionTo<DataTypeFloat64> { using Type = FunctionToFloat64; };
 
@@ -2972,6 +2975,7 @@ struct NameToInt32OrZero { static constexpr auto name = "toInt32OrZero"; };
 struct NameToInt64OrZero { static constexpr auto name = "toInt64OrZero"; };
 struct NameToInt128OrZero { static constexpr auto name = "toInt128OrZero"; };
 struct NameToInt256OrZero { static constexpr auto name = "toInt256OrZero"; };
+struct NameToBFloat16OrZero { static constexpr auto name = "toBFloat16OrZero"; };
 struct NameToFloat32OrZero { static constexpr auto name = "toFloat32OrZero"; };
 struct NameToFloat64OrZero { static constexpr auto name = "toFloat64OrZero"; };
 struct NameToDateOrZero { static constexpr auto name = "toDateOrZero"; };
@@ -2998,6 +3002,7 @@ using FunctionToInt32OrZero = FunctionConvertFromString<DataTypeInt32, NameToInt
 using FunctionToInt64OrZero = FunctionConvertFromString<DataTypeInt64, NameToInt64OrZero, ConvertFromStringExceptionMode::Zero>;
 using FunctionToInt128OrZero = FunctionConvertFromString<DataTypeInt128, NameToInt128OrZero, ConvertFromStringExceptionMode::Zero>;
 using FunctionToInt256OrZero = FunctionConvertFromString<DataTypeInt256, NameToInt256OrZero, ConvertFromStringExceptionMode::Zero>;
+using FunctionToBFloat16OrZero = FunctionConvertFromString<DataTypeBFloat16, NameToBFloat16OrZero, ConvertFromStringExceptionMode::Zero>;
 using FunctionToFloat32OrZero = FunctionConvertFromString<DataTypeFloat32, NameToFloat32OrZero, ConvertFromStringExceptionMode::Zero>;
 using FunctionToFloat64OrZero = FunctionConvertFromString<DataTypeFloat64, NameToFloat64OrZero, ConvertFromStringExceptionMode::Zero>;
 using FunctionToDateOrZero = FunctionConvertFromString<DataTypeDate, NameToDateOrZero, ConvertFromStringExceptionMode::Zero>;
@@ -3024,6 +3029,7 @@ struct NameToInt32OrNull { static constexpr auto name = "toInt32OrNull"; };
 struct NameToInt64OrNull { static constexpr auto name = "toInt64OrNull"; };
 struct NameToInt128OrNull { static constexpr auto name = "toInt128OrNull"; };
 struct NameToInt256OrNull { static constexpr auto name = "toInt256OrNull"; };
+struct NameToBFloat16OrNull { static constexpr auto name = "toBFloat16OrNull"; };
 struct NameToFloat32OrNull { static constexpr auto name = "toFloat32OrNull"; };
 struct NameToFloat64OrNull { static constexpr auto name = "toFloat64OrNull"; };
 struct NameToDateOrNull { static constexpr auto name = "toDateOrNull"; };
@@ -3050,6 +3056,7 @@ using FunctionToInt32OrNull = FunctionConvertFromString<DataTypeInt32, NameToInt
 using FunctionToInt64OrNull = FunctionConvertFromString<DataTypeInt64, NameToInt64OrNull, ConvertFromStringExceptionMode::Null>;
 using FunctionToInt128OrNull = FunctionConvertFromString<DataTypeInt128, NameToInt128OrNull, ConvertFromStringExceptionMode::Null>;
 using FunctionToInt256OrNull = FunctionConvertFromString<DataTypeInt256, NameToInt256OrNull, ConvertFromStringExceptionMode::Null>;
+using FunctionToBFloat16OrNull = FunctionConvertFromString<DataTypeBFloat16, NameToBFloat16OrNull, ConvertFromStringExceptionMode::Null>;
 using FunctionToFloat32OrNull = FunctionConvertFromString<DataTypeFloat32, NameToFloat32OrNull, ConvertFromStringExceptionMode::Null>;
 using FunctionToFloat64OrNull = FunctionConvertFromString<DataTypeFloat64, NameToFloat64OrNull, ConvertFromStringExceptionMode::Null>;
 using FunctionToDateOrNull = FunctionConvertFromString<DataTypeDate, NameToDateOrNull, ConvertFromStringExceptionMode::Null>;
@@ -5193,7 +5200,7 @@ private:
             if constexpr (is_any_of<ToDataType,
                 DataTypeUInt16, DataTypeUInt32, DataTypeUInt64, DataTypeUInt128, DataTypeUInt256,
                 DataTypeInt8, DataTypeInt16, DataTypeInt32, DataTypeInt64, DataTypeInt128, DataTypeInt256,
-                DataTypeFloat32, DataTypeFloat64,
+                DataTypeBFloat16, DataTypeFloat32, DataTypeFloat64,
                 DataTypeDate, DataTypeDate32, DataTypeDateTime,
                 DataTypeUUID, DataTypeIPv4, DataTypeIPv6>)
             {
@@ -5446,6 +5453,17 @@ REGISTER_FUNCTION(Conversion)
     factory.registerFunction<FunctionToInt64>();
     factory.registerFunction<FunctionToInt128>();
     factory.registerFunction<FunctionToInt256>();
+
+    factory.registerFunction<FunctionToBFloat16>(FunctionDocumentation{.description=R"(
+Converts Float32 to BFloat16 with losing the precision.
+
+Example:
+[example:typical]
+)",
+        .examples{
+            {"typical", "SELECT toBFloat16(12.3::Float32);", "12.3125"}},
+        .categories{"Conversion"}});
+
     factory.registerFunction<FunctionToFloat32>();
     factory.registerFunction<FunctionToFloat64>();
 
@@ -5484,6 +5502,31 @@ REGISTER_FUNCTION(Conversion)
     factory.registerFunction<FunctionToInt64OrZero>();
     factory.registerFunction<FunctionToInt128OrZero>();
     factory.registerFunction<FunctionToInt256OrZero>();
+
+    factory.registerFunction<FunctionToBFloat16OrZero>(FunctionDocumentation{.description=R"(
+Converts String to BFloat16.
+
+If the string does not represent a floating point value, the function returns zero.
+
+The function allows a silent loss of precision while converting from the string representation. In that case, it will return the truncated result.
+
+Example of successful conversion:
+[example:typical]
+
+Examples of not successful conversion:
+[example:invalid1]
+[example:invalid2]
+
+Example of a loss of precision:
+[example:precision]
+)",
+        .examples{
+            {"typical", "SELECT toBFloat16OrZero('12.3');", "12.3125"},
+            {"invalid1", "SELECT toBFloat16OrZero('abc');", "0"},
+            {"invalid2", "SELECT toBFloat16OrZero(' 1');", "0"},
+            {"precision", "SELECT toBFloat16OrZero('12.3456789');", "12.375"}},
+        .categories{"Conversion"}});
+
     factory.registerFunction<FunctionToFloat32OrZero>();
     factory.registerFunction<FunctionToFloat64OrZero>();
     factory.registerFunction<FunctionToDateOrZero>();
@@ -5512,6 +5555,31 @@ REGISTER_FUNCTION(Conversion)
     factory.registerFunction<FunctionToInt64OrNull>();
     factory.registerFunction<FunctionToInt128OrNull>();
     factory.registerFunction<FunctionToInt256OrNull>();
+
+    factory.registerFunction<FunctionToBFloat16OrNull>(FunctionDocumentation{.description=R"(
+Converts String to Nullable(BFloat16).
+
+If the string does not represent a floating point value, the function returns NULL.
+
+The function allows a silent loss of precision while converting from the string representation. In that case, it will return the truncated result.
+
+Example of successful conversion:
+[example:typical]
+
+Examples of not successful conversion:
+[example:invalid1]
+[example:invalid2]
+
+Example of a loss of precision:
+[example:precision]
+)",
+    .examples{
+        {"typical", "SELECT toBFloat16OrNull('12.3');", "12.3125"},
+        {"invalid1", "SELECT toBFloat16OrNull('abc');", "NULL"},
+        {"invalid2", "SELECT toBFloat16OrNull(' 1');", "NULL"},
+        {"precision", "SELECT toBFloat16OrNull('12.3456789');", "12.375"}},
+    .categories{"Conversion"}});
+
     factory.registerFunction<FunctionToFloat32OrNull>();
     factory.registerFunction<FunctionToFloat64OrNull>();
     factory.registerFunction<FunctionToDateOrNull>();
diff --git a/src/Functions/FunctionsRound.h b/src/Functions/FunctionsRound.h
index ed7fe1a5de1..6c9cc8a37b3 100644
--- a/src/Functions/FunctionsRound.h
+++ b/src/Functions/FunctionsRound.h
@@ -268,6 +268,19 @@ inline double roundWithMode(double x, RoundingMode mode)
     std::unreachable();
 }
 
+inline BFloat16 roundWithMode(BFloat16 x, RoundingMode mode)
+{
+    switch (mode)
+    {
+        case RoundingMode::Round: return BFloat16(nearbyintf(Float32(x)));
+        case RoundingMode::Floor: return BFloat16(floorf(Float32(x)));
+        case RoundingMode::Ceil: return BFloat16(ceilf(Float32(x)));
+        case RoundingMode::Trunc: return BFloat16(truncf(Float32(x)));
+    }
+
+    std::unreachable();
+}
+
 template <typename T>
 class FloatRoundingComputationBase<T, Vectorize::No>
 {
@@ -285,10 +298,15 @@ public:
 
     static VectorType prepare(size_t scale)
     {
-        return load1(scale);
+        return load1(ScalarType(scale));
     }
 };
 
+template <>
+class FloatRoundingComputationBase<BFloat16, Vectorize::Yes> : public FloatRoundingComputationBase<BFloat16, Vectorize::No>
+{
+};
+
 
 /** Implementation of low-level round-off functions for floating-point values.
   */
@@ -511,7 +529,7 @@ template <typename T, RoundingMode rounding_mode, TieBreakingMode tie_breaking_m
 struct Dispatcher
 {
     template <ScaleMode scale_mode>
-    using FunctionRoundingImpl = std::conditional_t<std::is_floating_point_v<T>,
+    using FunctionRoundingImpl = std::conditional_t<is_floating_point<T>,
         FloatRoundingImpl<T, rounding_mode, scale_mode>,
         IntegerRoundingImpl<T, rounding_mode, scale_mode, tie_breaking_mode>>;
 
diff --git a/src/Functions/FunctionsVisitParam.h b/src/Functions/FunctionsVisitParam.h
index a77fa740f9c..76d4a13687c 100644
--- a/src/Functions/FunctionsVisitParam.h
+++ b/src/Functions/FunctionsVisitParam.h
@@ -57,7 +57,7 @@ struct ExtractNumericType
         ResultType x = 0;
         if (!in.eof())
         {
-            if constexpr (std::is_floating_point_v<NumericType>)
+            if constexpr (is_floating_point<NumericType>)
                 tryReadFloatText(x, in);
             else
                 tryReadIntText(x, in);
diff --git a/src/Functions/PolygonUtils.h b/src/Functions/PolygonUtils.h
index bf8241774a6..601ffcb00b4 100644
--- a/src/Functions/PolygonUtils.h
+++ b/src/Functions/PolygonUtils.h
@@ -583,7 +583,7 @@ struct CallPointInPolygon<Type, Types ...>
     template <typename PointInPolygonImpl>
     static ColumnPtr call(const IColumn & x, const IColumn & y, PointInPolygonImpl && impl)
     {
-        using Impl = TypeListChangeRoot<CallPointInPolygon, TypeListIntAndFloat>;
+        using Impl = TypeListChangeRoot<CallPointInPolygon, TypeListNativeNumber>;
         if (auto column = typeid_cast<const ColumnVector<Type> *>(&x))
             return Impl::template call<Type>(*column, y, impl);
         return CallPointInPolygon<Types ...>::call(x, y, impl);
@@ -609,7 +609,7 @@ struct CallPointInPolygon<>
 template <typename PointInPolygonImpl>
 NO_INLINE ColumnPtr pointInPolygon(const IColumn & x, const IColumn & y, PointInPolygonImpl && impl)
 {
-    using Impl = TypeListChangeRoot<CallPointInPolygon, TypeListIntAndFloat>;
+    using Impl = TypeListChangeRoot<CallPointInPolygon, TypeListNativeNumber>;
     return Impl::call(x, y, impl);
 }
 
diff --git a/src/Functions/abs.cpp b/src/Functions/abs.cpp
index 742d3b85619..0b701bb7108 100644
--- a/src/Functions/abs.cpp
+++ b/src/Functions/abs.cpp
@@ -22,7 +22,7 @@ struct AbsImpl
             return a < 0 ? static_cast<ResultType>(~a) + 1 : static_cast<ResultType>(a);
         else if constexpr (is_integer<A> && is_unsigned_v<A>)
             return static_cast<ResultType>(a);
-        else if constexpr (std::is_floating_point_v<A>)
+        else if constexpr (is_floating_point<A>)
             return static_cast<ResultType>(std::abs(a));
     }
 
diff --git a/src/Functions/array/arrayAggregation.cpp b/src/Functions/array/arrayAggregation.cpp
index bb2503886f1..d073c79b832 100644
--- a/src/Functions/array/arrayAggregation.cpp
+++ b/src/Functions/array/arrayAggregation.cpp
@@ -87,7 +87,7 @@ struct ArrayAggregateResultImpl<ArrayElement, AggregateOperation::sum>
         std::conditional_t<std::is_same_v<ArrayElement, Decimal128>, Decimal128,
         std::conditional_t<std::is_same_v<ArrayElement, Decimal256>, Decimal256,
         std::conditional_t<std::is_same_v<ArrayElement, DateTime64>, Decimal128,
-        std::conditional_t<std::is_floating_point_v<ArrayElement>, Float64,
+        std::conditional_t<is_floating_point<ArrayElement>, Float64,
         std::conditional_t<std::is_signed_v<ArrayElement>, Int64,
             UInt64>>>>>>>>>>>;
 };
diff --git a/src/Functions/array/arrayDistance.cpp b/src/Functions/array/arrayDistance.cpp
index a1f48747eb6..8a6c2649898 100644
--- a/src/Functions/array/arrayDistance.cpp
+++ b/src/Functions/array/arrayDistance.cpp
@@ -14,6 +14,7 @@
 #include <immintrin.h>
 #endif
 
+
 namespace DB
 {
 namespace ErrorCodes
@@ -34,7 +35,7 @@ struct L1Distance
     template <typename FloatType>
     struct State
     {
-        FloatType sum = 0;
+        FloatType sum{};
     };
 
     template <typename ResultType>
@@ -65,7 +66,7 @@ struct L2Distance
     template <typename FloatType>
     struct State
     {
-        FloatType sum = 0;
+        FloatType sum{};
     };
 
     template <typename ResultType>
@@ -90,19 +91,17 @@ struct L2Distance
         size_t & i_y,
         State<ResultType> & state)
     {
-        static constexpr bool is_float32 = std::is_same_v<ResultType, Float32>;
-
         __m512 sums;
-        if constexpr (is_float32)
+        if constexpr (sizeof(ResultType) <= 4)
             sums = _mm512_setzero_ps();
         else
             sums = _mm512_setzero_pd();
 
-        constexpr size_t n = is_float32 ? 16 : 8;
+        constexpr size_t n = sizeof(__m512) / sizeof(ResultType);
 
         for (; i_x + n < i_max; i_x += n, i_y += n)
         {
-            if constexpr (is_float32)
+            if constexpr (sizeof(ResultType) == 4)
             {
                 __m512 x = _mm512_loadu_ps(data_x + i_x);
                 __m512 y = _mm512_loadu_ps(data_y + i_y);
@@ -118,11 +117,38 @@ struct L2Distance
             }
         }
 
-        if constexpr (is_float32)
+        if constexpr (sizeof(ResultType) <= 4)
             state.sum = _mm512_reduce_add_ps(sums);
         else
             state.sum = _mm512_reduce_add_pd(sums);
     }
+
+    AVX512BF16_FUNCTION_SPECIFIC_ATTRIBUTE static void accumulateCombineBF16(
+        const BFloat16 * __restrict data_x,
+        const BFloat16 * __restrict data_y,
+        size_t i_max,
+        size_t & i_x,
+        size_t & i_y,
+        State<Float32> & state)
+    {
+        __m512 sums = _mm512_setzero_ps();
+        constexpr size_t n = sizeof(__m512) / sizeof(BFloat16);
+
+        for (; i_x + n < i_max; i_x += n, i_y += n)
+        {
+            __m512 x_1 = _mm512_cvtpbh_ps(_mm256_loadu_ps(reinterpret_cast<const Float32 *>(data_x + i_x)));
+            __m512 x_2 = _mm512_cvtpbh_ps(_mm256_loadu_ps(reinterpret_cast<const Float32 *>(data_x + i_x + n / 2)));
+            __m512 y_1 = _mm512_cvtpbh_ps(_mm256_loadu_ps(reinterpret_cast<const Float32 *>(data_y + i_y)));
+            __m512 y_2 = _mm512_cvtpbh_ps(_mm256_loadu_ps(reinterpret_cast<const Float32 *>(data_y + i_y + n / 2)));
+
+            __m512 differences_1 = _mm512_sub_ps(x_1, y_1);
+            __m512 differences_2 = _mm512_sub_ps(x_2, y_2);
+            sums = _mm512_fmadd_ps(differences_1, differences_1, sums);
+            sums = _mm512_fmadd_ps(differences_2, differences_2, sums);
+        }
+
+        state.sum = _mm512_reduce_add_ps(sums);
+    }
 #endif
 
     template <typename ResultType>
@@ -156,13 +182,13 @@ struct LpDistance
     template <typename FloatType>
     struct State
     {
-        FloatType sum = 0;
+        FloatType sum{};
     };
 
     template <typename ResultType>
     static void accumulate(State<ResultType> & state, ResultType x, ResultType y, const ConstParams & params)
     {
-        state.sum += static_cast<ResultType>(std::pow(fabs(x - y), params.power));
+        state.sum += static_cast<ResultType>(pow(fabs(x - y), params.power));
     }
 
     template <typename ResultType>
@@ -174,7 +200,7 @@ struct LpDistance
     template <typename ResultType>
     static ResultType finalize(const State<ResultType> & state, const ConstParams & params)
     {
-        return static_cast<ResultType>(std::pow(state.sum, params.inverted_power));
+        return static_cast<ResultType>(pow(state.sum, params.inverted_power));
     }
 };
 
@@ -187,7 +213,7 @@ struct LinfDistance
     template <typename FloatType>
     struct State
     {
-        FloatType dist = 0;
+        FloatType dist{};
     };
 
     template <typename ResultType>
@@ -218,9 +244,9 @@ struct CosineDistance
     template <typename FloatType>
     struct State
     {
-        FloatType dot_prod = 0;
-        FloatType x_squared = 0;
-        FloatType y_squared = 0;
+        FloatType dot_prod{};
+        FloatType x_squared{};
+        FloatType y_squared{};
     };
 
     template <typename ResultType>
@@ -249,13 +275,11 @@ struct CosineDistance
         size_t & i_y,
         State<ResultType> & state)
     {
-        static constexpr bool is_float32 = std::is_same_v<ResultType, Float32>;
-
         __m512 dot_products;
         __m512 x_squareds;
         __m512 y_squareds;
 
-        if constexpr (is_float32)
+        if constexpr (sizeof(ResultType) <= 4)
         {
             dot_products = _mm512_setzero_ps();
             x_squareds = _mm512_setzero_ps();
@@ -268,11 +292,11 @@ struct CosineDistance
             y_squareds = _mm512_setzero_pd();
         }
 
-        constexpr size_t n = is_float32 ? 16 : 8;
+        constexpr size_t n = sizeof(__m512) / sizeof(ResultType);
 
         for (; i_x + n < i_max; i_x += n, i_y += n)
         {
-            if constexpr (is_float32)
+            if constexpr (sizeof(ResultType) == 4)
             {
                 __m512 x = _mm512_loadu_ps(data_x + i_x);
                 __m512 y = _mm512_loadu_ps(data_y + i_y);
@@ -290,7 +314,7 @@ struct CosineDistance
             }
         }
 
-        if constexpr (is_float32)
+        if constexpr (sizeof(ResultType) == 4)
         {
             state.dot_prod = _mm512_reduce_add_ps(dot_products);
             state.x_squared = _mm512_reduce_add_ps(x_squareds);
@@ -303,16 +327,48 @@ struct CosineDistance
             state.y_squared = _mm512_reduce_add_pd(y_squareds);
         }
     }
+
+    AVX512BF16_FUNCTION_SPECIFIC_ATTRIBUTE static void accumulateCombineBF16(
+        const BFloat16 * __restrict data_x,
+        const BFloat16 * __restrict data_y,
+        size_t i_max,
+        size_t & i_x,
+        size_t & i_y,
+        State<Float32> & state)
+    {
+        __m512 dot_products;
+        __m512 x_squareds;
+        __m512 y_squareds;
+
+        dot_products = _mm512_setzero_ps();
+        x_squareds = _mm512_setzero_ps();
+        y_squareds = _mm512_setzero_ps();
+
+        constexpr size_t n = sizeof(__m512) / sizeof(BFloat16);
+
+        for (; i_x + n < i_max; i_x += n, i_y += n)
+        {
+            __m512 x = _mm512_loadu_ps(data_x + i_x);
+            __m512 y = _mm512_loadu_ps(data_y + i_y);
+            dot_products = _mm512_dpbf16_ps(dot_products, x, y);
+            x_squareds = _mm512_dpbf16_ps(x_squareds, x, x);
+            y_squareds = _mm512_dpbf16_ps(y_squareds, y, y);
+        }
+
+        state.dot_prod = _mm512_reduce_add_ps(dot_products);
+        state.x_squared = _mm512_reduce_add_ps(x_squareds);
+        state.y_squared = _mm512_reduce_add_ps(y_squareds);
+    }
 #endif
 
     template <typename ResultType>
     static ResultType finalize(const State<ResultType> & state, const ConstParams &)
     {
-        return 1 - state.dot_prod / sqrt(state.x_squared * state.y_squared);
+        return 1.0f - state.dot_prod / sqrt(state.x_squared * state.y_squared);
     }
 };
 
-template <class Kernel>
+template <typename Kernel>
 class FunctionArrayDistance : public IFunction
 {
 public:
@@ -352,12 +408,13 @@ public:
             case TypeIndex::Float64:
                 return std::make_shared<DataTypeFloat64>();
             case TypeIndex::Float32:
+            case TypeIndex::BFloat16:
                 return std::make_shared<DataTypeFloat32>();
             default:
                 throw Exception(
                     ErrorCodes::ILLEGAL_TYPE_OF_ARGUMENT,
                     "Arguments of function {} has nested type {}. "
-                    "Supported types: UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64, Float32, Float64.",
+                    "Supported types: UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64, BFloat16, Float32, Float64.",
                     getName(),
                     common_type->getName());
         }
@@ -369,10 +426,8 @@ public:
         {
             case TypeIndex::Float32:
                 return executeWithResultType<Float32>(arguments, input_rows_count);
-                break;
             case TypeIndex::Float64:
                 return executeWithResultType<Float64>(arguments, input_rows_count);
-                break;
             default:
                 throw Exception(ErrorCodes::LOGICAL_ERROR, "Unexpected result type {}", result_type->getName());
         }
@@ -388,6 +443,7 @@ public:
     ACTION(Int16)   \
     ACTION(Int32)   \
     ACTION(Int64)   \
+    ACTION(BFloat16) \
     ACTION(Float32) \
     ACTION(Float64)
 
@@ -412,7 +468,7 @@ private:
                 throw Exception(
                     ErrorCodes::ILLEGAL_TYPE_OF_ARGUMENT,
                     "Arguments of function {} has nested type {}. "
-                    "Supported types: UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64, Float32, Float64.",
+                    "Supported types: UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64, BFloat16, Float32, Float64.",
                     getName(),
                     type_x->getName());
         }
@@ -437,7 +493,7 @@ private:
                 throw Exception(
                     ErrorCodes::ILLEGAL_TYPE_OF_ARGUMENT,
                     "Arguments of function {} has nested type {}. "
-                    "Supported types: UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64, Float32, Float64.",
+                    "Supported types: UInt8, UInt16, UInt32, UInt64, Int8, Int16, Int32, Int64, BFloat16, Float32, Float64.",
                     getName(),
                     type_y->getName());
         }
@@ -446,14 +502,10 @@ private:
     template <typename ResultType, typename LeftType, typename RightType>
     ColumnPtr executeWithResultTypeAndLeftTypeAndRightType(ColumnPtr col_x, ColumnPtr col_y, size_t input_rows_count, const ColumnsWithTypeAndName & arguments) const
     {
-        if (typeid_cast<const ColumnConst *>(col_x.get()))
-        {
+        if (col_x->isConst())
             return executeWithLeftArgConst<ResultType, LeftType, RightType>(col_x, col_y, input_rows_count, arguments);
-        }
-        if (typeid_cast<const ColumnConst *>(col_y.get()))
-        {
+        if (col_y->isConst())
             return executeWithLeftArgConst<ResultType, RightType, LeftType>(col_y, col_x, input_rows_count, arguments);
-        }
 
         const auto & array_x = *assert_cast<const ColumnArray *>(col_x.get());
         const auto & array_y = *assert_cast<const ColumnArray *>(col_y.get());
@@ -497,7 +549,7 @@ private:
                     state, static_cast<ResultType>(data_x[prev]), static_cast<ResultType>(data_y[prev]), kernel_params);
             }
             result_data[row] = Kernel::finalize(state, kernel_params);
-            row++;
+            ++row;
         }
         return col_res;
     }
@@ -548,36 +600,52 @@ private:
 
             /// SIMD optimization: process multiple elements in both input arrays at once.
             /// To avoid combinatorial explosion of SIMD kernels, focus on
-            /// - the two most common input/output types (Float32 x Float32) --> Float32 and (Float64 x Float64) --> Float64 instead of 10 x
-            ///   10 input types x 2 output types,
+            /// - the three most common input/output types (BFloat16 x BFloat16) --> Float32,
+            ///   (Float32 x Float32) --> Float32 and (Float64 x Float64) --> Float64
+            ///   instead of 10 x 10 input types x 2 output types,
             /// - const/non-const inputs instead of non-const/non-const inputs
             /// - the two most common metrics L2 and cosine distance,
             /// - the most powerful SIMD instruction set (AVX-512F).
+            bool processed = false;
 #if USE_MULTITARGET_CODE
-            if constexpr (std::is_same_v<ResultType, LeftType> && std::is_same_v<ResultType, RightType>) /// ResultType is Float32 or Float64
+            /// ResultType is Float32 or Float64
+            if constexpr (std::is_same_v<Kernel, L2Distance> || std::is_same_v<Kernel, CosineDistance>)
             {
-                if constexpr (std::is_same_v<Kernel, L2Distance>
-                           || std::is_same_v<Kernel, CosineDistance>)
+                if constexpr (std::is_same_v<ResultType, LeftType> && std::is_same_v<ResultType, RightType>)
                 {
                     if (isArchSupported(TargetArch::AVX512F))
+                    {
                         Kernel::template accumulateCombine<ResultType>(data_x.data(), data_y.data(), i + offsets_x[0], i, prev, state);
+                        processed = true;
+                    }
+                }
+                else if constexpr (std::is_same_v<Float32, ResultType> && std::is_same_v<BFloat16, LeftType> && std::is_same_v<BFloat16, RightType>)
+                {
+                    if (isArchSupported(TargetArch::AVX512BF16))
+                    {
+                        Kernel::accumulateCombineBF16(data_x.data(), data_y.data(), i + offsets_x[0], i, prev, state);
+                        processed = true;
+                    }
                 }
             }
-#else
-            /// Process chunks in vectorized manner
-            static constexpr size_t VEC_SIZE = 4;
-            typename Kernel::template State<ResultType> states[VEC_SIZE];
-            for (; prev + VEC_SIZE < off; i += VEC_SIZE, prev += VEC_SIZE)
+#endif
+            if (!processed)
             {
-                for (size_t s = 0; s < VEC_SIZE; ++s)
-                    Kernel::template accumulate<ResultType>(
-                        states[s], static_cast<ResultType>(data_x[i + s]), static_cast<ResultType>(data_y[prev + s]), kernel_params);
+                /// Process chunks in a vectorized manner.
+                static constexpr size_t VEC_SIZE = 32;
+                typename Kernel::template State<ResultType> states[VEC_SIZE];
+                for (; prev + VEC_SIZE < off; i += VEC_SIZE, prev += VEC_SIZE)
+                {
+                    for (size_t s = 0; s < VEC_SIZE; ++s)
+                        Kernel::template accumulate<ResultType>(
+                            states[s], static_cast<ResultType>(data_x[i + s]), static_cast<ResultType>(data_y[prev + s]), kernel_params);
+                }
+
+                for (const auto & other_state : states)
+                    Kernel::template combine<ResultType>(state, other_state, kernel_params);
             }
 
-            for (const auto & other_state : states)
-                Kernel::template combine<ResultType>(state, other_state, kernel_params);
-#endif
-            /// Process the tail
+            /// Process the tail.
             for (; prev < off; ++i, ++prev)
             {
                 Kernel::template accumulate<ResultType>(
@@ -638,4 +706,5 @@ FunctionPtr createFunctionArrayL2SquaredDistance(ContextPtr context_) { return F
 FunctionPtr createFunctionArrayLpDistance(ContextPtr context_) { return FunctionArrayDistance<LpDistance>::create(context_); }
 FunctionPtr createFunctionArrayLinfDistance(ContextPtr context_) { return FunctionArrayDistance<LinfDistance>::create(context_); }
 FunctionPtr createFunctionArrayCosineDistance(ContextPtr context_) { return FunctionArrayDistance<CosineDistance>::create(context_); }
+
 }
diff --git a/src/Functions/array/mapPopulateSeries.cpp b/src/Functions/array/mapPopulateSeries.cpp
index 4083678b209..d2ca48ff186 100644
--- a/src/Functions/array/mapPopulateSeries.cpp
+++ b/src/Functions/array/mapPopulateSeries.cpp
@@ -452,23 +452,29 @@ private:
             using ValueType = typename Types::RightType;
 
             static constexpr bool key_and_value_are_numbers = IsDataTypeNumber<KeyType> && IsDataTypeNumber<ValueType>;
-            static constexpr bool key_is_float = std::is_same_v<KeyType, DataTypeFloat32> || std::is_same_v<KeyType, DataTypeFloat64>;
 
-            if constexpr (key_and_value_are_numbers && !key_is_float)
+            if constexpr (key_and_value_are_numbers)
             {
-                using KeyFieldType = typename KeyType::FieldType;
-                using ValueFieldType = typename ValueType::FieldType;
+                if constexpr (is_floating_point<typename KeyType::FieldType>)
+                {
+                    return false;
+                }
+                else
+                {
+                    using KeyFieldType = typename KeyType::FieldType;
+                    using ValueFieldType = typename ValueType::FieldType;
 
-                executeImplTyped<KeyFieldType, ValueFieldType>(
-                    input.key_column,
-                    input.value_column,
-                    input.offsets_column,
-                    input.max_key_column,
-                    std::move(result_columns.result_key_column),
-                    std::move(result_columns.result_value_column),
-                    std::move(result_columns.result_offset_column));
+                    executeImplTyped<KeyFieldType, ValueFieldType>(
+                        input.key_column,
+                        input.value_column,
+                        input.offsets_column,
+                        input.max_key_column,
+                        std::move(result_columns.result_key_column),
+                        std::move(result_columns.result_value_column),
+                        std::move(result_columns.result_offset_column));
 
-                return true;
+                    return true;
+                }
             }
 
             return false;
diff --git a/src/Functions/divide.cpp b/src/Functions/divide.cpp
index 7c67245c382..3947ba2d142 100644
--- a/src/Functions/divide.cpp
+++ b/src/Functions/divide.cpp
@@ -18,7 +18,7 @@ struct DivideFloatingImpl
     template <typename Result = ResultType>
     static NO_SANITIZE_UNDEFINED Result apply(A a [[maybe_unused]], B b [[maybe_unused]])
     {
-        return static_cast<Result>(a) / b;
+        return static_cast<Result>(a) / static_cast<Result>(b);
     }
 
 #if USE_EMBEDDED_COMPILER
diff --git a/src/Functions/exp.cpp b/src/Functions/exp.cpp
index e67cbd6d819..24f1d313831 100644
--- a/src/Functions/exp.cpp
+++ b/src/Functions/exp.cpp
@@ -3,6 +3,12 @@
 
 namespace DB
 {
+
+namespace ErrorCodes
+{
+extern const int NOT_IMPLEMENTED;
+}
+
 namespace
 {
 
@@ -21,7 +27,14 @@ namespace
         template <typename T>
         static void execute(const T * src, size_t size, T * dst)
         {
-            NFastOps::Exp<true>(src, size, dst);
+            if constexpr (std::is_same_v<T, BFloat16>)
+            {
+                throw Exception(ErrorCodes::NOT_IMPLEMENTED, "Function `{}` is not implemented for BFloat16", name);
+            }
+            else
+            {
+                NFastOps::Exp<true>(src, size, dst);
+            }
         }
     };
 }
diff --git a/src/Functions/factorial.cpp b/src/Functions/factorial.cpp
index 9b319caad63..83c184f50cc 100644
--- a/src/Functions/factorial.cpp
+++ b/src/Functions/factorial.cpp
@@ -21,7 +21,7 @@ struct FactorialImpl
 
     static NO_SANITIZE_UNDEFINED ResultType apply(A a)
     {
-        if constexpr (std::is_floating_point_v<A> || is_over_big_int<A>)
+        if constexpr (is_floating_point<A> || is_over_big_int<A>)
             throw Exception(
                 ErrorCodes::ILLEGAL_TYPE_OF_ARGUMENT,
                 "Illegal type of argument of function factorial, should not be floating point or big int");
diff --git a/src/Functions/log.cpp b/src/Functions/log.cpp
index 8bebdb8d7bd..49fc509634b 100644
--- a/src/Functions/log.cpp
+++ b/src/Functions/log.cpp
@@ -4,6 +4,11 @@
 namespace DB
 {
 
+namespace ErrorCodes
+{
+extern const int NOT_IMPLEMENTED;
+}
+
 namespace
 {
 
@@ -20,7 +25,14 @@ struct LogName { static constexpr auto name = "log"; };
         template <typename T>
         static void execute(const T * src, size_t size, T * dst)
         {
-            NFastOps::Log<true>(src, size, dst);
+            if constexpr (std::is_same_v<T, BFloat16>)
+            {
+                throw Exception(ErrorCodes::NOT_IMPLEMENTED, "Function `{}` is not implemented for BFloat16", name);
+            }
+            else
+            {
+                NFastOps::Log<true>(src, size, dst);
+            }
         }
     };
 
diff --git a/src/Functions/minus.cpp b/src/Functions/minus.cpp
index f3b9b8a7bcb..a372e8d5d78 100644
--- a/src/Functions/minus.cpp
+++ b/src/Functions/minus.cpp
@@ -17,8 +17,8 @@ struct MinusImpl
     {
         if constexpr (is_big_int_v<A> || is_big_int_v<B>)
         {
-            using CastA = std::conditional_t<std::is_floating_point_v<B>, B, A>;
-            using CastB = std::conditional_t<std::is_floating_point_v<A>, A, B>;
+            using CastA = std::conditional_t<is_floating_point<B>, B, A>;
+            using CastB = std::conditional_t<is_floating_point<A>, A, B>;
 
             return static_cast<Result>(static_cast<CastA>(a)) - static_cast<Result>(static_cast<CastB>(b));
         }
diff --git a/src/Functions/moduloOrZero.cpp b/src/Functions/moduloOrZero.cpp
index cd7873b3b9e..5a4d1539345 100644
--- a/src/Functions/moduloOrZero.cpp
+++ b/src/Functions/moduloOrZero.cpp
@@ -17,7 +17,7 @@ struct ModuloOrZeroImpl
     template <typename Result = ResultType>
     static Result apply(A a, B b)
     {
-        if constexpr (std::is_floating_point_v<ResultType>)
+        if constexpr (is_floating_point<ResultType>)
         {
             /// This computation is similar to `fmod` but the latter is not inlined and has 40 times worse performance.
             return ResultType(a) - trunc(ResultType(a) / ResultType(b)) * ResultType(b);
diff --git a/src/Functions/multiply.cpp b/src/Functions/multiply.cpp
index 67b6fff6b58..740ab81d0d9 100644
--- a/src/Functions/multiply.cpp
+++ b/src/Functions/multiply.cpp
@@ -18,8 +18,8 @@ struct MultiplyImpl
     {
         if constexpr (is_big_int_v<A> || is_big_int_v<B>)
         {
-            using CastA = std::conditional_t<std::is_floating_point_v<B>, B, A>;
-            using CastB = std::conditional_t<std::is_floating_point_v<A>, A, B>;
+            using CastA = std::conditional_t<is_floating_point<B>, B, A>;
+            using CastB = std::conditional_t<is_floating_point<A>, A, B>;
 
             return static_cast<Result>(static_cast<CastA>(a)) * static_cast<Result>(static_cast<CastB>(b));
         }
diff --git a/src/Functions/plus.cpp b/src/Functions/plus.cpp
index ffb0fe2ade7..26921713f78 100644
--- a/src/Functions/plus.cpp
+++ b/src/Functions/plus.cpp
@@ -19,8 +19,8 @@ struct PlusImpl
         /// Next everywhere, static_cast - so that there is no wrong result in expressions of the form Int64 c = UInt32(a) * Int32(-1).
         if constexpr (is_big_int_v<A> || is_big_int_v<B>)
         {
-            using CastA = std::conditional_t<std::is_floating_point_v<B>, B, A>;
-            using CastB = std::conditional_t<std::is_floating_point_v<A>, A, B>;
+            using CastA = std::conditional_t<is_floating_point<B>, B, A>;
+            using CastB = std::conditional_t<is_floating_point<A>, A, B>;
 
             return static_cast<Result>(static_cast<CastA>(a)) + static_cast<Result>(static_cast<CastB>(b));
         }
diff --git a/src/Functions/sigmoid.cpp b/src/Functions/sigmoid.cpp
index d121bdc7389..bb9710a15fe 100644
--- a/src/Functions/sigmoid.cpp
+++ b/src/Functions/sigmoid.cpp
@@ -3,6 +3,12 @@
 
 namespace DB
 {
+
+namespace ErrorCodes
+{
+extern const int NOT_IMPLEMENTED;
+}
+
 namespace
 {
 
@@ -21,7 +27,14 @@ namespace
         template <typename T>
         static void execute(const T * src, size_t size, T * dst)
         {
-            NFastOps::Sigmoid<>(src, size, dst);
+            if constexpr (std::is_same_v<T, BFloat16>)
+            {
+                throw Exception(ErrorCodes::NOT_IMPLEMENTED, "Function `{}` is not implemented for BFloat16", name);
+            }
+            else
+            {
+                NFastOps::Sigmoid<>(src, size, dst);
+            }
         }
     };
 }
@@ -47,4 +60,3 @@ REGISTER_FUNCTION(Sigmoid)
 }
 
 }
-
diff --git a/src/Functions/sign.cpp b/src/Functions/sign.cpp
index 914e1ad9e1f..63d59e3e9da 100644
--- a/src/Functions/sign.cpp
+++ b/src/Functions/sign.cpp
@@ -13,7 +13,7 @@ struct SignImpl
 
     static NO_SANITIZE_UNDEFINED ResultType apply(A a)
     {
-        if constexpr (is_decimal<A> || std::is_floating_point_v<A>)
+        if constexpr (is_decimal<A> || is_floating_point<A>)
             return a < A(0) ? -1 : a == A(0) ? 0 : 1;
         else if constexpr (is_signed_v<A>)
             return a < 0 ? -1 : a == 0 ? 0 : 1;
diff --git a/src/Functions/tanh.cpp b/src/Functions/tanh.cpp
index 62755737f70..d0e1440485b 100644
--- a/src/Functions/tanh.cpp
+++ b/src/Functions/tanh.cpp
@@ -3,6 +3,12 @@
 
 namespace DB
 {
+
+namespace ErrorCodes
+{
+extern const int NOT_IMPLEMENTED;
+}
+
 namespace
 {
 
@@ -19,7 +25,14 @@ struct TanhName { static constexpr auto name = "tanh"; };
         template <typename T>
         static void execute(const T * src, size_t size, T * dst)
         {
-            NFastOps::Tanh<>(src, size, dst);
+            if constexpr (std::is_same_v<T, BFloat16>)
+            {
+                throw Exception(ErrorCodes::NOT_IMPLEMENTED, "Function `{}` is not implemented for BFloat16", name);
+            }
+            else
+            {
+                NFastOps::Tanh<>(src, size, dst);
+            }
         }
     };
 
diff --git a/src/IO/ReadHelpers.h b/src/IO/ReadHelpers.h
index 50131781e77..65808f76751 100644
--- a/src/IO/ReadHelpers.h
+++ b/src/IO/ReadHelpers.h
@@ -1420,7 +1420,9 @@ inline bool tryReadText(UUID & x, ReadBuffer & buf) { return tryReadUUIDText(x,
 inline bool tryReadText(IPv4 & x, ReadBuffer & buf) { return tryReadIPv4Text(x, buf); }
 inline bool tryReadText(IPv6 & x, ReadBuffer & buf) { return tryReadIPv6Text(x, buf); }
 
-inline void readText(is_floating_point auto & x, ReadBuffer & buf) { readFloatText(x, buf); }
+template <typename T>
+requires is_floating_point<T>
+inline void readText(T & x, ReadBuffer & buf) { readFloatText(x, buf); }
 
 inline void readText(String & x, ReadBuffer & buf) { readEscapedString(x, buf); }
 
diff --git a/src/IO/WriteHelpers.h b/src/IO/WriteHelpers.h
index 27363093549..f7b2504f664 100644
--- a/src/IO/WriteHelpers.h
+++ b/src/IO/WriteHelpers.h
@@ -150,11 +150,12 @@ inline void writeBoolText(bool x, WriteBuffer & buf)
 
 
 template <typename T>
+requires is_floating_point<T>
 inline size_t writeFloatTextFastPath(T x, char * buffer)
 {
     Int64 result = 0;
 
-    if constexpr (std::is_same_v<T, double>)
+    if constexpr (std::is_same_v<T, Float64>)
     {
         /// The library Ryu has low performance on integers.
         /// This workaround improves performance 6..10 times.
@@ -164,13 +165,22 @@ inline size_t writeFloatTextFastPath(T x, char * buffer)
         else
             result = jkj::dragonbox::to_chars_n(x, buffer) - buffer;
     }
-    else
+    else if constexpr (std::is_same_v<T, Float32>)
     {
         if (DecomposedFloat32(x).isIntegerInRepresentableRange())
             result = itoa(Int32(x), buffer) - buffer;
         else
             result = jkj::dragonbox::to_chars_n(x, buffer) - buffer;
     }
+    else if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 f32 = Float32(x);
+
+        if (DecomposedFloat32(f32).isIntegerInRepresentableRange())
+            result = itoa(Int32(f32), buffer) - buffer;
+        else
+            result = jkj::dragonbox::to_chars_n(f32, buffer) - buffer;
+    }
 
     if (result <= 0)
         throw Exception(ErrorCodes::CANNOT_PRINT_FLOAT_OR_DOUBLE_NUMBER, "Cannot print floating point number");
@@ -178,10 +188,9 @@ inline size_t writeFloatTextFastPath(T x, char * buffer)
 }
 
 template <typename T>
+requires is_floating_point<T>
 inline void writeFloatText(T x, WriteBuffer & buf)
 {
-    static_assert(std::is_same_v<T, double> || std::is_same_v<T, float>, "Argument for writeFloatText must be float or double");
-
     using Converter = DoubleConverter<false>;
     if (likely(buf.available() >= Converter::MAX_REPRESENTATION_LENGTH))
     {
@@ -540,9 +549,9 @@ void writeJSONNumber(T x, WriteBuffer & ostr, const FormatSettings & settings)
         writeCString("null", ostr);
     else
     {
-        if constexpr (std::is_floating_point_v<T>)
+        if constexpr (is_floating_point<T>)
         {
-            if (std::signbit(x))
+            if (signBit(x))
             {
                 if (isNaN(x))
                     writeCString("-nan", ostr);
@@ -798,7 +807,6 @@ inline void writeXMLStringForTextElement(std::string_view s, WriteBuffer & buf)
 }
 
 /// @brief Serialize `uuid` into an array of characters in big-endian byte order.
-/// @param uuid UUID to serialize.
 /// @return Array of characters in big-endian byte order.
 std::array<char, 36> formatUUID(const UUID & uuid);
 
@@ -1099,7 +1107,9 @@ inline void writeText(is_integer auto x, WriteBuffer & buf)
         writeIntText(x, buf);
 }
 
-inline void writeText(is_floating_point auto x, WriteBuffer & buf) { writeFloatText(x, buf); }
+template <typename T>
+requires is_floating_point<T>
+inline void writeText(T x, WriteBuffer & buf) { writeFloatText(x, buf); }
 
 inline void writeText(is_enum auto x, WriteBuffer & buf) { writeText(magic_enum::enum_name(x), buf); }
 
diff --git a/src/IO/readFloatText.cpp b/src/IO/readFloatText.cpp
index 17ccc1b25b7..fb3c86fd7b6 100644
--- a/src/IO/readFloatText.cpp
+++ b/src/IO/readFloatText.cpp
@@ -47,26 +47,35 @@ void assertNaN(ReadBuffer & buf)
 }
 
 
+template void readFloatTextPrecise<BFloat16>(BFloat16 &, ReadBuffer &);
 template void readFloatTextPrecise<Float32>(Float32 &, ReadBuffer &);
 template void readFloatTextPrecise<Float64>(Float64 &, ReadBuffer &);
+template bool tryReadFloatTextPrecise<BFloat16>(BFloat16 &, ReadBuffer &);
 template bool tryReadFloatTextPrecise<Float32>(Float32 &, ReadBuffer &);
 template bool tryReadFloatTextPrecise<Float64>(Float64 &, ReadBuffer &);
 
+template void readFloatTextFast<BFloat16>(BFloat16 &, ReadBuffer &);
 template void readFloatTextFast<Float32>(Float32 &, ReadBuffer &);
 template void readFloatTextFast<Float64>(Float64 &, ReadBuffer &);
+template bool tryReadFloatTextFast<BFloat16>(BFloat16 &, ReadBuffer &);
 template bool tryReadFloatTextFast<Float32>(Float32 &, ReadBuffer &);
 template bool tryReadFloatTextFast<Float64>(Float64 &, ReadBuffer &);
 
+template void readFloatTextSimple<BFloat16>(BFloat16 &, ReadBuffer &);
 template void readFloatTextSimple<Float32>(Float32 &, ReadBuffer &);
 template void readFloatTextSimple<Float64>(Float64 &, ReadBuffer &);
+template bool tryReadFloatTextSimple<BFloat16>(BFloat16 &, ReadBuffer &);
 template bool tryReadFloatTextSimple<Float32>(Float32 &, ReadBuffer &);
 template bool tryReadFloatTextSimple<Float64>(Float64 &, ReadBuffer &);
 
+template void readFloatText<BFloat16>(BFloat16 &, ReadBuffer &);
 template void readFloatText<Float32>(Float32 &, ReadBuffer &);
 template void readFloatText<Float64>(Float64 &, ReadBuffer &);
+template bool tryReadFloatText<BFloat16>(BFloat16 &, ReadBuffer &);
 template bool tryReadFloatText<Float32>(Float32 &, ReadBuffer &);
 template bool tryReadFloatText<Float64>(Float64 &, ReadBuffer &);
 
+template bool tryReadFloatTextNoExponent<BFloat16>(BFloat16 &, ReadBuffer &);
 template bool tryReadFloatTextNoExponent<Float32>(Float32 &, ReadBuffer &);
 template bool tryReadFloatTextNoExponent<Float64>(Float64 &, ReadBuffer &);
 
diff --git a/src/IO/readFloatText.h b/src/IO/readFloatText.h
index c2fec9d4b0b..a7fd6058dd9 100644
--- a/src/IO/readFloatText.h
+++ b/src/IO/readFloatText.h
@@ -222,7 +222,6 @@ ReturnType readFloatTextPreciseImpl(T & x, ReadBuffer & buf)
             break;
         }
 
-
         char tmp_buf[MAX_LENGTH];
         int num_copied_chars = 0;
 
@@ -597,22 +596,85 @@ ReturnType readFloatTextSimpleImpl(T & x, ReadBuffer & buf)
     return ReturnType(true);
 }
 
-template <typename T> void readFloatTextPrecise(T & x, ReadBuffer & in) { readFloatTextPreciseImpl<T, void>(x, in); }
-template <typename T> bool tryReadFloatTextPrecise(T & x, ReadBuffer & in) { return readFloatTextPreciseImpl<T, bool>(x, in); }
+template <typename T> void readFloatTextPrecise(T & x, ReadBuffer & in)
+{
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        readFloatTextPreciseImpl<Float32, void>(tmp, in);
+        x = BFloat16(tmp);
+    }
+    else
+        readFloatTextPreciseImpl<T, void>(x, in);
+}
+
+template <typename T> bool tryReadFloatTextPrecise(T & x, ReadBuffer & in)
+{
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        bool res = readFloatTextPreciseImpl<Float32, bool>(tmp, in);
+        if (res)
+            x = BFloat16(tmp);
+        return res;
+    }
+    else
+        return readFloatTextPreciseImpl<T, bool>(x, in);
+}
 
 template <typename T> void readFloatTextFast(T & x, ReadBuffer & in)
 {
     bool has_fractional;
-    readFloatTextFastImpl<T, void>(x, in, has_fractional);
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        readFloatTextFastImpl<Float32, void>(tmp, in, has_fractional);
+        x = BFloat16(tmp);
+    }
+    else
+        readFloatTextFastImpl<T, void>(x, in, has_fractional);
 }
+
 template <typename T> bool tryReadFloatTextFast(T & x, ReadBuffer & in)
 {
     bool has_fractional;
-    return readFloatTextFastImpl<T, bool>(x, in, has_fractional);
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        bool res = readFloatTextFastImpl<Float32, bool>(tmp, in, has_fractional);
+        if (res)
+            x = BFloat16(tmp);
+        return res;
+    }
+    else
+        return readFloatTextFastImpl<T, bool>(x, in, has_fractional);
 }
 
-template <typename T> void readFloatTextSimple(T & x, ReadBuffer & in) { readFloatTextSimpleImpl<T, void>(x, in); }
-template <typename T> bool tryReadFloatTextSimple(T & x, ReadBuffer & in) { return readFloatTextSimpleImpl<T, bool>(x, in); }
+template <typename T> void readFloatTextSimple(T & x, ReadBuffer & in)
+{
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        readFloatTextSimpleImpl<Float32, void>(tmp, in);
+        x = BFloat16(tmp);
+    }
+    else
+        readFloatTextSimpleImpl<T, void>(x, in);
+}
+
+template <typename T> bool tryReadFloatTextSimple(T & x, ReadBuffer & in)
+{
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        bool res = readFloatTextSimpleImpl<Float32, bool>(tmp, in);
+        if (res)
+            x = BFloat16(tmp);
+        return res;
+    }
+    else
+        return readFloatTextSimpleImpl<T, bool>(x, in);
+}
 
 
 /// Implementation that is selected as default.
@@ -624,18 +686,47 @@ template <typename T> bool tryReadFloatText(T & x, ReadBuffer & in) { return try
 template <typename T> bool tryReadFloatTextNoExponent(T & x, ReadBuffer & in)
 {
     bool has_fractional;
-    return readFloatTextFastImpl<T, bool, false>(x, in, has_fractional);
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        bool res = readFloatTextFastImpl<Float32, bool, false>(tmp, in, has_fractional);
+        if (res)
+            x = BFloat16(tmp);
+        return res;
+
+    }
+    else
+        return readFloatTextFastImpl<T, bool, false>(x, in, has_fractional);
 }
 
 /// With a @has_fractional flag
 /// Used for input_format_try_infer_integers
 template <typename T> bool tryReadFloatTextExt(T & x, ReadBuffer & in, bool & has_fractional)
 {
-    return readFloatTextFastImpl<T, bool>(x, in, has_fractional);
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        bool res = readFloatTextFastImpl<Float32, bool>(tmp, in, has_fractional);
+        if (res)
+            x = BFloat16(tmp);
+        return res;
+    }
+    else
+        return readFloatTextFastImpl<T, bool>(x, in, has_fractional);
 }
+
 template <typename T> bool tryReadFloatTextExtNoExponent(T & x, ReadBuffer & in, bool & has_fractional)
 {
-    return readFloatTextFastImpl<T, bool, false>(x, in, has_fractional);
+    if constexpr (std::is_same_v<T, BFloat16>)
+    {
+        Float32 tmp;
+        bool res = readFloatTextFastImpl<Float32, bool, false>(tmp, in, has_fractional);
+        if (res)
+            x = BFloat16(tmp);
+        return res;
+    }
+    else
+        return readFloatTextFastImpl<T, bool, false>(x, in, has_fractional);
 }
 
 }
diff --git a/src/Interpreters/AggregatedData.h b/src/Interpreters/AggregatedData.h
index 4b581c682ca..1a804afd33e 100644
--- a/src/Interpreters/AggregatedData.h
+++ b/src/Interpreters/AggregatedData.h
@@ -1,12 +1,16 @@
 #pragma once
+
 #include <AggregateFunctions/IAggregateFunction_fwd.h>
 
 #include <Common/HashTable/FixedHashMap.h>
 #include <Common/HashTable/StringHashMap.h>
 #include <Common/HashTable/TwoLevelHashMap.h>
 #include <Common/HashTable/TwoLevelStringHashMap.h>
+
+
 namespace DB
 {
+
 /** Different data structures that can be used for aggregation
   * For efficiency, the aggregation data itself is put into the pool.
   * Data and pool ownership (states of aggregate functions)
diff --git a/src/Interpreters/RowRefs.cpp b/src/Interpreters/RowRefs.cpp
index 1b397ab56ef..26524218b2b 100644
--- a/src/Interpreters/RowRefs.cpp
+++ b/src/Interpreters/RowRefs.cpp
@@ -183,7 +183,7 @@ private:
         if (sorted.load(std::memory_order_relaxed))
             return;
 
-        if constexpr (std::is_arithmetic_v<TKey> && !std::is_floating_point_v<TKey>)
+        if constexpr (std::is_arithmetic_v<TKey> && !is_floating_point<TKey>)
         {
             if (likely(entries.size() > 256))
             {
diff --git a/src/Interpreters/parseColumnsListForTableFunction.cpp b/src/Interpreters/parseColumnsListForTableFunction.cpp
index d8f526db592..d3bf6f860f3 100644
--- a/src/Interpreters/parseColumnsListForTableFunction.cpp
+++ b/src/Interpreters/parseColumnsListForTableFunction.cpp
@@ -20,6 +20,7 @@ namespace Setting
     extern const SettingsBool allow_experimental_json_type;
     extern const SettingsBool allow_experimental_object_type;
     extern const SettingsBool allow_experimental_variant_type;
+    extern const SettingsBool allow_experimental_bfloat16_type;
     extern const SettingsBool allow_suspicious_fixed_string_types;
     extern const SettingsBool allow_suspicious_low_cardinality_types;
     extern const SettingsBool allow_suspicious_variant_types;
@@ -42,6 +43,7 @@ DataTypeValidationSettings::DataTypeValidationSettings(const DB::Settings & sett
     , allow_experimental_object_type(settings[Setting::allow_experimental_object_type])
     , allow_suspicious_fixed_string_types(settings[Setting::allow_suspicious_fixed_string_types])
     , allow_experimental_variant_type(settings[Setting::allow_experimental_variant_type])
+    , allow_experimental_bfloat16_type(settings[Setting::allow_experimental_bfloat16_type])
     , allow_suspicious_variant_types(settings[Setting::allow_suspicious_variant_types])
     , validate_nested_types(settings[Setting::validate_experimental_and_suspicious_types_inside_nested_types])
     , allow_experimental_dynamic_type(settings[Setting::allow_experimental_dynamic_type])
@@ -105,6 +107,18 @@ void validateDataType(const DataTypePtr & type_to_check, const DataTypeValidatio
             }
         }
 
+        if (!settings.allow_experimental_bfloat16_type)
+        {
+            if (WhichDataType(data_type).isBFloat16())
+            {
+                throw Exception(
+                    ErrorCodes::ILLEGAL_COLUMN,
+                    "Cannot create column with type '{}' because experimental BFloat16 type is not allowed. "
+                    "Set setting allow_experimental_bfloat16_type = 1 in order to allow it",
+                    data_type.getName());
+            }
+        }
+
         if (!settings.allow_suspicious_variant_types)
         {
             if (const auto * variant_type = typeid_cast<const DataTypeVariant *>(&data_type))
diff --git a/src/Interpreters/parseColumnsListForTableFunction.h b/src/Interpreters/parseColumnsListForTableFunction.h
index 6e00492c0ad..39b9f092d89 100644
--- a/src/Interpreters/parseColumnsListForTableFunction.h
+++ b/src/Interpreters/parseColumnsListForTableFunction.h
@@ -20,6 +20,7 @@ struct DataTypeValidationSettings
     bool allow_experimental_object_type = true;
     bool allow_suspicious_fixed_string_types = true;
     bool allow_experimental_variant_type = true;
+    bool allow_experimental_bfloat16_type = true;
     bool allow_suspicious_variant_types = true;
     bool validate_nested_types = true;
     bool allow_experimental_dynamic_type = true;
diff --git a/src/Processors/Formats/Impl/Parquet/ParquetDataValuesReader.cpp b/src/Processors/Formats/Impl/Parquet/ParquetDataValuesReader.cpp
index b471989076b..4b79be98810 100644
--- a/src/Processors/Formats/Impl/Parquet/ParquetDataValuesReader.cpp
+++ b/src/Processors/Formats/Impl/Parquet/ParquetDataValuesReader.cpp
@@ -580,6 +580,7 @@ template class ParquetPlainValuesReader<ColumnInt32>;
 template class ParquetPlainValuesReader<ColumnUInt32>;
 template class ParquetPlainValuesReader<ColumnInt64>;
 template class ParquetPlainValuesReader<ColumnUInt64>;
+template class ParquetPlainValuesReader<ColumnBFloat16>;
 template class ParquetPlainValuesReader<ColumnFloat32>;
 template class ParquetPlainValuesReader<ColumnFloat64>;
 template class ParquetPlainValuesReader<ColumnDecimal<Decimal32>>;
@@ -602,6 +603,7 @@ template class ParquetRleDictReader<ColumnInt32>;
 template class ParquetRleDictReader<ColumnUInt32>;
 template class ParquetRleDictReader<ColumnInt64>;
 template class ParquetRleDictReader<ColumnUInt64>;
+template class ParquetRleDictReader<ColumnBFloat16>;
 template class ParquetRleDictReader<ColumnFloat32>;
 template class ParquetRleDictReader<ColumnFloat64>;
 template class ParquetRleDictReader<ColumnDecimal<Decimal32>>;
diff --git a/src/Processors/Formats/Impl/Parquet/ParquetLeafColReader.cpp b/src/Processors/Formats/Impl/Parquet/ParquetLeafColReader.cpp
index c3c7db510ed..328dd37107e 100644
--- a/src/Processors/Formats/Impl/Parquet/ParquetLeafColReader.cpp
+++ b/src/Processors/Formats/Impl/Parquet/ParquetLeafColReader.cpp
@@ -644,6 +644,7 @@ template class ParquetLeafColReader<ColumnInt32>;
 template class ParquetLeafColReader<ColumnUInt32>;
 template class ParquetLeafColReader<ColumnInt64>;
 template class ParquetLeafColReader<ColumnUInt64>;
+template class ParquetLeafColReader<ColumnBFloat16>;
 template class ParquetLeafColReader<ColumnFloat32>;
 template class ParquetLeafColReader<ColumnFloat64>;
 template class ParquetLeafColReader<ColumnString>;
diff --git a/tests/performance/avg_weighted.xml b/tests/performance/avg_weighted.xml
index edf3c19fdfa..ec1b7aae5c2 100644
--- a/tests/performance/avg_weighted.xml
+++ b/tests/performance/avg_weighted.xml
@@ -1,6 +1,5 @@
 <test>
     <settings>
-        <allow_experimental_bigint_types>1</allow_experimental_bigint_types>
         <max_threads>1</max_threads>
         <max_insert_threads>8</max_insert_threads>
     </settings>
diff --git a/tests/performance/reinterpret_as.xml b/tests/performance/reinterpret_as.xml
index d05ef3bb038..2e0fa0571c3 100644
--- a/tests/performance/reinterpret_as.xml
+++ b/tests/performance/reinterpret_as.xml
@@ -1,6 +1,5 @@
 <test>
     <settings>
-        <allow_experimental_bigint_types>1</allow_experimental_bigint_types>
         <max_memory_usage>15G</max_memory_usage>
     </settings>
 
diff --git a/tests/queries/0_stateless/01035_avg.sql b/tests/queries/0_stateless/01035_avg.sql
index a3cb35a80ec..0f7baddaec5 100644
--- a/tests/queries/0_stateless/01035_avg.sql
+++ b/tests/queries/0_stateless/01035_avg.sql
@@ -1,5 +1,3 @@
-SET allow_experimental_bigint_types=1;
-
 CREATE TABLE IF NOT EXISTS test_01035_avg (
     i8 Int8         DEFAULT i64,
     i16 Int16       DEFAULT i64,
diff --git a/tests/queries/0_stateless/01182_materialized_view_different_structure.sql b/tests/queries/0_stateless/01182_materialized_view_different_structure.sql
index 485f9985974..7e41172bd0c 100644
--- a/tests/queries/0_stateless/01182_materialized_view_different_structure.sql
+++ b/tests/queries/0_stateless/01182_materialized_view_different_structure.sql
@@ -20,7 +20,6 @@ SELECT sum(value) FROM (SELECT number, sum(number) AS value FROM (SELECT *, toDe
 CREATE TABLE src (n UInt64, s FixedString(16)) ENGINE=Memory;
 CREATE TABLE dst (n UInt8, s String) ENGINE = Memory;
 CREATE MATERIALIZED VIEW mv TO dst (n String) AS SELECT * FROM src;
-SET allow_experimental_bigint_types=1;
 CREATE TABLE dist (n Int128) ENGINE=Distributed(test_cluster_two_shards, currentDatabase(), mv);
 
 INSERT INTO src SELECT number, toString(number) FROM numbers(1000);
diff --git a/tests/queries/0_stateless/01440_big_int_exotic_casts.sql b/tests/queries/0_stateless/01440_big_int_exotic_casts.sql
index 42fde9da01b..f411af897e8 100644
--- a/tests/queries/0_stateless/01440_big_int_exotic_casts.sql
+++ b/tests/queries/0_stateless/01440_big_int_exotic_casts.sql
@@ -32,8 +32,6 @@ SELECT number y, toInt128(number) - y, toInt256(number) - y, toUInt256(number) -
 SELECT -number y, toInt128(number) + y, toInt256(number) + y, toUInt256(number) + y FROM numbers_mt(10) ORDER BY number;
 
 
-SET allow_experimental_bigint_types = 1;
-
 DROP TABLE IF EXISTS t;
 CREATE TABLE t (x UInt64, i256 Int256, u256 UInt256, d256 Decimal256(2)) ENGINE = Memory;
 
diff --git a/tests/queries/0_stateless/01554_bloom_filter_index_big_integer_uuid.sql b/tests/queries/0_stateless/01554_bloom_filter_index_big_integer_uuid.sql
index 3472f41092d..f82fe39f439 100644
--- a/tests/queries/0_stateless/01554_bloom_filter_index_big_integer_uuid.sql
+++ b/tests/queries/0_stateless/01554_bloom_filter_index_big_integer_uuid.sql
@@ -1,5 +1,3 @@
-SET allow_experimental_bigint_types = 1;
-
 CREATE TABLE 01154_test (x Int128, INDEX ix_x x TYPE bloom_filter(0.01) GRANULARITY 1) ENGINE = MergeTree() ORDER BY x SETTINGS index_granularity=8192;
 INSERT INTO 01154_test VALUES (1), (2), (3);
 SELECT x FROM 01154_test WHERE x = 1;
diff --git a/tests/queries/0_stateless/01622_byte_size.sql b/tests/queries/0_stateless/01622_byte_size.sql
index 9f9de4e58e9..f73011f4151 100644
--- a/tests/queries/0_stateless/01622_byte_size.sql
+++ b/tests/queries/0_stateless/01622_byte_size.sql
@@ -4,8 +4,6 @@
 select '';
 select '# byteSize';
 
-set allow_experimental_bigint_types = 1;
-
 -- numbers #0 --
 select '';
 select 'byteSize for numbers #0';
diff --git a/tests/queries/0_stateless/01721_dictionary_decimal_p_s.sql b/tests/queries/0_stateless/01721_dictionary_decimal_p_s.sql
index 272bd2d7104..57483430cc0 100644
--- a/tests/queries/0_stateless/01721_dictionary_decimal_p_s.sql
+++ b/tests/queries/0_stateless/01721_dictionary_decimal_p_s.sql
@@ -1,6 +1,5 @@
 -- Tags: no-parallel
 
-set allow_experimental_bigint_types=1;
 drop database if exists db_01721;
 drop table if exists db_01721.table_decimal_dict;
 drop dictionary if exists db_01721.decimal_dict;
@@ -77,4 +76,3 @@ SELECT dictGet('db_01721.decimal_dict', 'Decimal32_', toUInt64(5000)),
 drop table if exists table_decimal_dict;
 drop dictionary if exists cache_dict;
 drop database if exists db_01721;
-
diff --git a/tests/queries/0_stateless/01804_dictionary_decimal256_type.sql b/tests/queries/0_stateless/01804_dictionary_decimal256_type.sql
index 08a8d0feb27..32b029442b9 100644
--- a/tests/queries/0_stateless/01804_dictionary_decimal256_type.sql
+++ b/tests/queries/0_stateless/01804_dictionary_decimal256_type.sql
@@ -1,7 +1,5 @@
 -- Tags: no-parallel
 
-SET allow_experimental_bigint_types = 1;
-
 DROP TABLE IF EXISTS dictionary_decimal_source_table;
 CREATE TABLE dictionary_decimal_source_table
 (
diff --git a/tests/queries/0_stateless/01875_ssd_cache_dictionary_decimal256_type.sh b/tests/queries/0_stateless/01875_ssd_cache_dictionary_decimal256_type.sh
index 1294ba53e82..2a24a931696 100755
--- a/tests/queries/0_stateless/01875_ssd_cache_dictionary_decimal256_type.sh
+++ b/tests/queries/0_stateless/01875_ssd_cache_dictionary_decimal256_type.sh
@@ -6,8 +6,6 @@ CURDIR=$(cd "$(dirname "${BASH_SOURCE[0]}")" && pwd)
 . "$CURDIR"/../shell_config.sh
 
 $CLICKHOUSE_CLIENT --query="
-    SET allow_experimental_bigint_types = 1;
-
     DROP TABLE IF EXISTS dictionary_decimal_source_table;
     CREATE TABLE dictionary_decimal_source_table
     (
diff --git a/tests/queries/0_stateless/03269_bf16.reference b/tests/queries/0_stateless/03269_bf16.reference
new file mode 100644
index 00000000000..896cc307623
--- /dev/null
+++ b/tests/queries/0_stateless/03269_bf16.reference
@@ -0,0 +1,47 @@
+1	-1	1.09375	-1.09375	1	-1	1.09375	-1.09375	18446744000000000000	-0	inf	-inf	nan
+1.09375	1.09375	1.09375	1
+1	1	0	1	1
+0	2.1875	1.1962891	1	Float32	Float32	Float32	Float64
+-0.006250000000000089	2.19375	1.203125	1.0057142857142858	Float64	Float64	Float64	Float64
+0	0	1	0
+1000	1000	1	0
+2000	2000	1	0
+3000	2992	0	8
+4000	4000	1	0
+5000	4992	0	8
+6000	5984	0	16
+7000	6976	0	24
+8000	8000	1	0
+9000	8960	0	40
+49995000	49855104	4999.5	4985.5104	0	0	9999	9984	10000	925	10000	925
+0	0	1	0
+1000	1000	1	0
+2000	2000	1	0
+3000	2992	0	8
+4000	4000	1	0
+5000	4992	0	8
+6000	5984	0	16
+7000	6976	0	24
+8000	8000	1	0
+9000	8960	0	40
+49995000	49855104	4999.5	4985.5104	0	0	9999	9984	10000	925	10000	925
+Row 1:
+──────
+a32:                        [0,0.5,1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8,8.5,9,9.5,10,10.5,11,11.5,12,12.5,13,13.5,14,14.5,15,15.5,16,16.5,17,17.5,18,18.5,19,19.5,20,20.5,21,21.5,22,22.5,23,23.5,24,24.5,25,25.5,26,26.5,27,27.5,28,28.5,29,29.5,30,30.5,31,31.5,32,32.5,33,33.5,34,34.5,35,35.5,36,36.5,37,37.5,38,38.5,39,39.5,40,40.5,41,41.5,42,42.5,43,43.5,44,44.5,45,45.5,46,46.5,47,47.5,48,48.5,49,49.5,50,50.5,51,51.5,52,52.5,53,53.5,54,54.5,55,55.5,56,56.5,57,57.5,58,58.5,59,59.5,60,60.5,61,61.5,62,62.5,63,63.5,64,64.5,65,65.5,66,66.5,67,67.5,68,68.5,69,69.5,70,70.5,71,71.5,72,72.5,73,73.5,74,74.5,75,75.5,76,76.5,77,77.5,78,78.5,79,79.5,80,80.5,81,81.5,82,82.5,83,83.5,84,84.5,85,85.5,86,86.5,87,87.5,88,88.5,89,89.5,90,90.5,91,91.5,92,92.5,93,93.5,94,94.5,95,95.5,96,96.5,97,97.5,98,98.5,99,99.5,100,100.5,101,101.5,102,102.5,103,103.5,104,104.5,105,105.5,106,106.5,107,107.5,108,108.5,109,109.5,110,110.5,111,111.5,112,112.5,113,113.5,114,114.5,115,115.5,116,116.5,117,117.5,118,118.5,119,119.5,120,120.5,121,121.5,122,122.5,123,123.5,124,124.5,125,125.5,126,126.5,127,127.5,128,128.5,129,129.5,130,130.5,131,131.5,132,132.5,133,133.5,134,134.5,135,135.5,136,136.5,137,137.5,138,138.5,139,139.5,140,140.5,141,141.5,142,142.5,143,143.5,144,144.5,145,145.5,146,146.5,147,147.5,148,148.5,149,149.5,150,150.5,151,151.5,152,152.5,153,153.5,154,154.5,155,155.5,156,156.5,157,157.5,158,158.5,159,159.5,160,160.5,161,161.5,162,162.5,163,163.5,164,164.5,165,165.5,166,166.5,167,167.5,168,168.5,169,169.5,170,170.5,171,171.5,172,172.5,173,173.5,174,174.5,175,175.5,176,176.5,177,177.5,178,178.5,179,179.5,180,180.5,181,181.5,182,182.5,183,183.5,184,184.5,185,185.5,186,186.5,187,187.5,188,188.5,189,189.5,190,190.5,191,191.5]
+a16:                        [0,0.5,1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8,8.5,9,9.5,10,10.5,11,11.5,12,12.5,13,13.5,14,14.5,15,15.5,16,16.5,17,17.5,18,18.5,19,19.5,20,20.5,21,21.5,22,22.5,23,23.5,24,24.5,25,25.5,26,26.5,27,27.5,28,28.5,29,29.5,30,30.5,31,31.5,32,32.5,33,33.5,34,34.5,35,35.5,36,36.5,37,37.5,38,38.5,39,39.5,40,40.5,41,41.5,42,42.5,43,43.5,44,44.5,45,45.5,46,46.5,47,47.5,48,48.5,49,49.5,50,50.5,51,51.5,52,52.5,53,53.5,54,54.5,55,55.5,56,56.5,57,57.5,58,58.5,59,59.5,60,60.5,61,61.5,62,62.5,63,63.5,64,64.5,65,65.5,66,66.5,67,67.5,68,68.5,69,69.5,70,70.5,71,71.5,72,72.5,73,73.5,74,74.5,75,75.5,76,76.5,77,77.5,78,78.5,79,79.5,80,80.5,81,81.5,82,82.5,83,83.5,84,84.5,85,85.5,86,86.5,87,87.5,88,88.5,89,89.5,90,90.5,91,91.5,92,92.5,93,93.5,94,94.5,95,95.5,96,96.5,97,97.5,98,98.5,99,99.5,100,100.5,101,101.5,102,102.5,103,103.5,104,104.5,105,105.5,106,106.5,107,107.5,108,108.5,109,109.5,110,110.5,111,111.5,112,112.5,113,113.5,114,114.5,115,115.5,116,116.5,117,117.5,118,118.5,119,119.5,120,120.5,121,121.5,122,122.5,123,123.5,124,124.5,125,125.5,126,126.5,127,127.5,128,128,129,129,130,130,131,131,132,132,133,133,134,134,135,135,136,136,137,137,138,138,139,139,140,140,141,141,142,142,143,143,144,144,145,145,146,146,147,147,148,148,149,149,150,150,151,151,152,152,153,153,154,154,155,155,156,156,157,157,158,158,159,159,160,160,161,161,162,162,163,163,164,164,165,165,166,166,167,167,168,168,169,169,170,170,171,171,172,172,173,173,174,174,175,175,176,176,177,177,178,178,179,179,180,180,181,181,182,182,183,183,184,184,185,185,186,186,187,187,188,188,189,189,190,190,191,191]
+a32_1:                      [1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8,8.5,9,9.5,10,10.5,11,11.5,12,12.5,13,13.5,14,14.5,15,15.5,16,16.5,17,17.5,18,18.5,19,19.5,20,20.5,21,21.5,22,22.5,23,23.5,24,24.5,25,25.5,26,26.5,27,27.5,28,28.5,29,29.5,30,30.5,31,31.5,32,32.5,33,33.5,34,34.5,35,35.5,36,36.5,37,37.5,38,38.5,39,39.5,40,40.5,41,41.5,42,42.5,43,43.5,44,44.5,45,45.5,46,46.5,47,47.5,48,48.5,49,49.5,50,50.5,51,51.5,52,52.5,53,53.5,54,54.5,55,55.5,56,56.5,57,57.5,58,58.5,59,59.5,60,60.5,61,61.5,62,62.5,63,63.5,64,64.5,65,65.5,66,66.5,67,67.5,68,68.5,69,69.5,70,70.5,71,71.5,72,72.5,73,73.5,74,74.5,75,75.5,76,76.5,77,77.5,78,78.5,79,79.5,80,80.5,81,81.5,82,82.5,83,83.5,84,84.5,85,85.5,86,86.5,87,87.5,88,88.5,89,89.5,90,90.5,91,91.5,92,92.5,93,93.5,94,94.5,95,95.5,96,96.5,97,97.5,98,98.5,99,99.5,100,100.5,101,101.5,102,102.5,103,103.5,104,104.5,105,105.5,106,106.5,107,107.5,108,108.5,109,109.5,110,110.5,111,111.5,112,112.5,113,113.5,114,114.5,115,115.5,116,116.5,117,117.5,118,118.5,119,119.5,120,120.5,121,121.5,122,122.5,123,123.5,124,124.5,125,125.5,126,126.5,127,127.5,128,128.5,129,129.5,130,130.5,131,131.5,132,132.5,133,133.5,134,134.5,135,135.5,136,136.5,137,137.5,138,138.5,139,139.5,140,140.5,141,141.5,142,142.5,143,143.5,144,144.5,145,145.5,146,146.5,147,147.5,148,148.5,149,149.5,150,150.5,151,151.5,152,152.5,153,153.5,154,154.5,155,155.5,156,156.5,157,157.5,158,158.5,159,159.5,160,160.5,161,161.5,162,162.5,163,163.5,164,164.5,165,165.5,166,166.5,167,167.5,168,168.5,169,169.5,170,170.5,171,171.5,172,172.5,173,173.5,174,174.5,175,175.5,176,176.5,177,177.5,178,178.5,179,179.5,180,180.5,181,181.5,182,182.5,183,183.5,184,184.5,185,185.5,186,186.5,187,187.5,188,188.5,189,189.5,190,190.5,191,191.5,192,192.5]
+a16_1:                      [1,1.5,2,2.5,3,3.5,4,4.5,5,5.5,6,6.5,7,7.5,8,8.5,9,9.5,10,10.5,11,11.5,12,12.5,13,13.5,14,14.5,15,15.5,16,16.5,17,17.5,18,18.5,19,19.5,20,20.5,21,21.5,22,22.5,23,23.5,24,24.5,25,25.5,26,26.5,27,27.5,28,28.5,29,29.5,30,30.5,31,31.5,32,32.5,33,33.5,34,34.5,35,35.5,36,36.5,37,37.5,38,38.5,39,39.5,40,40.5,41,41.5,42,42.5,43,43.5,44,44.5,45,45.5,46,46.5,47,47.5,48,48.5,49,49.5,50,50.5,51,51.5,52,52.5,53,53.5,54,54.5,55,55.5,56,56.5,57,57.5,58,58.5,59,59.5,60,60.5,61,61.5,62,62.5,63,63.5,64,64.5,65,65.5,66,66.5,67,67.5,68,68.5,69,69.5,70,70.5,71,71.5,72,72.5,73,73.5,74,74.5,75,75.5,76,76.5,77,77.5,78,78.5,79,79.5,80,80.5,81,81.5,82,82.5,83,83.5,84,84.5,85,85.5,86,86.5,87,87.5,88,88.5,89,89.5,90,90.5,91,91.5,92,92.5,93,93.5,94,94.5,95,95.5,96,96.5,97,97.5,98,98.5,99,99.5,100,100.5,101,101.5,102,102.5,103,103.5,104,104.5,105,105.5,106,106.5,107,107.5,108,108.5,109,109.5,110,110.5,111,111.5,112,112.5,113,113.5,114,114.5,115,115.5,116,116.5,117,117.5,118,118.5,119,119.5,120,120.5,121,121.5,122,122.5,123,123.5,124,124.5,125,125.5,126,126.5,127,127.5,128,128.5,129,129,130,130,131,131,132,132,133,133,134,134,135,135,136,136,137,137,138,138,139,139,140,140,141,141,142,142,143,143,144,144,145,145,146,146,147,147,148,148,149,149,150,150,151,151,152,152,153,153,154,154,155,155,156,156,157,157,158,158,159,159,160,160,161,161,162,162,163,163,164,164,165,165,166,166,167,167,168,168,169,169,170,170,171,171,172,172,173,173,174,174,175,175,176,176,177,177,178,178,179,179,180,180,181,181,182,182,183,183,184,184,185,185,186,186,187,187,188,188,189,189,190,190,191,191,192,192]
+dotProduct(a32, a32_1):     4736944 -- 4.74 million
+dotProduct(a16, a16_1):     4726688 -- 4.73 million
+cosineDistance(a32, a32_1): 0.000010093636084174129
+cosineDistance(a16, a16_1): 0.00001010226319664298
+L2Distance(a32, a32_1):     19.595917942265423
+L2Distance(a16, a16_1):     19.595917942265423
+L1Distance(a32, a32_1):     384
+L1Distance(a16, a16_1):     384
+LinfDistance(a32, a32_1):   1
+LinfDistance(a16, a16_1):   1
+LpDistance(a32, a32_1, 5):  3.2875036590344515
+LpDistance(a16, a16_1, 5):  3.2875036590344515
+1.09375	8C3F	1000110000111111	2	16268	8C3F
+1.09375	1	1.09375	1.0859375	0
diff --git a/tests/queries/0_stateless/03269_bf16.sql b/tests/queries/0_stateless/03269_bf16.sql
new file mode 100644
index 00000000000..b332a6e3119
--- /dev/null
+++ b/tests/queries/0_stateless/03269_bf16.sql
@@ -0,0 +1,100 @@
+SET allow_experimental_bfloat16_type = 1;
+
+-- This is a smoke test, non exhaustive.
+
+-- Conversions
+
+SELECT
+    1::BFloat16,
+    -1::BFloat16,
+    1.1::BFloat16,
+    -1.1::BFloat16,
+    CAST(1 AS BFloat16),
+    CAST(-1 AS BFloat16),
+    CAST(1.1 AS BFloat16),
+    CAST(-1.1 AS BFloat16),
+    CAST(0xFFFFFFFFFFFFFFFF AS BFloat16),
+    CAST(-0.0 AS BFloat16),
+    CAST(inf AS BFloat16),
+    CAST(-inf AS BFloat16),
+    CAST(nan AS BFloat16);
+
+-- Conversions back
+
+SELECT
+    CAST(1.1::BFloat16 AS BFloat16),
+    CAST(1.1::BFloat16 AS Float32),
+    CAST(1.1::BFloat16 AS Float64),
+    CAST(1.1::BFloat16 AS Int8);
+
+-- Comparisons
+
+SELECT
+    1.1::BFloat16 = 1.1::BFloat16,
+    1.1::BFloat16 < 1.1,
+    1.1::BFloat16 > 1.1,
+    1.1::BFloat16 > 1,
+    1.1::BFloat16 = 1.09375;
+
+-- Arithmetic
+
+SELECT
+    1.1::BFloat16 - 1.1::BFloat16 AS a,
+    1.1::BFloat16 + 1.1::BFloat16 AS b,
+    1.1::BFloat16 * 1.1::BFloat16 AS c,
+    1.1::BFloat16 / 1.1::BFloat16 AS d,
+    toTypeName(a), toTypeName(b), toTypeName(c), toTypeName(d);
+
+SELECT
+    1.1::BFloat16 - 1.1 AS a,
+    1.1 + 1.1::BFloat16 AS b,
+    1.1::BFloat16 * 1.1 AS c,
+    1.1 / 1.1::BFloat16 AS d,
+    toTypeName(a), toTypeName(b), toTypeName(c), toTypeName(d);
+
+-- Tables
+
+DROP TABLE IF EXISTS t;
+CREATE TEMPORARY TABLE t (n UInt64, x BFloat16);
+INSERT INTO t SELECT number, number FROM numbers(10000);
+SELECT *, n = x, n - x FROM t WHERE n % 1000 = 0 ORDER BY n;
+
+-- Aggregate functions
+
+SELECT sum(n), sum(x), avg(n), avg(x), min(n), min(x), max(n), max(x), uniq(n), uniq(x), uniqExact(n), uniqExact(x) FROM t;
+
+-- MergeTree
+
+DROP TABLE t;
+CREATE TABLE t (n UInt64, x BFloat16) ENGINE = MergeTree ORDER BY n;
+INSERT INTO t SELECT number, number FROM numbers(10000);
+SELECT *, n = x, n - x FROM t WHERE n % 1000 = 0 ORDER BY n;
+SELECT sum(n), sum(x), avg(n), avg(x), min(n), min(x), max(n), max(x), uniq(n), uniq(x), uniqExact(n), uniqExact(x) FROM t;
+
+-- Distances
+
+WITH
+    arrayMap(x -> toFloat32(x) / 2, range(384)) AS a32,
+    arrayMap(x -> toBFloat16(x) / 2, range(384)) AS a16,
+    arrayMap(x -> x + 1, a32) AS a32_1,
+    arrayMap(x -> x + 1, a16) AS a16_1
+SELECT a32, a16, a32_1, a16_1,
+    dotProduct(a32, a32_1), dotProduct(a16, a16_1),
+    cosineDistance(a32, a32_1), cosineDistance(a16, a16_1),
+    L2Distance(a32, a32_1), L2Distance(a16, a16_1),
+    L1Distance(a32, a32_1), L1Distance(a16, a16_1),
+    LinfDistance(a32, a32_1), LinfDistance(a16, a16_1),
+    LpDistance(a32, a32_1, 5), LpDistance(a16, a16_1, 5)
+FORMAT Vertical;
+
+-- Introspection
+
+SELECT 1.1::BFloat16 AS x,
+    hex(x), bin(x),
+    byteSize(x),
+    reinterpretAsUInt16(x), hex(reinterpretAsString(x));
+
+-- Rounding (this could be not towards the nearest)
+
+SELECT 1.1::BFloat16 AS x,
+    round(x), round(x, 1), round(x, 2), round(x, -1);
diff --git a/utils/check-style/aspell-ignore/en/aspell-dict.txt b/utils/check-style/aspell-ignore/en/aspell-dict.txt
index a0d4d1d349e..fa581c705f2 100644
--- a/utils/check-style/aspell-ignore/en/aspell-dict.txt
+++ b/utils/check-style/aspell-ignore/en/aspell-dict.txt
@@ -3156,3 +3156,4 @@ znode
 znodes
 zookeeperSessionUptime
 zstd
+BFloat