From ce91b7a868f9684c75ad260d8c719ecd1c72350c Mon Sep 17 00:00:00 2001
From: kradchen <chenjuijun@realcan.cn>
Date: Mon, 29 May 2023 11:09:09 +0800
Subject: [PATCH] Add Matrix method setBlockValue

---
 src/Matrix.cpp       | 63 ++++++++++++++++++++++++++++++++++++++++++++
 src/Matrix.h         |  2 ++
 test/Matrix_Test.cpp | 32 +++++++++++++++++++---
 3 files changed, 94 insertions(+), 3 deletions(-)

diff --git a/src/Matrix.cpp b/src/Matrix.cpp
index 8043598..aec01ea 100644
--- a/src/Matrix.cpp
+++ b/src/Matrix.cpp
@@ -616,6 +616,69 @@ namespace Aurora {
         }
 
         
+    }
+
+    bool Matrix::setBlockValue(int aDim,int aBeginIndex, int aEndIndex, double value) {
+        if(aDim>2 ){
+            std::cerr<<"block only support 1D-3D data!"<<std::endl;
+            return false;
+        }
+        //横向vector，切面为0，为1，都强制设置aDim为1来处理
+        if (isVector() && aDim == 0 && getDimSize(1)>1){
+            aDim = 1;
+        }
+        if (aBeginIndex>=getDimSize(aDim) || aBeginIndex<0){
+            std::cerr<<"block BeginIndx error!BeginIndx:"<<aBeginIndex<<std::endl;
+            return false;
+        }
+        if (aEndIndex>=getDimSize(aDim) || aEndIndex<0){
+            std::cerr<<"block EndIndex error!EndIndex:"<<aEndIndex<<std::endl;
+            return false;
+        }
+        int dimLength = std::abs(aEndIndex-aBeginIndex)+1;
+        int dataSize = getDataSize()/getDimSize(aDim)*dimLength;
+        double * dataOutput = malloc(dataSize);
+        int colStride = getDimSize(0);
+        int sliceStride = getDimSize(0)*getDimSize(1);
+        switch (aDim) {
+            case 0:{
+                int colStride2 = dimLength;
+                int sliceStride2 = dimLength*getDimSize(1);
+                for (size_t i = 0; i < getDimSize(2); i++)
+                {
+                    for (size_t j = 0; j < getDimSize(1); j++)
+                    {
+                        cblas_dcopy(
+                            dimLength,
+                            &value,0,
+                            getData() + aBeginIndex + j * colStride +
+                                i * sliceStride,
+                            1);
+                    }
+                }
+                return true;
+            }
+            case 1:{
+                int colStride2 = getDimSize(0);
+                int sliceStride2 = dimLength*getDimSize(0);
+                int copySize = sliceStride2;
+                for (size_t i = 0; i < getDimSize(2); i++)
+                {
+                    cblas_dcopy(copySize,
+                                &value, 0 ,
+                                getData() + aBeginIndex * colStride +
+                                    i * sliceStride, 1);
+                }
+                return true;
+            }
+            case 2:{
+                int copySize = dimLength*sliceStride;
+                cblas_dcopy(copySize, &value, 0 ,getData() + aBeginIndex * sliceStride ,1);
+                return true;
+            }
+        }
+
+        
     }
 
     void Matrix::printf() {
diff --git a/src/Matrix.h b/src/Matrix.h
index 38f38ed..964f390 100644
--- a/src/Matrix.h
+++ b/src/Matrix.h
@@ -183,6 +183,8 @@ namespace Aurora {
          */
         Matrix block(int aDim,int aBeginIndx, int aEndIndex) const;
 
+        bool setBlockValue(int aDim,int aBeginIndx, int aEndIndex,double value);
+
         /**
          * 矩阵乘法
          * @attention 目前只支持矩阵乘向量
diff --git a/test/Matrix_Test.cpp b/test/Matrix_Test.cpp
index 6875a30..f8a6e2f 100644
--- a/test/Matrix_Test.cpp
+++ b/test/Matrix_Test.cpp
@@ -5,11 +5,9 @@
 #include "Function1D.h"
 #include "Function2D.h"
 #include "Function3D.h"
+#include "TestUtility.h"
 
 
-
-#define DISPLAY_MATRIX
-
 class Matrix_Test : public ::testing::Test {
 protected:
     static void SetUpMatrixTester() {
@@ -401,4 +399,32 @@ TEST_F(Matrix_Test, matrixfunction){
     EXPECT_DOUBLE_EQ(C.getData()[6], 24);
     EXPECT_DOUBLE_EQ(C.getData()[7], 306);
     EXPECT_DOUBLE_EQ(C.getData()[8], 460);
+
+    double *dataE = Aurora::random(60);
+    Aurora::Matrix E = Aurora::Matrix::New(dataE, 3, 4, 5);
+    Aurora::Matrix block1 = E.block(0, 1, 2);
+    for (size_t i = 0; i < block1.getDimSize(2); i++)
+    {   
+        for (size_t j = 0; j < block1.getDimSize(1); j++)
+        {
+            for (size_t k = 1; k < E.getDimSize(0); k++)
+            {
+                auto index1 = k+j*E.getDimSize(0)+i*E.getDimSize(1)*E.getDimSize(0);
+                auto index2 = k-1+j*block1.getDimSize(0)+i*block1.getDimSize(1)*block1.getDimSize(0);
+                EXPECT_DOUBLE_AE(E[index1],block1[index2]);
+            }
+        }
+    }
+    E.setBlockValue(0, 1, 2,-1);
+    for (size_t i = 0; i < block1.getDimSize(2); i++)
+    {   
+        for (size_t j = 0; j < block1.getDimSize(1); j++)
+        {
+            for (size_t k = 1; k < E.getDimSize(0); k++)
+            {
+                auto index1 = k+j*E.getDimSize(0)+i*E.getDimSize(1)*E.getDimSize(0);
+                EXPECT_DOUBLE_AE(E[index1],-1);
+            }
+        }
+    }
 }
\ No newline at end of file