diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index efa9f6a203a5479c0440e2ae89c4a133ecf19b3a..65fc903097c0021e3a6378c1780fee6a407566d3 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -18,6 +18,7 @@ include:
   - .gitlab/test_parsec.yml
   - .gitlab/test_quark.yml
   - .gitlab/test_openmp.yml
+  - .gitlab/test_hmat.yml
   - .gitlab/coverage.yml
   - .gitlab/coverity.yml
   - .gitlab/sonarqube.yml
diff --git a/.gitlab/build.yml b/.gitlab/build.yml
index 5e0fd37742189ae74db6f46926958ab0394e0e68..d311506836a79e08faaa1a9d1108ceb9bab56fc4 100644
--- a/.gitlab/build.yml
+++ b/.gitlab/build.yml
@@ -74,4 +74,15 @@ build_starpu_macosx:
   script:
     - tools/build_macosx.sh
 
+build_hmat:
+  extends: .build_script_template
+  artifacts:
+    name: build_hmat
+    expire_in: 48 hours
+    paths:
+      - build-hmat
+      - chameleon_hmat.log
+  variables:
+    BUILD_OPTIONS: "-DCHAMELEON_USE_CUDA=OFF -DCHAMELEON_USE_MPI=ON -DCHAMELEON_USE_HMAT=ON"
+    VERSION: hmat
 
diff --git a/.gitlab/test_hmat.yml b/.gitlab/test_hmat.yml
new file mode 100644
index 0000000000000000000000000000000000000000..97e2d3702c96ad57116fd37c76ab6e6949686dc2
--- /dev/null
+++ b/.gitlab/test_hmat.yml
@@ -0,0 +1,34 @@
+---
+.test_hmat_template:
+  extends: .test_script_template
+  needs: [build_hmat]
+  dependencies:
+    - build_hmat
+  artifacts:
+    name: test_hmat
+    expire_in: 48 hours
+    paths:
+      - build_hmat
+      - chameleon_hmat_*.log
+      - chameleon_hmat_*.lcov
+      - bis_chameleon_hmat_*.lcov
+    reports:
+      junit: report.xml
+
+.test_hmat_template_master: &test_hmat_master
+  extends:
+    - .test_hmat_template
+    - .only-master
+
+.test_hmat_template_branches: &test_hmat_branches
+  extends:
+    - .test_hmat_template
+    - .only-branches
+
+test_hmat:
+  <<: *test_hmat_branches
+  variables:
+    VERSION: hmat
+    CATEGORY: none
+    PRECISION: all
+    TESTS_RESTRICTION: "-R \"testBEM_\""
diff --git a/.gitmodules b/.gitmodules
index 75281289e3d509740ac3d1bbaacafe260a704dcc..0ab797b12027a56786a5464f6376ccec9350719c 100644
--- a/.gitmodules
+++ b/.gitmodules
@@ -7,3 +7,9 @@
 [submodule "doc/orgmode/org-html-themes"]
 	path = doc/orgmode/org-html-themes
 	url = https://github.com/fniessen/org-html-themes.git
+[submodule "coreblas/hmat-oss"]
+	path = coreblas/hmat-oss
+	url = https://github.com/jeromerobert/hmat-oss.git
+[submodule "testing/test_fembem"]
+	path = testing/test_fembem
+	url = https://gitlab.inria.fr/solverstack/test_fembem.git
diff --git a/CMakeLists.txt b/CMakeLists.txt
index f1f938142e16fb9804510842f0bd7397c255433a..770b404d86abd6731d769cfbfb0b75adb59f6aff 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -199,6 +199,9 @@ endif()
 
 # Enable Hmat-OSS kernels
 option(CHAMELEON_USE_HMAT "Enable hmat kernels" OFF)
+cmake_dependent_option(CHAMELEON_HMAT_EXTERNAL
+                    "Use external hmat-oss library" OFF
+                    "CHAMELEON_USE_HMAT" OFF)
 
 option(CHAMELEON_RUNTIME_SYNC "Enable synchronous task submission when available to debug the code without parallelism" OFF)
 if (CHAMELEON_RUNTIME_SYNC)
@@ -269,17 +272,23 @@ add_subdirectory(hqr)
 # Build dependency HMAT-OSS library #
 #####################################
 if ( CHAMELEON_USE_HMAT )
