FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
solve_LINEQ_mkl.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!-------------------------------------------------------------------------------
7
9 use hecmw_util
10 use m_fstr
16
17 implicit none
18
19 private
22
23 logical, save :: NEED_ANALYSIS = .true.
24 type (sparse_matrix), save :: spMAT
25
26contains
27
28 subroutine solve_lineq_mkl_contact_init(hecMESH,is_sym)
29 type (hecmwst_local_mesh), intent(in) :: hecmesh
30 logical, intent(in) :: is_sym
31
32 integer(kind=kint) :: spmat_type
33 integer(kind=kint) :: spmat_symtype
34 integer(kind=kint) :: i
35
36 call sparse_matrix_finalize(spmat)
37
38 if (is_sym) then
39 spmat_symtype = sparse_matrix_symtype_sym
40 else
41 spmat_symtype = sparse_matrix_symtype_asym
42 end if
43 if(hecmesh%PETOT.GT.1) then
44 spmat_type = sparse_matrix_type_coo
45 else
46 spmat_type = sparse_matrix_type_csr
47 endif
48 call sparse_matrix_set_type(spmat, spmat_type, spmat_symtype)
49
50 need_analysis = .true.
51 end subroutine solve_lineq_mkl_contact_init
52
54 subroutine solve_lineq_mkl_contact(hecMESH,hecMAT,fstrMAT,istat,conMAT)
55 type (hecmwst_local_mesh), intent(in) :: hecmesh
56 type (hecmwst_matrix ), intent(inout) :: hecmat
57 type (fstrst_matrix_contact_lagrange), intent(inout) :: fstrmat
58 integer(kind=kint), intent(out) :: istat
59 type (hecmwst_matrix), intent(in),optional :: conmat
60
61 integer(kind=kint) :: phase_start
62 real(kind=kreal) :: t1,t2
63
64 t1=hecmw_wtime()
65 call hecmw_mat_dump(hecmat, hecmesh)
66
67 if (need_analysis) then
68 !constrtuct new structure
69 call sparse_matrix_contact_init_prof(spmat, hecmat, fstrmat, hecmesh)
70 endif
71
72 ! ---- For Parallel Contact with Multi-Partition Domains
73 if(paracontactflag.and.present(conmat)) then
74 call sparse_matrix_para_contact_set_vals(spmat, hecmat, fstrmat, conmat)
75 call sparse_matrix_para_contact_set_rhs(spmat, hecmat, fstrmat, conmat)
76 else
77 call sparse_matrix_contact_set_vals(spmat, hecmat, fstrmat)
78 !call sparse_matrix_dump(spMAT)
79 call sparse_matrix_contact_set_rhs(spmat, hecmat, fstrmat)
80 endif
81
82 t2=hecmw_wtime()
83 if (myrank==0 .and. spmat%timelog > 0) then
84 if( hecmesh%PETOT .GT. 1 ) then
85 write(*,'(A,f10.3)') ' [Cluster Pardiso]: Setup completed. time(sec)=',t2-t1
86 else
87 write(*,'(A,f10.3)') ' [Pardiso]: Setup completed. time(sec)=',t2-t1
88 end if
89 endif
90
91 phase_start = 2
92 if (need_analysis) then
93 phase_start = 1
94 need_analysis = .false.
95 endif
96
97 ! SOLVE
98 if( hecmesh%PETOT.GT.1 ) then
99 call hecmw_clustermkl_wrapper(spmat, phase_start, hecmat%X, istat)
100 call sparse_matrix_contact_get_rhs(spmat, hecmat, fstrmat)
101 else
102 call hecmw_mkl_wrapper(spmat, phase_start, hecmat%X, istat)
103 deallocate(spmat%rhs)
104 endif
105
106 call hecmw_mat_dump_solution(hecmat)
107 end subroutine solve_lineq_mkl_contact
108
This module provides functions of reconstructing.
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint) myrank
PARALLEL EXECUTION.
Definition: m_fstr.f90:80
logical paracontactflag
PARALLEL CONTACT FLAG.
Definition: m_fstr.f90:84
This module provides linear equation solver interface for Cluster Pardiso.
subroutine, public hecmw_clustermkl_wrapper(spmat, phase_start, solution, istat)
This module provides linear equation solver interface for Pardiso.
subroutine, public hecmw_mkl_wrapper(spmat, phase_start, solx, istat)
This module provides functions to solve sparse system of \linear equitions using intel MKL direct spa...
subroutine, public solve_lineq_mkl_contact_init(hecmesh, is_sym)
subroutine, public solve_lineq_mkl_contact(hecmesh, hecmat, fstrmat, istat, conmat)
This subroutine executes the MKL solver.
This module provides conversion routines between HEC and FISTR data structures and DOF based sparse m...
subroutine, public sparse_matrix_contact_init_prof(spmat, hecmat, fstrmat, hecmesh)
subroutine, public sparse_matrix_para_contact_set_rhs(spmat, hecmat, fstrmat, conmat)
subroutine, public sparse_matrix_contact_set_vals(spmat, hecmat, fstrmat)
subroutine, public sparse_matrix_contact_get_rhs(spmat, hecmat, fstrmat)
subroutine, public sparse_matrix_contact_set_rhs(spmat, hecmat, fstrmat)
subroutine, public sparse_matrix_para_contact_set_vals(spmat, hecmat, fstrmat, conmat)
This module provides DOF based sparse matrix data structure (CSR and COO)
integer(kind=kint), parameter, public sparse_matrix_type_coo
integer(kind=kint), parameter, public sparse_matrix_type_csr
subroutine, public sparse_matrix_set_type(spmat, type, symtype)
integer(kind=kint), parameter, public sparse_matrix_symtype_sym
subroutine, public sparse_matrix_finalize(spmat)
integer(kind=kint), parameter, public sparse_matrix_symtype_asym
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...