Actual source code: htool.hpp
1: #include <petsc/private/matimpl.h>
2: #include <htool/misc/petsc.hpp>
4: class WrapperHtool : public htool::VirtualGenerator<PetscScalar> {
5: PetscInt dim;
6: MatHtoolKernel &kernel;
7: void *ctx;
9: public:
10: WrapperHtool(PetscInt M, PetscInt N, PetscInt sdim, MatHtoolKernel &g, void *kernelctx) : VirtualGenerator(M, N), dim(sdim), kernel(g), ctx(kernelctx) { }
11: void copy_submatrix(PetscInt M, PetscInt N, const PetscInt *rows, const PetscInt *cols, PetscScalar *ptr) const
12: {
13: #if !PetscDefined(HAVE_OPENMP)
14: #endif
15: PETSC_COMM_SELF, kernel(dim, M, N, rows, cols, ptr, ctx);
16: #if !PetscDefined(HAVE_OPENMP)
17: return;
18: #endif
19: }
20: };
22: struct Mat_Htool {
23: PetscInt dim;
24: PetscReal *gcoords_target;
25: PetscReal *gcoords_source;
26: PetscScalar *work_target;
27: PetscScalar *work_source;
28: PetscScalar s;
29: PetscInt bs[2];
30: PetscReal epsilon;
31: PetscReal eta;
32: PetscInt depth[2];
33: MatHtoolCompressorType compressor;
34: MatHtoolClusteringType clustering;
35: MatHtoolKernel kernel;
36: void *kernelctx;
37: WrapperHtool *wrapper;
38: htool::VirtualHMatrix<PetscScalar> *hmatrix;
39: };
41: struct MatHtoolKernelTranspose {
42: Mat A;
43: MatHtoolKernel kernel;
44: void *kernelctx;
45: };