FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_MKL_wrapper.F90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
6#ifdef WITH_MKL
7include 'mkl_pardiso.f90'
8#endif
9
11 use hecmw_util
13
14#ifdef WITH_MKL
15 use mkl_pardiso
16#endif
17
18 implicit none
19
20 private ! default
21 public :: hecmw_mkl_wrapper ! only entry point of Parallel Direct Solver is public
22
23 logical, save :: INITIALIZED = .false.
24#ifdef WITH_MKL
25 type(MKL_PARDISO_HANDLE) :: pt(64)
26#endif
27 integer maxfct, mnum, mtype, nrhs, msglvl
28 integer idum(1), i
29 data nrhs /1/, maxfct /1/, mnum /1/
30
31 integer, save :: iparm(64)
32
33contains
34
35 subroutine hecmw_mkl_wrapper(spMAT, phase_start, solx, istat)
36 implicit none
37 type (sparse_matrix), intent(inout) :: spmat
38 integer(kind=kint), intent(in) :: phase_start
39 integer(kind=kint), intent(out) :: istat
40 real(kind=kreal), pointer, intent(inout) :: solx(:)
41
42 integer(kind=kint) :: myrank, phase
43 real(kind=kreal) :: t2,t3,t4,t5
44
45#ifdef WITH_MKL
46
47 myrank=hecmw_comm_get_rank()
48
49 if (.not. initialized) then
50 do i=1,64
51 pt(i)%dummy = 0
52 enddo
53 iparm(:) = 0
54 iparm(1) = 1 ! no solver default
55 iparm(2) = 3 ! fill-in reordering from METIS
56 iparm(3) = 1
57 iparm(10) = 13 ! perturbe the pivot elements with 1E-13
58 iparm(11) = 1 ! use nonsymmetric permutation and scaling MPS
59 iparm(13) = 1 ! maximum weighted matching algorithm is switched-off
60 iparm(18) = -1
61 iparm(19) = -1
62 msglvl = 0 ! print statistical information
63 initialized = .true.
64 endif
65
66 if ( phase_start == 1 ) then
67 phase = -1
68 call pardiso(pt, maxfct, mnum, mtype, phase, spmat%N, spmat%A, spmat%IRN, spmat%JCN, idum, &
69 nrhs, iparm, msglvl, spmat%RHS, solx, istat)
70 if (istat < 0) then
71 write(*,*) 'ERROR: MKL returned with error', istat
72 return
73 endif
74
75 if (spmat%symtype == sparse_matrix_symtype_spd) then
76 mtype = 2
77 else if (spmat%symtype == sparse_matrix_symtype_sym) then
78 mtype = -2
79 else
80 mtype = 11
81 endif
82
83 ! ANALYSIS
84 t2=hecmw_wtime()
85
86 phase = 11
87 call pardiso(pt, maxfct, mnum, mtype, phase, spmat%N, spmat%A, spmat%IRN, spmat%JCN, idum, &
88 nrhs, iparm, msglvl, spmat%RHS, solx, istat)
89 if (istat < 0) then
90 write(*,*) 'ERROR: MKL returned with error', istat
91 return
92 endif
93
94 t3=hecmw_wtime()
95 if (myrank==0 .and. (spmat%iterlog > 0 .or. spmat%timelog > 0)) &
96 write(*,'(A,f10.3)') ' [Pardiso]: Analysis completed. time(sec)=',t3-t2
97 endif
98
99 ! FACTORIZATION
100 if ( phase_start .le. 2 ) then
101 t3=hecmw_wtime()
102
103 phase = 22
104 call pardiso(pt, maxfct, mnum, mtype, phase, spmat%N, spmat%A, spmat%IRN, spmat%JCN, idum, &
105 nrhs, iparm, msglvl, spmat%RHS, solx, istat)
106 if (istat < 0) then
107 write(*,*) 'ERROR: MKL returned with error', istat
108 return
109 endif
110
111 t4=hecmw_wtime()
112 if (myrank==0 .and. (spmat%iterlog > 0 .or. spmat%timelog > 0)) &
113 write(*,'(A,f10.3)') ' [Pardiso]: Factorization completed. time(sec)=',t4-t3
114 endif
115
116 ! SOLUTION
117 t4=hecmw_wtime()
118 phase = 33
119 call pardiso(pt, maxfct, mnum, mtype, phase, spmat%N, spmat%A, spmat%IRN, spmat%JCN, idum, &
120 nrhs, iparm, msglvl, spmat%RHS, solx, istat)
121 if (istat < 0) then
122 write(*,*) 'ERROR: MKL returned with error', istat
123 return
124 endif
125
126 t5=hecmw_wtime()
127 if (myrank==0 .and. (spmat%iterlog > 0 .or. spmat%timelog > 0)) &
128 write(*,'(A,f10.3)') ' [Pardiso]: Solution completed. time(sec)=',t5-t4
129
130#else
131 stop "MKL Pardiso not available"
132#endif
133 end subroutine hecmw_mkl_wrapper
134
135end module m_hecmw_mkl_wrapper
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
integer(kind=kint) function hecmw_comm_get_rank()
real(kind=kreal) function hecmw_wtime()
This module provides linear equation solver interface for Pardiso.
subroutine, public hecmw_mkl_wrapper(spmat, phase_start, solx, istat)
This module provides DOF based sparse matrix data structure (CSR and COO)
integer(kind=kint), parameter, public sparse_matrix_symtype_sym
integer(kind=kint), parameter, public sparse_matrix_symtype_spd