-  find_package(HMAT REQUIRED)
-  if(HMAT_FOUND)
-    set_target_properties(HMAT::hmat PROPERTIES
-      INTERFACE_INCLUDE_DIRECTORIES "${HMAT_INCLUDE_DIRS}"
-      INTERFACE_COMPILE_DEFINITIONS "${HMAT_DEFINITIONS}"
-      INTERFACE_LINK_DIRECTORIES "${HMAT_LIBRARY_DIRS}"
-      INTERFACE_LINK_LIBRARIES "hmat"
-      )
-    message( STATUS "HMAT: Found" )
+  if ( CHAMELEON_HMAT_EXTERNAL )
+    find_package(HMAT REQUIRED)
+    if(HMAT_FOUND)
+      if( NOT TARGET HMAT::hmat )
+        set_target_properties(HMAT::hmat PROPERTIES
+          INTERFACE_INCLUDE_DIRECTORIES "${HMAT_INCLUDE_DIRS}"
+          INTERFACE_COMPILE_DEFINITIONS "${HMAT_DEFINITIONS}"
+          INTERFACE_LINK_DIRECTORIES "${HMAT_LIBRARY_DIRS}"
+          INTERFACE_LINK_LIBRARIES "hmat"
+          )
+      endif()
+      message( STATUS "HMAT: Found" )
+    else()
+      message( ERROR "HMAT: Not found" )
+    endif()
   else()
-    message( ERROR "HMAT: Not found" )
+    message( STATUS "HMAT: Use internal HMat-OSS" )
   endif()
 endif()
 
diff --git a/cmake_modules/CHAMELEONConfig.cmake.in b/cmake_modules/CHAMELEONConfig.cmake.in
index 89ae753b1d741dd90c70b31a36aaa2daaa0fcfe6..4817a7aa022d75f6459e8c1a93a1430c046468de 100644
--- a/cmake_modules/CHAMELEONConfig.cmake.in
+++ b/cmake_modules/CHAMELEONConfig.cmake.in
@@ -57,6 +57,15 @@ if (@CHAMELEON_USE_CUDA@ AND NOT @CHAMELEON_SIMULATION@)
 endif()
 
 # add the targets file
+if (@CHAMELEON_USE_HMAT@)
+  if ( @CHAMELEON_HMAT_EXTERNAL@ )
+    # Should be a find_dependency() but we fix it to make sure we use
+    # the correct path in higher level projects.
+    include("@HMAT_DIR@/HMATTargets.cmake")
+  else()
+    include("${CMAKE_CURRENT_LIST_DIR}/../hmat/HMATTargets.cmake")
+  endif()
+endif()
 include("${CMAKE_CURRENT_LIST_DIR}/../hqr/hqrTargets.cmake")
 if (NOT @CHAMELEON_SIMULATION@)
   include("${CMAKE_CURRENT_LIST_DIR}/coreblasTargets.cmake")
