FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
solve_LINEQ_MUMPS_contact.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 !-------------------------------------------------------------------------------
8  use m_fstr
13 
14  private
17 
18  logical, save :: INITIALIZED = .false.
19  logical, save :: NEED_ANALYSIS = .true.
20  type (sparse_matrix), save :: spMAT
21 
22 contains
23 
24  subroutine solve_lineq_mumps_contact_init(hecMESH,hecMAT,fstrMAT,is_sym)
25  implicit none
26  type (hecmwst_local_mesh), intent(in) :: hecmesh
27  type (hecmwst_matrix ), intent(inout) :: hecmat
28  type (fstrst_matrix_contact_lagrange), intent(in) :: fstrmat
29  logical, intent(in) :: is_sym
30  integer(kind=kint) :: spmat_type
31  integer(kind=kint) :: spmat_symtype
32  integer(kind=kint) :: mumps_job
33  integer(kind=kint) :: istat
34 
35  if (initialized) then
36  mumps_job=-2
37  call hecmw_mumps_wrapper(spmat, mumps_job, istat)
38  if (istat < 0) then
39  write(*,*) 'ERROR: MUMPS returned with error', istat
40  stop
41  endif
42  call sparse_matrix_finalize(spmat)
43  initialized = .false.
44  endif
45 
46  if (is_sym) then
47  !spmat_symtype = SPARSE_MATRIX_SYMTYPE_SPD
48  spmat_symtype = sparse_matrix_symtype_sym
49  else
50  spmat_symtype = sparse_matrix_symtype_asym
51  endif
52  spmat_type = sparse_matrix_type_coo
53  call sparse_matrix_set_type(spmat, spmat_type, spmat_symtype)
54  mumps_job=-1
55  call hecmw_mumps_wrapper(spmat, mumps_job, istat)
56  if (istat < 0) then
57  write(*,*) 'ERROR: MUMPS returned with error', istat
58  stop
59  endif
60 
61  need_analysis = .true.
62  initialized = .true.
63  end subroutine solve_lineq_mumps_contact_init
64 
65  subroutine solve_lineq_mumps_contact(hecMESH,hecMAT,fstrMAT,istat,conMAT)
66  implicit none
67  type (hecmwst_local_mesh), intent(in) :: hecmesh
68  type (hecmwst_matrix ), intent(inout) :: hecmat
69  type (fstrst_matrix_contact_lagrange), intent(inout) :: fstrmat
70  integer(kind=kint), intent(out) :: istat
71  type (hecmwst_matrix), intent(in),optional :: conmat
72 
73  integer(kind=kint) :: mumps_job
74 
75  call hecmw_mat_dump(hecmat, hecmesh)
76 
77  ! ANALYSIS
78  if (need_analysis) then
79  call sparse_matrix_contact_init_prof(spmat, hecmat, fstrmat, hecmesh)
80  mumps_job=1
81  call hecmw_mumps_wrapper(spmat, mumps_job, istat)
82  if (istat < 0) then
83  write(*,*) 'ERROR: MUMPS returned with error', istat
84  stop
85  endif
86  if (myrank==0) write(*,*) ' [MUMPS]: Analysis completed.'
87  need_analysis = .false.
88  endif
89 
90  ! FACTORIZATION and SOLUTION
91  ! ---- For Parallel Contact with Multi-Partition Domains
92  if(paracontactflag.and.present(conmat)) then
93  call sparse_matrix_para_contact_set_vals(spmat, hecmat, fstrmat, conmat)
94  call sparse_matrix_para_contact_set_rhs(spmat, hecmat, fstrmat, conmat)
95  else
96  call sparse_matrix_contact_set_vals(spmat, hecmat, fstrmat)
97  !call sparse_matrix_dump(spMAT)
98  call sparse_matrix_contact_set_rhs(spmat, hecmat, fstrmat)
99  endif
100  mumps_job=5
101  call hecmw_mumps_wrapper(spmat, mumps_job, istat)
102  if (istat < 0) then
103  write(*,*) 'ERROR: MUMPS returned with error', istat
104  return
105  endif
106  call sparse_matrix_contact_get_rhs(spmat, hecmat, fstrmat)
107  if (myrank==0) write(*,*) ' [MUMPS]: Factorization and Solution completed.'
108 
109  call hecmw_mat_dump_solution(hecmat)
110 
111  !call sparse_matrix_finalize(spMAT)
112  end subroutine solve_lineq_mumps_contact
113 
This module provides functions of reconstructing.
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 wrapper for parallel sparse direct solver MUMPS.
subroutine, public hecmw_mumps_wrapper(spMAT, job, istat)
This module provides linear equation solver interface of MUMPS for contact problems using Lagrange mu...
subroutine, public solve_lineq_mumps_contact_init(hecMESH, hecMAT, fstrMAT, is_sym)
subroutine, public solve_lineq_mumps_contact(hecMESH, hecMAT, fstrMAT, istat, conMAT)
This module provides conversion routines between HEC and FISTR data structures and DOF based sparse m...
subroutine, public sparse_matrix_contact_get_rhs(spMAT, hecMAT, fstrMAT)
subroutine, public sparse_matrix_contact_set_rhs(spMAT, hecMAT, fstrMAT)
subroutine, public sparse_matrix_contact_set_vals(spMAT, hecMAT, fstrMAT)
subroutine, public sparse_matrix_para_contact_set_vals(spMAT, hecMAT, fstrMAT, conMAT)
subroutine, public sparse_matrix_para_contact_set_rhs(spMAT, hecMAT, fstrMAT, conMAT)
subroutine, public sparse_matrix_contact_init_prof(spMAT, hecMAT, fstrMAT, hecMESH)
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_symtype_sym
subroutine, public sparse_matrix_finalize(spMAT)
integer(kind=kint), parameter, public sparse_matrix_symtype_asym
subroutine, public sparse_matrix_set_type(spMAT, type, symtype)
Structure for Lagrange multiplier-related part of stiffness matrix (Lagrange multiplier-related matri...