1!------------------------------------------------------------------------------- 2! Copyright (c) 2019 FrontISTR Commons 3! This software is released under the MIT License, see LICENSE.txt 4!------------------------------------------------------------------------------- 5!> This module provides linear equation solver interface for Pardiso 6 7module hecmw_solver_direct_MKL 8 use hecmw_util 9 use m_sparse_matrix 10 use m_sparse_matrix_hec 11 use m_hecmw_MKL_wrapper 12 use hecmw_matrix_ass 13 use hecmw_matrix_dump 14 15 private 16 public :: hecmw_solve_direct_MKL 17 18 logical, save :: INITIALIZED = .false. 19 type (sparse_matrix), save :: spMAT 20 21contains 22 23 subroutine hecmw_solve_direct_MKL(hecMESH,hecMAT) 24 implicit none 25 type (hecmwST_local_mesh), intent(in) :: hecMESH 26 type (hecmwST_matrix ), intent(inout) :: hecMAT 27 integer(kind=kint) :: spmat_type 28 integer(kind=kint) :: spmat_symtype 29 integer(kind=kint) :: phase_start 30 integer(kind=kint) :: istat,myrank 31 real(kind=kreal) :: t1,t2 32 33 t1=hecmw_wtime() 34 35 call hecmw_mat_dump(hecMAT, hecMESH) 36 37 myrank=hecmw_comm_get_rank() 38 39 if (INITIALIZED .and. hecMAT%Iarray(98) .eq. 1) then 40 call sparse_matrix_finalize(spMAT) 41 INITIALIZED = .false. 42 endif 43 44 if (.not. INITIALIZED) then 45 spmat_type = SPARSE_MATRIX_TYPE_CSR 46 if (hecMAT%symmetric) then 47 spmat_symtype = SPARSE_MATRIX_SYMTYPE_SYM 48 else 49 spmat_symtype = SPARSE_MATRIX_SYMTYPE_ASYM 50 end if 51 call sparse_matrix_set_type(spMAT, spmat_type, spmat_symtype) 52 53 INITIALIZED = .true. 54 hecMAT%Iarray(98) = 1 55 endif 56 57 !* Flag to activate symbolic factorization: 1(yes) 0(no) hecMESH%Iarray(98) 58 !* Flag to activate numeric factorization: 1(yes) 0(no) hecMESH%Iarray(97) 59 60 if (hecMAT%Iarray(98) .eq. 1) then 61 ! ANALYSIS and FACTORIZATION 62 call sparse_matrix_hec_init_prof(spMAT, hecMAT, hecMESH) 63 call sparse_matrix_hec_set_vals(spMAT, hecMAT) 64 !call sparse_matrix_dump(spMAT) 65 phase_start = 1 66 hecMAT%Iarray(98) = 0 67 hecMAT%Iarray(97) = 0 68 endif 69 if (hecMAT%Iarray(97) .eq. 1) then 70 ! FACTORIZATION 71 call sparse_matrix_hec_set_vals(spMAT, hecMAT) 72 !call sparse_matrix_dump(spMAT) 73 phase_start = 2 74 hecMAT%Iarray(97) = 0 75 endif 76 call sparse_matrix_hec_set_rhs(spMAT, hecMAT) 77 78 t2=hecmw_wtime() 79 if (myrank==0 .and. (spMAT%iterlog > 0 .or. spMAT%timelog > 0)) & 80 write(*,'(A,f10.3)') ' [Pardiso]: Setup completed. time(sec)=',t2-t1 81 82 ! SOLVE 83 call hecmw_mkl_wrapper(spMAT, phase_start, hecMAT%X, istat) 84 85 call hecmw_mat_dump_solution(hecMAT) 86 87 end subroutine hecmw_solve_direct_MKL 88 89end module hecmw_solver_direct_MKL 90