diff --git a/bindings/Fortran/modules/adios2_parameters_mod.f90 b/bindings/Fortran/modules/adios2_parameters_mod.f90
index ac89b4f1ed..dc27a9c210 100644
--- a/bindings/Fortran/modules/adios2_parameters_mod.f90
+++ b/bindings/Fortran/modules/adios2_parameters_mod.f90
@@ -55,6 +55,8 @@ module adios2_parameters_mod
     integer, parameter :: adios2_mode_write = 1
     integer, parameter :: adios2_mode_read = 2
     integer, parameter :: adios2_mode_append = 3
+    integer, parameter :: adios2_mode_readRandomAccess = 6,
+
     integer, parameter :: adios2_mode_deferred = 4
     integer, parameter :: adios2_mode_sync = 5
 
diff --git a/source/adios2/core/Operator.cpp b/source/adios2/core/Operator.cpp
index ad36b92785..23c36cb64f 100644
--- a/source/adios2/core/Operator.cpp
+++ b/source/adios2/core/Operator.cpp
@@ -98,6 +98,8 @@ Dims Operator::ConvertDims(const Dims &dimensions, const DataType type,
     return ret;
 }
 
+size_t Operator::GetHeaderSize() const { return 0; }
+
 // PRIVATE
 void Operator::CheckCallbackType(const std::string type) const
 {
diff --git a/source/adios2/core/Operator.h b/source/adios2/core/Operator.h
index 0be396dc5e..63a33e749f 100644
--- a/source/adios2/core/Operator.h
+++ b/source/adios2/core/Operator.h
@@ -68,6 +68,8 @@ class Operator
                               const std::string &, const size_t, const Dims &,
                               const Dims &, const Dims &) const;
 
