1 #pragma once 2 3 #include <petsc/private/matimpl.h> 4 #include <petscmathtool.h> 5 6 PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wsign-compare") 7 #include <htool/clustering/tree_builder/tree_builder.hpp> 8 #include <htool/hmatrix/lrmat/SVD.hpp> 9 #include <htool/hmatrix/lrmat/fullACA.hpp> 10 #include <htool/hmatrix/lrmat/recompressed_low_rank_generator.hpp> 11 #include <htool/hmatrix/hmatrix_distributed_output.hpp> 12 #include <htool/hmatrix/linalg/factorization.hpp> 13 #include <htool/distributed_operator/utility.hpp> 14 #include <htool/distributed_operator/linalg/add_distributed_operator_vector_product_local_to_local.hpp> 15 #include <htool/distributed_operator/linalg/add_distributed_operator_matrix_product_local_to_local.hpp> 16 PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END() 17 18 class WrapperHtool : public htool::VirtualGenerator<PetscScalar> { 19 MatHtoolKernelFn *&kernel; 20 PetscInt sdim; 21 void *ctx; 22 23 public: WrapperHtool(PetscInt dim,MatHtoolKernelFn * & g,void * kernelctx)24 WrapperHtool(PetscInt dim, MatHtoolKernelFn *&g, void *kernelctx) : VirtualGenerator<PetscScalar>(), kernel(g), sdim(dim), ctx(kernelctx) { } copy_submatrix(PetscInt M,PetscInt N,const PetscInt * rows,const PetscInt * cols,PetscScalar * ptr) const25 void copy_submatrix(PetscInt M, PetscInt N, const PetscInt *rows, const PetscInt *cols, PetscScalar *ptr) const 26 { 27 #if !PetscDefined(HAVE_OPENMP) 28 PetscFunctionBegin; 29 #endif 30 PetscCallAbort(PETSC_COMM_SELF, kernel(sdim, M, N, rows, cols, ptr, ctx)); 31 #if !PetscDefined(HAVE_OPENMP) 32 PetscFunctionReturnVoid(); 33 #endif 34 } 35 }; 36 37 struct Mat_Htool { 38 PetscInt dim; 39 PetscReal *gcoords_target; 40 PetscReal *gcoords_source; 41 PetscInt max_cluster_leaf_size; 42 PetscReal epsilon; 43 PetscReal eta; 44 PetscInt depth[2]; 45 PetscBool block_tree_consistency; 46 PetscBool permutation; 47 PetscBool recompression; 48 MatHtoolCompressorType compressor; 49 MatHtoolClusteringType clustering; 50 MatHtoolKernelFn *kernel; 51 void *kernelctx; 52 WrapperHtool *wrapper; 53 std::unique_ptr<htool::Cluster<PetscReal>> target_cluster; 54 std::unique_ptr<htool::Cluster<PetscReal>> source_cluster; 55 std::unique_ptr<htool::DefaultApproximationBuilder<PetscScalar>> distributed_operator_holder; 56 }; 57 58 struct MatHtoolKernelTranspose { 59 Mat A; 60 MatHtoolKernelFn *kernel; 61 void *kernelctx; 62 }; 63