1a4963045SJacob Faibussowitsch #pragma once 2586621ddSJed Brown 3586621ddSJed Brown #include <petscsys.h> 4586621ddSJed Brown 5586621ddSJed Brown #if defined(PETSC_USE_COMPLEX) 6586621ddSJed Brown #define CHOLMOD_SCALAR_TYPE CHOLMOD_COMPLEX 7586621ddSJed Brown #else 8586621ddSJed Brown #define CHOLMOD_SCALAR_TYPE CHOLMOD_REAL 9586621ddSJed Brown #endif 10586621ddSJed Brown 11586621ddSJed Brown #if defined(PETSC_USE_64BIT_INDICES) 12586621ddSJed Brown #define CHOLMOD_INT_TYPE CHOLMOD_LONG 13586621ddSJed Brown #define cholmod_X_start cholmod_l_start 14586621ddSJed Brown #define cholmod_X_analyze cholmod_l_analyze 15d755ee67SBarry Smith /* the type casts are needed because PetscInt is long long while SuiteSparse_long is long and compilers warn even when they are identical */ 16d755ee67SBarry Smith #define cholmod_X_analyze_p(a, b, c, d, f) cholmod_l_analyze_p(a, (SuiteSparse_long *)b, (SuiteSparse_long *)c, d, f) 17586621ddSJed Brown #define cholmod_X_copy cholmod_l_copy 18586621ddSJed Brown #define cholmod_X_factorize cholmod_l_factorize 19586621ddSJed Brown #define cholmod_X_finish cholmod_l_finish 20586621ddSJed Brown #define cholmod_X_free_factor cholmod_l_free_factor 21586621ddSJed Brown #define cholmod_X_free_dense cholmod_l_free_dense 22d755ee67SBarry Smith #define cholmod_X_resymbol(a, b, c, d, f, e) cholmod_l_resymbol(a, (SuiteSparse_long *)b, c, d, f, e) 23586621ddSJed Brown #define cholmod_X_solve cholmod_l_solve 24218db3c1SStefano Zampini #define cholmod_X_solve2 cholmod_l_solve2 25a2fc1e05SToby Isaac #define cholmod_X_check_sparse cholmod_l_check_sparse 26586621ddSJed Brown #else 27586621ddSJed Brown #define CHOLMOD_INT_TYPE CHOLMOD_INT 28586621ddSJed Brown #define cholmod_X_start cholmod_start 29586621ddSJed Brown #define cholmod_X_analyze cholmod_analyze 30586621ddSJed Brown #define cholmod_X_analyze_p cholmod_analyze_p 31586621ddSJed Brown #define cholmod_X_copy cholmod_copy 32586621ddSJed Brown #define cholmod_X_factorize cholmod_factorize 33586621ddSJed Brown #define cholmod_X_finish cholmod_finish 34586621ddSJed Brown #define cholmod_X_free_factor cholmod_free_factor 35586621ddSJed Brown #define cholmod_X_free_dense cholmod_free_dense 36586621ddSJed Brown #define cholmod_X_resymbol cholmod_resymbol 37586621ddSJed Brown #define cholmod_X_solve cholmod_solve 38218db3c1SStefano Zampini #define cholmod_X_solve2 cholmod_solve2 39a2fc1e05SToby Isaac #define cholmod_X_check_sparse cholmod_check_sparse 40586621ddSJed Brown #endif 41586621ddSJed Brown 42586621ddSJed Brown EXTERN_C_BEGIN 43586621ddSJed Brown #include <cholmod.h> 44a2fc1e05SToby Isaac #include <SuiteSparseQR_C.h> 45586621ddSJed Brown EXTERN_C_END 46586621ddSJed Brown 47586621ddSJed Brown typedef struct { 48218db3c1SStefano Zampini PetscErrorCode (*Wrap)(Mat, PetscBool, cholmod_sparse *, PetscBool *, PetscBool *); 49586621ddSJed Brown cholmod_sparse *matrix; 50586621ddSJed Brown cholmod_factor *factor; 51586621ddSJed Brown cholmod_common *common; 52a2fc1e05SToby Isaac SuiteSparseQR_C_factorization *spqrfact; 53*e2b46ddfSPierre Jolivet PetscScalar scale; 54ace3abfcSBarry Smith PetscBool pack; 5502ef7dfaSPierre Jolivet PetscBool normal; 56586621ddSJed Brown } Mat_CHOLMOD; 57586621ddSJed Brown 585a576424SJed Brown PETSC_INTERN PetscErrorCode CholmodStart(Mat); 595a576424SJed Brown PETSC_INTERN PetscErrorCode MatView_CHOLMOD(Mat, PetscViewer); 605a576424SJed Brown PETSC_INTERN PetscErrorCode MatCholeskyFactorSymbolic_CHOLMOD(Mat, Mat, IS, const MatFactorInfo *); 61218db3c1SStefano Zampini PETSC_INTERN PetscErrorCode MatGetInfo_CHOLMOD(Mat, MatInfoType, MatInfo *); 625a576424SJed Brown PETSC_INTERN PetscErrorCode MatDestroy_CHOLMOD(Mat); 63586621ddSJed Brown 64a2fc1e05SToby Isaac PETSC_INTERN PetscErrorCode VecWrapCholmod(Vec, PetscInt, cholmod_dense *); 65a2fc1e05SToby Isaac PETSC_INTERN PetscErrorCode VecUnWrapCholmod(Vec, PetscInt, cholmod_dense *); 66a2fc1e05SToby Isaac PETSC_INTERN PetscErrorCode MatDenseWrapCholmod(Mat, PetscInt, cholmod_dense *); 67a2fc1e05SToby Isaac PETSC_INTERN PetscErrorCode MatDenseUnWrapCholmod(Mat, PetscInt, cholmod_dense *); 68