+    virtual size_t GetHeaderSize() const;
+
     /**
      * @param dataIn
      * @param blockStart
diff --git a/source/adios2/helper/adiosMemory.cpp b/source/adios2/helper/adiosMemory.cpp
index ad7928eaef..779cd68cef 100644
--- a/source/adios2/helper/adiosMemory.cpp
+++ b/source/adios2/helper/adiosMemory.cpp
@@ -20,6 +20,16 @@ namespace adios2
 namespace helper
 {
 
+size_t CopyMemoryWithOpHeader(const char *src, const Dims &blockCount,
+                              const DataType type, char *dest,
+                              size_t destOffset, const MemorySpace memSpace)
+{
+    const size_t sizeIn = GetTotalSize(blockCount, GetDataTypeSize(type));
+    CopyContiguousMemory(src, sizeIn, dest + destOffset,
+                         /* endianReverse */ false, memSpace);
+    return destOffset + sizeIn;
+}
+
 namespace
 {
 
diff --git a/source/adios2/helper/adiosMemory.h b/source/adios2/helper/adiosMemory.h
index 09261d565d..02a1248825 100644
--- a/source/adios2/helper/adiosMemory.h
+++ b/source/adios2/helper/adiosMemory.h
@@ -191,6 +191,10 @@ void ClipContiguousMemory(T *dest, const Dims &destStart, const Dims &destCount,
                           const bool endianReverse = false,
                           const MemorySpace memSpace = MemorySpace::Host);
 
+size_t CopyMemoryWithOpHeader(const char *src, const Dims &blockCount,
+                              const DataType type, char *dest,
+                              size_t headerSize, const MemorySpace memSpace);
+
 template <class T>
 void CopyContiguousMemory(const char *src, const size_t stride, T *dest,
                           const bool endianReverse = false,
diff --git a/source/adios2/helper/adiosMemory.inl b/source/adios2/helper/adiosMemory.inl
index fdf8654364..7854056c37 100644
--- a/source/adios2/helper/adiosMemory.inl
+++ b/source/adios2/helper/adiosMemory.inl
@@ -650,6 +650,11 @@ template <class T>
 void CopyContiguousMemory(const char *src, const size_t payloadStride, T *dest,
                           const bool endianReverse, const MemorySpace memSpace)
 {
+    if (payloadStride == 0)
+    {
+        return;
+    }
+
 #ifdef ADIOS2_HAVE_GPU_SUPPORT
     if (memSpace == MemorySpace::GPU)
     {
diff --git a/source/adios2/operator/OperatorFactory.cpp b/source/adios2/operator/OperatorFactory.cpp
index 1979aeeb02..f6a6ab7037 100644
--- a/source/adios2/operator/OperatorFactory.cpp
+++ b/source/adios2/operator/OperatorFactory.cpp
@@ -169,7 +169,7 @@ std::shared_ptr<Operator> MakeOperator(const std::string &type,
 }
 
 size_t Decompress(const char *bufferIn, const size_t sizeIn, char *dataOut,
-                  std::shared_ptr<Operator> op)
+                  MemorySpace memSpace, std::shared_ptr<Operator> op)
 {
     Operator::OperatorType compressorType;
     std::memcpy(&compressorType, bufferIn, 1);
@@ -177,7 +177,15 @@ size_t Decompress(const char *bufferIn, const size_t sizeIn, char *dataOut,
     {
         op = MakeOperator(OperatorTypeToString(compressorType), {});
     }
-    return op->InverseOperate(bufferIn, sizeIn, dataOut);
+    size_t sizeOut = op->InverseOperate(bufferIn, sizeIn, dataOut);
+    if (sizeOut == 0) // the inverse operator was not applied
+    {
+        size_t headerSize = op->GetHeaderSize();
+        sizeOut = sizeIn - headerSize;
+        helper::CopyContiguousMemory(bufferIn + headerSize, sizeOut, dataOut,
+                                     /*endianReverse*/ false, memSpace);
+    }
+    return sizeOut;
 }
 
 } // end namespace core
diff --git a/source/adios2/operator/OperatorFactory.h b/source/adios2/operator/OperatorFactory.h
index ccac962c38..ba93b35178 100644
--- a/source/adios2/operator/OperatorFactory.h
+++ b/source/adios2/operator/OperatorFactory.h
@@ -21,7 +21,7 @@ std::shared_ptr<Operator> MakeOperator(const std::string &type,
                                        const Params &parameters);
 
 size_t Decompress(const char *bufferIn, const size_t sizeIn, char *dataOut,
-                  std::shared_ptr<Operator> op = nullptr);
+                  MemorySpace memSpace, std::shared_ptr<Operator> op = nullptr);
 
 } // end namespace core
 } // end namespace adios2
diff --git a/source/adios2/operator/compress/CompressBlosc.cpp b/source/adios2/operator/compress/CompressBlosc.cpp
index 1e70bf878e..5b243ccdcf 100644
--- a/source/adios2/operator/compress/CompressBlosc.cpp
+++ b/source/adios2/operator/compress/CompressBlosc.cpp
@@ -212,17 +212,19 @@ size_t CompressBlosc::Operate(const char *dataIn, const Dims &blockStart,
         }
     }
 
+    headerSize = bufferOutOffset;
     if (useMemcpy)
     {
-        std::memcpy(bufferOut + bufferOutOffset, dataIn + inputOffset, sizeIn);
-        bufferOutOffset += sizeIn;
         headerPtr->SetNumChunks(0u);
+        bufferOutOffset = 0; // indicate that the operator was not applied
     }
 
     blosc2_destroy();
     return bufferOutOffset;
 }
 