diff --git a/compute/CMakeLists.txt b/compute/CMakeLists.txt
index 27a2b93fa58ce6787aac0334be50214d955c7027..db91b1d3b5cfed66eafc42fdcc0e16049fbaea4d 100644
--- a/compute/CMakeLists.txt
+++ b/compute/CMakeLists.txt
@@ -289,6 +289,8 @@ set(CHAMELEON_SOURCES_TARGETS "${CHAMELEON_SOURCES_TARGETS};chameleon_sources" C
 # Compile step
 # ------------
 add_library(chameleon ${CHAMELEON_SRCS} ${CHAMELEON_SRCSF})
+add_library(CHAMELEON::chameleon ALIAS chameleon)
+
 set_target_properties(chameleon PROPERTIES VERSION ${CHAMELEON_VERSION})
 set_target_properties(chameleon PROPERTIES SOVERSION ${CHAMELEON_VERSION_MAJOR})
 
diff --git a/coreblas/CMakeLists.txt b/coreblas/CMakeLists.txt
index 91eeb09f45a6b22caf41d1aa3b1d611d7b8817ed..9b660293997dede6e9ee6da546760a6ae3a2db5f 100644
--- a/coreblas/CMakeLists.txt
+++ b/coreblas/CMakeLists.txt
@@ -26,6 +26,9 @@
 #
 ###
 
+if (CHAMELEON_USE_HMAT AND NOT CHAMELEON_HMAT_EXTERNAL)
+  add_subdirectory(hmat-oss)
+endif()
 add_subdirectory(include)
 add_subdirectory(compute)
 add_subdirectory(eztrace_module)
diff --git a/coreblas/compute/hmat_z.c b/coreblas/compute/hmat_z.c
index b232bf71bd580a3be395cae60ce3c7ca866cb21b..aa7f0d219fd206fa40d6fc2502d8945d12bd2766 100644
--- a/coreblas/compute/hmat_z.c
+++ b/coreblas/compute/hmat_z.c
@@ -1,12 +1,12 @@
 /**
  *
- * @file core_zhmat.c
+ * @file hmat_z.c
  *
  * @copyright 2019-2019 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
  *                      Univ. Bordeaux. All rights reserved.
  * @copyright 2019-2019 Universidad Jaume I. All rights reserved.
  *
- * @brief Chameleon CPU kernel interface from CHAM_tile_t layout to the real one.
+ * @brief Chameleon interface for H-Mat kernels
  *
  * @version 1.0.0
  * @author Rocio Carratala-Saez
@@ -22,16 +22,31 @@
  */
 hmat_interface_t hmat_zinterface;
 
+/**
+ * @brief Constructor of the C++ interface
+ */
 void __hmat_zinit() __attribute__(( constructor ));
 void __hmat_zinit() {
     hmat_init_default_interface( &hmat_zinterface, HMAT_DOUBLE_COMPLEX );
 }
 
+/**
+ * @brief Destructor of the C++ interface
+ */
 void __hmat_zfini() __attribute__(( destructor ));
 void __hmat_zfini() {
     hmat_zinterface.finalize();
 }
 
+/**
+ * @brief Cholesky factorization of an H-Matrix
+ *
+ * @param[in,out] A
+ *      On entry, the input matrix A.
+ *      On exit, the factorized matrix A.
+ *
+ * @return 0 (Return an integer to match prototype of LAPACK)
+ */
 int hmat_zpotrf( hmat_matrix_t *A ) {
     hmat_factorization_context_t ctx_facto;
     hmat_factorization_context_init( &ctx_facto );
@@ -40,6 +55,15 @@ int hmat_zpotrf( hmat_matrix_t *A ) {
     return 0;
 }
 
+/**
+ * @brief LU factorization of an H-Matrix
+ *
+ * @param[in,out] A
+ *      On entry, the input matrix A.
+ *      On exit, the factorized matrix A.
+ *
+ * @return 0 (Return an integer to match prototype of LAPACK)
+ */
 int hmat_zgetrf( hmat_matrix_t *A ) {
     hmat_factorization_context_t ctx_facto;
     hmat_factorization_context_init( &ctx_facto );
@@ -48,21 +72,147 @@ int hmat_zgetrf( hmat_matrix_t *A ) {
     return 0;
 }
 
+/**
+ * @brief Matrix-Matrix product of H-Matrices
+ *
+ * Computes \f[ C = \alpha op(A) op(B) + \beta C \f]
+ *
+ * where op(A) = A, A^t or conj(A^t)
+ *
+ * @param[in] transA
+ *      'N' : op(A) = A
+ *      'T' : op(A) = A^t
+ *      'C' : op(A) = conj(A^t)
+ *
+ * @param[in] transB
+ *      'N' : op(B) = B
+ *      'T' : op(B) = B'
+ *      'C' : op(B) = conj(B')
+ *
+ * @param[in] alpha
+ *      The scalar alpha.
+ *
+ * @param[in] alpha
+ *      The scalar alpha.
+ *
+ * @param[in] A
+ *      On entry, the input matrix A.
+ *
+ * @param[in] B
+ *      On entry, the input matrix B.
+ *
+ * @param[in] beta
+ *      The scalar beta.
+ *
+ * @param[in,out] C
+ *      On entry, the input matrix C.
+ *      On exit, the updated matrix C.
+ *
+ * @return 0 (Return an integer to match prototype of LAPACK)
+ */
 int hmat_zgemm( char transA, char transB, void* alpha, hmat_matrix_t* A,
 		hmat_matrix_t* B, void* beta, hmat_matrix_t* C ) {
     return hmat_zinterface.gemm( transA, transB, alpha, A, B, beta, C );
 }
 
+/**
+ * @brief Matrix-Matrix product of an H-Matrix with a full-rank matrix.
+ *
+ * Computes \f[ C = \alpha op(A) B + \beta C \f]
+ *
+ * where op(A) = A, A^t or conj(A^t), and A is an H-Matrix, and B and C are
+ * full-rank matrices.
+ *
+ * @param[in] transA
+ *      'N' : op(A) = A
+ *      'T' : op(A) = A^t
+ *      'C' : op(A) = conj(A^t)
+ *
+ * @param[in] alpha
+ *      The scalar alpha.
+ *
+ * @param[in] alpha
+ *      The scalar alpha.
+ *
+ * @param[in] A
+ *      On entry, the input H-Matrix A.
+ *
+ * @param[in] B
+ *      On entry, the full-rank matrix B.
+ *
+ * @param[in] beta
+ *      The scalar beta.
+ *
+ * @param[in,out] C
+ *      On entry, the full-rank matrix C.
+ *      On exit, the updated full-rank matrix C.
+ *
+ * @return 0 on success, -1 otherwise
+ */
 int hmat_zgemv( char transA, void* alpha, hmat_matrix_t* A,
 		void* B, void* beta, void* C, int nrhs ) {
     return hmat_zinterface.gemm_scalar( transA, alpha, A, B, beta, C, nrhs );
 }
 
+/**
+ * @brief Solve a triangular system \f[ A x = b \f] with A being an H-Matrix
+ *
+ * Solves \f[ \alpha op(A) X = B \f] or  \f[ \alpha X op(A) = B \f]
+ *
+ * where op(A) = A, A^t or conj(A^t), and A is an H-Matrix, and B and C are
+ * full-rank matrices.
+ *
+ * @param[in] side
+ *      'L' : Solves \alpha A X = B
+ *      'R' : Solves \alpha X A = B
+ *
+ * @param[in] uplo
+ *      'L' : A is lower triangular, upper part is not referenced
+ *      'U' : A is upper triangular, lower part is not referenced
+ *
+ * @param[in] trans
+ *      'N' : op(A) = A
+ *      'T' : op(A) = A^t
+ *      'C' : op(A) = conj(A^t)
+ *
+ * @param[in] uplo
+ *      'U' : A has a unitary diagonal which is ot referenced
+ *      'N' : A has a non unitary diagonal
+ *
+ * @param[in] m
+ *      The number of rows of B
+ *
+ * @param[in] n
+ *      The number of columns of B
+ *
+ * @param[in] alpha
+ *      The scalar alpha.
+ *
+ * @param[in] A
+ *      On entry, the input H-Matrix A.
+ *
+ * @param[in] is_b_hmat
+ *      if true, B is an H-Matrix, otherwise B is full-rank and stored in lapack layout
+ *
+ * @param[in] B
+ *      On entry, the H-Mat of full-rank matrix B.
+ *      On entry, the solution of the trinagular system.
+ *
+ * @return 0 on success, -1 otherwise
+ */
 int hmat_ztrsm( char side, char uplo, char trans, char diag, int m, int n,
 		void* alpha, hmat_matrix_t* A, int is_b_hmat, void* B ) {
     return hmat_zinterface.trsm( side, uplo, trans, diag, m, n, alpha, A, is_b_hmat, B );
 }
 
+/**
+ * @brief Read/Unpack an H-Matrix from a packed format into a contiguous buffer
+ *
+ * @param[in] buffer
+ *      The buffer that contains the packed H-Matrix.
+ *
+ * @return The unpack H-Matrix
+ */
 hmat_matrix_t *hmat_zread( void *buffer ) {
     hmat_buffer_comm_t buffer_struct = { 0, (char*)buffer };
     hmat_matrix_t *hmat = hmat_zinterface.read_struct( (hmat_iostream)(&buffer_comm_read),
@@ -71,6 +221,14 @@ hmat_matrix_t *hmat_zread( void *buffer ) {
     return hmat;
 }
 
+/**
+ * @brief Compute the size required to store the H-Matrix in a packed format
+ *
+ * @param[in] hmat
+ *      The H-Matrix for which the size must be computed
+ *
+ * @return The size of the H-Matrix (structure+data)
+ */
 size_t hmat_zsize( hmat_matrix_t *hmat ) {
     size_t size = 0;
     hmat_zinterface.write_struct( hmat, (hmat_iostream)(&buffer_comm_size), (void *)(&size) );
@@ -78,12 +236,28 @@ size_t hmat_zsize( hmat_matrix_t *hmat ) {
     return size;
 }
 
+/**
+ * @brief Write/pack an H-Matrix in a packed format in a given buffer
+ *
+ * @param[in] hmat
+ *      The H-Matrix to pack
+ *
+ * @param[in,out] ptr
+ *      On entry the allocated buffer with sufficent space to store the H-Matrix.
+ *      On exit, contains the packed H-Matrix.
+ */
 void hmat_zwrite( hmat_matrix_t *hmat, char *ptr ) {
     hmat_buffer_comm_t buffer_struct = { 0, ptr };
     hmat_zinterface.write_struct( hmat, (hmat_iostream)(&buffer_comm_write), (void *)(&buffer_struct) );
     hmat_zinterface.write_data(   hmat, (hmat_iostream)(&buffer_comm_write), (void *)(&buffer_struct) );
 }
 
+/**
+ * @brief Free the data structure (structure+data) associated to an H-Matrix
+ *
+ * @param[in] hmat
+ *      The H-Matrix to destroy
+ */
 void hmat_zdestroy( hmat_matrix_t *hmat ) {
     hmat_zinterface.destroy( hmat );
 }
diff --git a/coreblas/hmat-oss b/coreblas/hmat-oss
new file mode 160000
index 0000000000000000000000000000000000000000..e0f34af5eaf36329a336742183e17755df174b4c
--- /dev/null
+++ b/coreblas/hmat-oss
@@ -0,0 +1 @@
+Subproject commit e0f34af5eaf36329a336742183e17755df174b4c
diff --git a/testing/CMakeLists.txt b/testing/CMakeLists.txt
index de990d69c2034a4968e3f1e8594d2278dbeaf438..ee99eb5983a823174d044869d947d45e4b0180b7 100644
--- a/testing/CMakeLists.txt
+++ b/testing/CMakeLists.txt
@@ -169,6 +169,23 @@ include(CTestLists.cmake)
 # copy input files
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/input DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
 
+if ( CHAMELEON_USE_HMAT AND (NOT CHAMELEON_SIMULATION) AND CHAMELEON_USE_MPI )
+
+  add_subdirectory( test_fembem )
+
+  add_test( testBEM_cham_gemv_c test_fembem/test_FEMBEM --bem -nbpts 1500 -c -gemvchameleon  )
+  add_test( testBEM_cham_gemv_d test_fembem/test_FEMBEM --bem -nbpts 1500 -d -gemvchameleon  )
+  add_test( testBEM_cham_posv_c test_fembem/test_FEMBEM --bem -nbpts  800 -c -solvechameleon )
+  add_test( testBEM_cham_gesv_s test_fembem/test_FEMBEM --bem -nbpts 1000 -s -solvechameleon --nosym)
+  add_test( testBEM_cham_gesv_z test_fembem/test_FEMBEM --bem -nbpts  900 -z -solvechameleon --nosym)
+
+  add_test( testBEM_hcham_gesv_s test_fembem/test_FEMBEM  --bem -nbpts 4000 -s -solvehchameleon --nosym -check_result)
+  add_test( testBEM_hcham_gesv_d test_fembem/test_FEMBEM  --bem -nbpts 3000 -d -solvehchameleon --nosym -check_result)
+  add_test( testBEM_hcham_gesv_c test_fembem/test_FEMBEM  --bem -nbpts 2000 -c -solvehchameleon --nosym -check_result)
+  add_test( testBEM_hcham_gesv_z test_fembem/test_FEMBEM  --bem -nbpts 1500 -z -solvehchameleon --nosym -check_result)
+
+endif()
+
 ###
 ### END CMakeLists.txt
 ###
diff --git a/testing/test_fembem b/testing/test_fembem
new file mode 160000
index 0000000000000000000000000000000000000000..59d1983493a95ec483b4d931afcc70b687474aa5
--- /dev/null
+++ b/testing/test_fembem
@@ -0,0 +1 @@
+Subproject commit 59d1983493a95ec483b4d931afcc70b687474aa5