feat(elmat): adds option for user allocated mat

labmec · May 28, 2024 · 0bb9c72 · 0bb9c72
1 parent 7688a2b
commit 0bb9c72
Show file tree

Hide file tree

Showing 2 changed files with 13 additions and 7 deletions.
diff --git a/Mesh/TPZElementMatrixT.cpp b/Mesh/TPZElementMatrixT.cpp
@@ -249,6 +249,10 @@ void TPZElementMatrixT<TVar>::ApplyConstraints(){
 
 	int numnodes_processed = 0;
 	int current_order = 0;
+  TPZFNMatrix<150,TVar> localmat;
+
+  TPZFMatrix<TVar> *fulldepmat = fUserAllocMat ? fUserAllocMat : &localmat;
+
 	while(numnodes_processed < totalnodes) {
 		int in;
 		for(in=0; in<totalnodes; in++) {
@@ -268,7 +272,7 @@ void TPZElementMatrixT<TVar>::ApplyConstraints(){
 			// loop over the nodes from which dfn depends
 			TPZConnect::TPZDependBase *dep = dfn->FirstDepend();
 
-      TPZFNMatrix<150,TVar> fulldepmat;
+
 			while(dep) {
         auto dept = dynamic_cast<TPZConnect::TPZDepend<TVar>*>(dep);
         if(!dept){DebugStop();}
@@ -304,15 +308,15 @@ void TPZElementMatrixT<TVar>::ApplyConstraints(){
           const int orig_col=depmat.Cols();
           const int full_row=numstate*orig_row;
           const int full_col=numstate*orig_col;
-          fulldepmat.Redim(full_row,full_col);
+          fulldepmat->Redim(full_row,full_col);
 
           for(int ic = 0; ic < orig_col; ic++){
             const auto first_c = ic*numstate;
             for(int ir = 0; ir < orig_row; ir++){
               const auto first_r = ir*numstate;
               const auto val = depmat.GetVal(ir,ic);
               for(int istate = 0; istate < numstate; istate++){
-                fulldepmat.PutVal(first_r+istate,first_c+istate,val);
+                fulldepmat->PutVal(first_r+istate,first_c+istate,val);
               }
             }
           }
@@ -327,12 +331,12 @@ void TPZElementMatrixT<TVar>::ApplyConstraints(){
           since there are no trial functions
         */
 
-        const auto deprows = fulldepmat.Rows();
-        const auto depcols = fulldepmat.Cols();
+        const auto deprows = fulldepmat->Rows();
+        const auto depcols = fulldepmat->Cols();
         //the window is the full matrix
 
         {
-          TPZMatrixWindow<TVar> dep_window(fulldepmat,0,0,deprows,depcols);
+          TPZMatrixWindow<TVar> dep_window(*fulldepmat,0,0,deprows,depcols);
           TPZMatrixWindow<TVar> send_window(this->fConstrMat,inpos,0,insize,this->fConstrMat.Cols());
           TPZMatrixWindow<TVar> receive_window(this->fConstrMat,deppos,0,depsize,this->fConstrMat.Cols());
           const TVar alpha{1};
@@ -346,7 +350,7 @@ void TPZElementMatrixT<TVar>::ApplyConstraints(){
         }
         if (this->fType == TPZElementMatrix::EK){
           //now we multiply it on the right side too
-          TPZMatrixWindow<TVar> dep_window(fulldepmat,0,0,deprows,depcols);
+          TPZMatrixWindow<TVar> dep_window(*fulldepmat,0,0,deprows,depcols);
           TPZMatrixWindow<TVar> send_window(this->fConstrMat,0,inpos,this->fConstrMat.Rows(),insize);
           TPZMatrixWindow<TVar> receive_window(this->fConstrMat,0,deppos,this->fConstrMat.Rows(),depsize);
           const TVar alpha{1};

diff --git a/Mesh/TPZElementMatrixT.h b/Mesh/TPZElementMatrixT.h
@@ -75,6 +75,7 @@ struct TPZElementMatrixT : public TPZElementMatrix {
     TPZBlock & ConstrBlock() override{
         return fConstrBlock;
     }
+  void SetUserAllocMat(TPZFMatrix<TVar> *mat){fUserAllocMat=mat;}
     /** @brief Pointer to a blocked matrix object*/
 	TPZFNMatrix<1000, TVar> fMat;
     /** @brief Block structure associated with fMat*/
@@ -83,6 +84,7 @@ struct TPZElementMatrixT : public TPZElementMatrix {
 	TPZFNMatrix<1000, TVar> fConstrMat;
     /** @brief Block structure associated with fConstrMat*/
 	TPZBlock fConstrBlock;
+  TPZFMatrix<TVar> *fUserAllocMat{nullptr};
 };
 
 extern template class TPZElementMatrixT<STATE>;