+size_t CompressBlosc::GetHeaderSize() const { return headerSize; }
+
 size_t CompressBlosc::InverseOperate(const char *bufferIn, const size_t sizeIn,
                                      char *dataOut)
 {
@@ -230,6 +232,7 @@ size_t CompressBlosc::InverseOperate(const char *bufferIn, const size_t sizeIn,
     const uint8_t bufferVersion =
         GetParameter<uint8_t>(bufferIn, bufferInOffset);
     bufferInOffset += 2; // skip two reserved bytes
+    headerSize = bufferInOffset;
 
     if (bufferVersion == 1)
     {
@@ -263,6 +266,7 @@ size_t CompressBlosc::DecompressV1(const char *bufferIn, const size_t sizeIn,
 
     size_t bufferInOffset = 0;
     size_t sizeOut = GetParameter<size_t, size_t>(bufferIn, bufferInOffset);
+    bool isCompressed = true;
 
     m_VersionInfo =
         " Data is compressed using BLOSC Version " +
@@ -294,18 +298,26 @@ size_t CompressBlosc::DecompressV1(const char *bufferIn, const size_t sizeIn,
             DecompressOldFormat(bufferIn + bufferInOffset,
                                 sizeIn - bufferInOffset, dataOut, sizeOut);
     }
+    if (decompressedSize == 0) // the decompression was not applied
+    {
+        isCompressed = false;
+        decompressedSize = bufferDecompressedSize;
+    }
     if (decompressedSize != sizeOut)
     {
         helper::Throw<std::runtime_error>("Operator", "CompressBlosc",
                                           "DecompressV1", m_VersionInfo);
     }
+    headerSize += sizeIn - sizeOut;
+    if (!isCompressed)
+        return 0; // indicate that the inverse operator was not applied
     return sizeOut;
 }
 
 size_t CompressBlosc::DecompressChunkedFormat(const char *bufferIn,
                                               const size_t sizeIn,
                                               char *dataOut,
-                                              const size_t sizeOut) const
+                                              const size_t sizeOut)
 {
     const DataHeader *dataPtr = reinterpret_cast<const DataHeader *>(bufferIn);
     uint32_t num_chunks = dataPtr->GetNumChunks();
@@ -386,9 +398,8 @@ size_t CompressBlosc::DecompressChunkedFormat(const char *bufferIn,
     }
     else
     {
-        std::memcpy(dataOut, inputDataBuff, inputDataSize);
-        currentOutputSize = inputDataSize;
-        inputOffset += inputDataSize;
+        bufferDecompressedSize = inputDataSize;
+        return 0; // the inverse operator was not applied
     }
 
     assert(currentOutputSize == uncompressedSize);
diff --git a/source/adios2/operator/compress/CompressBlosc.h b/source/adios2/operator/compress/CompressBlosc.h
index 6af032a271..3e1d8d273b 100644
--- a/source/adios2/operator/compress/CompressBlosc.h
+++ b/source/adios2/operator/compress/CompressBlosc.h
@@ -67,12 +67,16 @@ class CompressBlosc : public Operator
 
     bool IsDataTypeValid(const DataType type) const final;
 
+    size_t GetHeaderSize() const;
+
 private:
     using bloscSize_t = int32_t;
+    size_t headerSize = 0;
+    size_t bufferDecompressedSize = 0;
 
     /** Decompress chunked data */
     size_t DecompressChunkedFormat(const char *bufferIn, const size_t sizeIn,
-                                   char *dataOut, const size_t sizeOut) const;
+                                   char *dataOut, const size_t sizeOut);
 
     /** Decompress data written before ADIOS2 supported large variables larger
      * 2GiB. */
diff --git a/source/adios2/operator/compress/CompressMGARD.cpp b/source/adios2/operator/compress/CompressMGARD.cpp
index 2371a6354d..7b22456426 100644
--- a/source/adios2/operator/compress/CompressMGARD.cpp
+++ b/source/adios2/operator/compress/CompressMGARD.cpp
@@ -103,6 +103,14 @@ size_t CompressMGARD::Operate(const char *dataIn, const Dims &blockStart,
     double s = 0.0;
     auto errorBoundType = mgard_x::error_bound_type::REL;
 
+    // input size under this bound will not compress
+    size_t thresholdSize = 100000;
+
+    auto itThreshold = m_Parameters.find("threshold");
+    if (itThreshold != m_Parameters.end())
+    {
+        thresholdSize = std::stod(itThreshold->second);
+    }
     auto itAccuracy = m_Parameters.find("accuracy");
     if (itAccuracy != m_Parameters.end())
     {
@@ -142,15 +150,26 @@ size_t CompressMGARD::Operate(const char *dataIn, const Dims &blockStart,
     // let mgard know the output buffer size
     size_t sizeOut =
         helper::GetTotalSize(blockCount, helper::GetDataTypeSize(type));
+
+    if (sizeOut < thresholdSize)
+    {
+        /* disable compression and add marker in the header*/
+        PutParameter(bufferOut, bufferOutOffset, false);
+        headerSize = bufferOutOffset;
+        return 0;
+    }
+
+    PutParameter(bufferOut, bufferOutOffset, true);
     void *compressedData = bufferOut + bufferOutOffset;
     mgard_x::compress(mgardDim, mgardType, mgardCount, tolerance, s,
                       errorBoundType, dataIn, compressedData, sizeOut, true);
-
     bufferOutOffset += sizeOut;
 
     return bufferOutOffset;
 }
 
+size_t CompressMGARD::GetHeaderSize() const { return headerSize; }
+
 size_t CompressMGARD::DecompressV1(const char *bufferIn, const size_t sizeIn,
                                    char *dataOut)
 {
@@ -175,6 +194,8 @@ size_t CompressMGARD::DecompressV1(const char *bufferIn, const size_t sizeIn,
         std::to_string(GetParameter<uint8_t>(bufferIn, bufferInOffset)) +
         ". Please make sure a compatible version is used for decompression.";
 
+    const bool isCompressed = GetParameter<bool>(bufferIn, bufferInOffset);
+
     size_t sizeOut =
         helper::GetTotalSize(blockCount, helper::GetDataTypeSize(type));
 
@@ -183,19 +204,24 @@ size_t CompressMGARD::DecompressV1(const char *bufferIn, const size_t sizeIn,
         sizeOut /= 2;
     }
 
-    try
-    {
-        void *dataOutVoid = dataOut;
-        mgard_x::decompress(bufferIn + bufferInOffset, sizeIn - bufferInOffset,
-                            dataOutVoid, true);
-    }
-    catch (...)
+    if (isCompressed)
     {
-        helper::Throw<std::runtime_error>("Operator", "CompressMGARD",
-                                          "DecompressV1", m_VersionInfo);
+        try
+        {
+            void *dataOutVoid = dataOut;
+            mgard_x::decompress(bufferIn + bufferInOffset,
+                                sizeIn - bufferInOffset, dataOutVoid, true);
+        }
+        catch (...)
+        {
+            helper::Throw<std::runtime_error>("Operator", "CompressMGARD",
+                                              "DecompressV1", m_VersionInfo);
+        }
+        return sizeOut;
     }
 
-    return sizeOut;
+    headerSize += bufferInOffset;
+    return 0;
 }
 
 size_t CompressMGARD::InverseOperate(const char *bufferIn, const size_t sizeIn,
@@ -205,6 +231,7 @@ size_t CompressMGARD::InverseOperate(const char *bufferIn, const size_t sizeIn,
     const uint8_t bufferVersion =
         GetParameter<uint8_t>(bufferIn, bufferInOffset);
     bufferInOffset += 2; // skip two reserved bytes
+    headerSize = bufferInOffset;
 
     if (bufferVersion == 1)
     {
diff --git a/source/adios2/operator/compress/CompressMGARD.h b/source/adios2/operator/compress/CompressMGARD.h
index be6bc7e56e..a6824debff 100644
--- a/source/adios2/operator/compress/CompressMGARD.h
+++ b/source/adios2/operator/compress/CompressMGARD.h
@@ -51,7 +51,10 @@ class CompressMGARD : public Operator
 
     bool IsDataTypeValid(const DataType type) const final;
 
+    size_t GetHeaderSize() const;
+
 private:
+    size_t headerSize = 0;
     /**
      * Decompress function for V1 buffer. Do NOT remove even if the buffer
      * version is updated. Data might be still in lagacy formats. This function
diff --git a/source/adios2/operator/compress/CompressMGARDPlus.cpp b/source/adios2/operator/compress/CompressMGARDPlus.cpp
index 6f48bf165a..752091c87f 100644
--- a/source/adios2/operator/compress/CompressMGARDPlus.cpp
+++ b/source/adios2/operator/compress/CompressMGARDPlus.cpp
@@ -54,6 +54,11 @@ size_t CompressMGARDPlus::Operate(const char *dataIn, const Dims &blockStart,
     CompressMGARD mgard(m_Parameters);
     size_t mgardBufferSize = mgard.Operate(dataIn, blockStart, blockCount, type,
                                            bufferOut + bufferOutOffset);
+    if (mgardBufferSize == 0)
+    {
+        headerSize += (bufferOutOffset + mgard.GetHeaderSize());
+        return 0;
+    }
 
     if (*reinterpret_cast<OperatorType *>(bufferOut + bufferOutOffset) ==
         COMPRESS_MGARD)
@@ -90,6 +95,8 @@ size_t CompressMGARDPlus::Operate(const char *dataIn, const Dims &blockStart,
     return bufferOutOffset;
 }
 
+size_t CompressMGARDPlus::GetHeaderSize() const { return headerSize; }
+
 size_t CompressMGARDPlus::DecompressV1(const char *bufferIn,
                                        const size_t sizeIn, char *dataOut)
 {
@@ -114,6 +121,7 @@ size_t CompressMGARDPlus::DecompressV1(const char *bufferIn,
     // sizeOut. Here you may want to do your magic to change the decompressed
     // data somehow to improve its accuracy :)
 
+    headerSize += (bufferInOffset + mgard.GetHeaderSize());
     return sizeOut;
 }
 
@@ -124,6 +132,7 @@ size_t CompressMGARDPlus::InverseOperate(const char *bufferIn,
     const uint8_t bufferVersion =
         GetParameter<uint8_t>(bufferIn, bufferInOffset);
     bufferInOffset += 2; // skip two reserved bytes
+    headerSize = bufferInOffset;
 
     if (bufferVersion == 1)
     {
diff --git a/source/adios2/operator/compress/CompressMGARDPlus.h b/source/adios2/operator/compress/CompressMGARDPlus.h
index 83b54c3665..33dbb55ced 100644
--- a/source/adios2/operator/compress/CompressMGARDPlus.h
+++ b/source/adios2/operator/compress/CompressMGARDPlus.h
@@ -51,7 +51,11 @@ class CompressMGARDPlus : public Operator
 
     bool IsDataTypeValid(const DataType type) const final;
 
+    size_t GetHeaderSize() const;
+
 private:
+    size_t headerSize = 0;
+
     /**
      * Decompress function for V1 buffer. Do NOT remove even if the buffer
      * version is updated. Data might be still in lagacy formats. This function
diff --git a/source/adios2/toolkit/format/bp/BPSerializer.tcc b/source/adios2/toolkit/format/bp/BPSerializer.tcc
index 3471e8179e..ff4f6f9eb3 100644
--- a/source/adios2/toolkit/format/bp/BPSerializer.tcc
+++ b/source/adios2/toolkit/format/bp/BPSerializer.tcc
@@ -405,11 +405,17 @@ void BPSerializer::PutOperationPayloadInBuffer(
     const core::Variable<T> &variable,
     const typename core::Variable<T>::BPInfo &blockInfo)
 {
-    const size_t outputSize = blockInfo.Operations[0]->Operate(
+    size_t outputSize = blockInfo.Operations[0]->Operate(
         reinterpret_cast<char *>(blockInfo.Data), blockInfo.Start,
         blockInfo.Count, variable.m_Type,
         m_Data.m_Buffer.data() + m_Data.m_Position);
 
+    if (outputSize == 0) // the operator was not applied
+        outputSize = helper::CopyMemoryWithOpHeader(
+            reinterpret_cast<char *>(blockInfo.Data), blockInfo.Count,
+            variable.m_Type, m_Data.m_Buffer.data() + m_Data.m_Position,
+            blockInfo.Operations[0]->GetHeaderSize(), blockInfo.MemSpace);
+
     m_Data.m_Position += outputSize;
     m_Data.m_AbsolutePosition += outputSize;
 
diff --git a/source/adios2/toolkit/format/bp/bp3/BP3Deserializer.tcc b/source/adios2/toolkit/format/bp/bp3/BP3Deserializer.tcc
index c07829db6b..afbda7c8c8 100644
--- a/source/adios2/toolkit/format/bp/bp3/BP3Deserializer.tcc
+++ b/source/adios2/toolkit/format/bp/bp3/BP3Deserializer.tcc
@@ -566,7 +566,7 @@ void BP3Deserializer::PostDataRead(
                 }
             }
             core::Decompress(postOpData, blockOperationInfo.PayloadSize,
-                             preOpData, op);
+                             preOpData, blockInfo.MemSpace, op);
 
             // clip block to match selection
             helper::ClipVector(m_ThreadBuffers[threadID][0],
diff --git a/source/adios2/toolkit/format/bp/bp4/BP4Deserializer.tcc b/source/adios2/toolkit/format/bp/bp4/BP4Deserializer.tcc
index cace837be8..c1494afbfb 100644
--- a/source/adios2/toolkit/format/bp/bp4/BP4Deserializer.tcc
+++ b/source/adios2/toolkit/format/bp/bp4/BP4Deserializer.tcc
@@ -570,7 +570,7 @@ void BP4Deserializer::PostDataRead(
                 }
             }
             core::Decompress(postOpData, blockOperationInfo.PayloadSize,
-                             preOpData, op);
+                             preOpData, blockInfo.MemSpace, op);
 
             // clip block to match selection
             helper::ClipVector(m_ThreadBuffers[threadID][0],
diff --git a/source/adios2/toolkit/format/bp5/BP5Deserializer.cpp b/source/adios2/toolkit/format/bp5/BP5Deserializer.cpp
index 8519c50208..857e366d14 100644
--- a/source/adios2/toolkit/format/bp5/BP5Deserializer.cpp
+++ b/source/adios2/toolkit/format/bp5/BP5Deserializer.cpp
@@ -1773,7 +1773,7 @@ void BP5Deserializer::FinalizeGet(const ReadRequest &Read, const bool freeAddr)
             core::Decompress(IncomingData,
                              ((MetaArrayRecOperator *)writer_meta_base)
                                  ->DataBlockSize[Read.BlockID],
-                             decompressBuffer.data());
+                             decompressBuffer.data(), Req.MemSpace);
         }
         IncomingData = decompressBuffer.data();
         VirtualIncomingData = IncomingData;
diff --git a/source/adios2/toolkit/format/bp5/BP5Serializer.cpp b/source/adios2/toolkit/format/bp5/BP5Serializer.cpp
index 32f8e6508c..1100894b62 100644
--- a/source/adios2/toolkit/format/bp5/BP5Serializer.cpp
+++ b/source/adios2/toolkit/format/bp5/BP5Serializer.cpp
@@ -736,6 +736,12 @@ void BP5Serializer::Marshal(void *Variable, const char *Name,
             CompressedSize = VB->m_Operations[0]->Operate(
                 (const char *)Data, tmpOffsets, tmpCount, (DataType)Rec->Type,
                 CompressedData);
+            // if the operator was not applied
+            if (CompressedSize == 0)
+                CompressedSize = helper::CopyMemoryWithOpHeader(
+                    (const char *)Data, tmpCount, (DataType)Rec->Type,
+                    CompressedData, VB->m_Operations[0]->GetHeaderSize(),
+                    MemSpace);
             CurDataBuffer->DownsizeLastAlloc(AllocSize, CompressedSize);
         }
         else if (Span == nullptr)
diff --git a/source/adios2/toolkit/format/dataman/DataManSerializer.tcc b/source/adios2/toolkit/format/dataman/DataManSerializer.tcc
index 884c8e451b..535e0a00cc 100644
--- a/source/adios2/toolkit/format/dataman/DataManSerializer.tcc
+++ b/source/adios2/toolkit/format/dataman/DataManSerializer.tcc
@@ -143,6 +143,11 @@ void DataManSerializer::PutData(
         datasize = ops[0]->Operate(reinterpret_cast<const char *>(inputData),
                                    varStart, varCount, helper::GetDataType<T>(),
                                    m_CompressBuffer.data());
+        if (datasize == 0) // operator was not applied
+            datasize = helper::CopyMemoryWithOpHeader(
+                reinterpret_cast<const char *>(inputData), varCount,
+                helper::GetDataType<T>(), m_CompressBuffer.data(),
+                ops[0]->GetHeaderSize(), MemorySpace::Host);
         compressed = true;
     }
 
@@ -247,7 +252,7 @@ int DataManSerializer::GetData(T *outputData, const std::string &varName,
                 decompressBuffer.reserve(
                     helper::GetTotalSize(j.count, sizeof(T)));
                 core::Decompress(j.buffer->data() + j.position, j.size,
-                                 decompressBuffer.data());
+                                 decompressBuffer.data(), MemorySpace::Host);
                 decompressed = true;
                 input_data = decompressBuffer.data();
             }
diff --git a/testing/adios2/engine/bp/operations/TestBPWriteReadMGARDCuda.cpp b/testing/adios2/engine/bp/operations/TestBPWriteReadMGARDCuda.cpp
index b164132e8b..34ab31a230 100644
--- a/testing/adios2/engine/bp/operations/TestBPWriteReadMGARDCuda.cpp
+++ b/testing/adios2/engine/bp/operations/TestBPWriteReadMGARDCuda.cpp
@@ -70,7 +70,7 @@ void MGARDAccuracy2D(const std::string tolerance)
             bpWriter.BeginStep();
             cudaMemcpy(gpu64s, r64s.data(), Nx * Ny * sizeof(double),
                        cudaMemcpyHostToDevice);
-            var_r64.SetMemorySpace(adios2::MemorySpace::CUDA);
+            var_r64.SetMemorySpace(adios2::MemorySpace::GPU);
             bpWriter.Put<double>("r64", gpu64s);
             bpWriter.EndStep();
         }
@@ -140,6 +140,135 @@ void MGARDAccuracy2D(const std::string tolerance)
     }
 }
 
+void MGARDAccuracySmall(const std::string tolerance)
+{
+    // Each process would write a 1x8 array and all processes would
+    // form a mpiSize * Nx 1D array
+    const std::string fname("BPWRMGARD1D_" + tolerance + ".bp");
+
+    int mpiRank = 0, mpiSize = 1;
+    // Number of rows
+    const size_t Nx = 100;
+
+    // Number of steps
+    const size_t NSteps = 1;
+
+    std::vector<float> r32s(Nx);
+    std::iota(r32s.begin(), r32s.end(), 0.f);
+
+#if ADIOS2_USE_MPI
+    MPI_Comm_rank(MPI_COMM_WORLD, &mpiRank);
+    MPI_Comm_size(MPI_COMM_WORLD, &mpiSize);
+#endif
+
+#if ADIOS2_USE_MPI
+    adios2::ADIOS adios(MPI_COMM_WORLD);
+#else
+    adios2::ADIOS adios;
+#endif
+    {
+        adios2::IO io = adios.DeclareIO("TestIO");
+
+        if (!engineName.empty())
+        {
+            io.SetEngine(engineName);
+        }
+
+        const adios2::Dims shape{static_cast<size_t>(Nx * mpiSize)};
+        const adios2::Dims start{static_cast<size_t>(Nx * mpiRank)};
+        const adios2::Dims count{Nx};
+
+        auto var_r32 = io.DefineVariable<float>("r32", shape, start, count,
+                                                adios2::ConstantDims);
+
+        // add operations
+        adios2::Operator mgardOp =
+            adios.DefineOperator("mgardCompressor", adios2::ops::LossyMGARD);
+
+        var_r32.AddOperation(mgardOp,
+                             {{adios2::ops::mgard::key::tolerance, tolerance},
+                              {adios2::ops::mgard::key::s, "inf"}});
+
+        adios2::Engine bpWriter = io.Open(fname, adios2::Mode::Write);
+
+        float *gpu32s = nullptr;
+        cudaMalloc(&gpu32s, Nx * sizeof(float));
+        for (size_t step = 0; step < NSteps; ++step)
+        {
+            bpWriter.BeginStep();
+            cudaMemcpy(gpu32s, r32s.data(), Nx * sizeof(float),
+                       cudaMemcpyHostToDevice);
+            bpWriter.Put<float>("r32", gpu32s);
+            bpWriter.EndStep();
+        }
+
+        bpWriter.Close();
+    }
+
+    {
+        adios2::IO io = adios.DeclareIO("ReadIO");
+
+        if (!engineName.empty())
+        {
+            io.SetEngine(engineName);
+        }
+
+        adios2::Engine bpReader = io.Open(fname, adios2::Mode::Read);
+
+        unsigned int t = 0;
+        std::vector<float> decompressedR32s(Nx);
+
+        while (bpReader.BeginStep() == adios2::StepStatus::OK)
+        {
+            auto var_r32 = io.InquireVariable<float>("r32");
+            EXPECT_TRUE(var_r32);
+            ASSERT_EQ(var_r32.ShapeID(), adios2::ShapeID::GlobalArray);
+            ASSERT_EQ(var_r32.Steps(), NSteps);
+            ASSERT_EQ(var_r32.Shape()[0], mpiSize * Nx);
+
+            const adios2::Dims start{mpiRank * Nx};
+            const adios2::Dims count{Nx};
+            const adios2::Box<adios2::Dims> sel(start, count);
+            var_r32.SetSelection(sel);
+
+            float *gpu32s = nullptr;
+            cudaMalloc(&gpu32s, Nx * sizeof(float));
+            bpReader.Get(var_r32, gpu32s);
+            bpReader.EndStep();
+            cudaMemcpy(decompressedR32s.data(), gpu32s, Nx * sizeof(float),
+                       cudaMemcpyDeviceToHost);
+
+            double maxDiff = 0, relativeMaxDiff = 0;
+
+            for (size_t i = 0; i < Nx; ++i)
+            {
+                std::stringstream ss;
+                ss << "t=" << t << " i=" << i << " rank=" << mpiRank;
+                std::string msg = ss.str();
+
+                double diff = std::abs(r32s[i] - decompressedR32s[i]);
+
+                if (diff > maxDiff)
+                {
+                    maxDiff = diff;
+                }
+            }
+            ++t;
+
+            auto r32s_Max = std::max_element(r32s.begin(), r32s.end());
+            relativeMaxDiff = maxDiff / *r32s_Max;
+
+            ASSERT_LT(relativeMaxDiff, std::stod(tolerance));
+            std::cout << "Relative Max Diff " << relativeMaxDiff
+                      << " tolerance " << tolerance << "\n";
+        }
+
+        EXPECT_EQ(t, NSteps);
+
+        bpReader.Close();
+    }
+}
+
 class BPWriteReadMGARD : public ::testing::TestWithParam<std::string>
 {
 public:
@@ -149,6 +278,7 @@ class BPWriteReadMGARD : public ::testing::TestWithParam<std::string>
 };
 
 TEST_P(BPWriteReadMGARD, BPWRMGARDCU2D) { MGARDAccuracy2D(GetParam()); }
+TEST_P(BPWriteReadMGARD, BPWRMGARDCU1D) { MGARDAccuracySmall(GetParam()); }
 
 INSTANTIATE_TEST_SUITE_P(MGARDAccuracy, BPWriteReadMGARD,
                          ::testing::Values("0.01", "0.001", "0.0001",