FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
solve_LINEQ_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 
10 
11  use m_fstr
17 
18  implicit none
19 
20  private
21  public :: solve_lineq_contact_init
22  public :: solve_lineq_contact
23 
24 contains
25 
27  subroutine solve_lineq_contact_init(hecMESH,hecMAT,fstrMAT,is_sym)
28  type (hecmwst_local_mesh) :: hecmesh
29  type (hecmwst_matrix) :: hecmat
30  type (fstrst_matrix_contact_lagrange) :: fstrmat
31  logical :: is_sym
32 
33  if( hecmat%Iarray(99)==1 )then
34  call solve_lineq_iter_contact_init(hecmesh,hecmat,fstrmat,is_sym)
35  elseif( hecmat%Iarray(99)==2 )then
36  call solve_lineq_serial_lag_hecmw_init(hecmat,fstrmat,is_sym)
37  else if( hecmat%Iarray(99)==3 )then
38  call solve_lineq_mkl_contact_init(hecmesh,is_sym)
39  elseif( hecmat%Iarray(99)==5 ) then
40  call solve_lineq_mumps_contact_init(hecmesh,hecmat,fstrmat,is_sym)
41  endif
42  end subroutine solve_lineq_contact_init
43 
44 
46  subroutine solve_lineq_contact(hecMESH,hecMAT,fstrMAT,istat,rf,conMAT)
47 
48  type (hecmwst_local_mesh) :: hecmesh
49  type (hecmwst_matrix) :: hecmat
50  type (fstrst_matrix_contact_lagrange) :: fstrmat
51  integer(kind=kint), intent(out) :: istat
52  real(kind=kreal), optional :: rf
53  type (hecmwst_matrix),optional :: conmat
54 
55  real(kind=kreal) :: factor
56  real(kind=kreal) :: resid
57  real(kind=kreal) :: t1, t2
58  integer(kind=kint) :: ndof
59 
60  factor = 1.0d0
61  if( present(rf) )factor = rf
62 
63  t1 = hecmw_wtime()
64 
65  istat = 0
66  if( hecmat%Iarray(99)==1 )then
67  if(paracontactflag.and.present(conmat)) then
68  call solve_lineq_iter_contact(hecmesh,hecmat,fstrmat,istat,conmat)
69  else
70  call solve_lineq_iter_contact(hecmesh,hecmat,fstrmat,istat)
71  endif
72  elseif( hecmat%Iarray(99)==2 )then
73  call solve_lineq_serial_lag_hecmw(hecmesh,hecmat,fstrmat)
74  elseif( hecmat%Iarray(99)==3 )then
75  if(paracontactflag.and.present(conmat)) then
76  call solve_lineq_mkl_contact(hecmesh,hecmat,fstrmat,istat,conmat)
77  else
78  call solve_lineq_mkl_contact(hecmesh,hecmat,fstrmat,istat)
79  endif
80  elseif( hecmat%Iarray(99)==5 ) then
81  ! ---- For Parallel Contact with Multi-Partition Domains
82  if(paracontactflag.and.present(conmat)) then
83  call solve_lineq_mumps_contact(hecmesh,hecmat,fstrmat,istat,conmat)
84  else
85  call solve_lineq_mumps_contact(hecmesh,hecmat,fstrmat,istat)
86  endif
87  endif
88 
89  ndof = hecmat%NDOF
90  if( ndof==3 ) then
91  call hecmw_update_3_r(hecmesh,hecmat%X,hecmesh%n_node)
92  else if( ndof==2 ) then
93  call hecmw_update_2_r(hecmesh,hecmat%X,hecmesh%n_node)
94  else if( ndof==4 ) then
95  call hecmw_update_4_r(hecmesh,hecmat%X,hecmesh%n_node)
96  else if( ndof==6 ) then
97  call hecmw_update_m_r(hecmesh,hecmat%X,hecmesh%n_node,6)
98  endif
99 
100  t2 = hecmw_wtime()
101  if (hecmw_mat_get_timelog(hecmat) .ge. 1) then
102  if (myrank==0) write(*,*) ' solve time :', t2 - t1
103  endif
104 
105  if(paracontactflag.and.present(conmat)) then
106  else
107  resid=fstr_get_resid_max_contact(hecmesh,hecmat,fstrmat)
108  if (myrank==0) then
109  write(*,*) ' maximum residual = ', resid
110  if( hecmw_mat_get_solver_type(hecmat) /= 1 .and. resid >= 1.0d-8) then
111  write(*,*) ' ###Maximum residual exceeded 1.0d-8---Direct Solver### '
112  ! stop
113  endif
114  endif
115  endif
116 
117  hecmat%X=factor*hecmat%X
118 
119  end subroutine solve_lineq_contact
120 
121 end module m_solve_lineq_contact
real(kind=kreal) function, public fstr_get_resid_max_contact(hecMESH, hecMAT, fstrMAT)
This function calculates maximum residual.
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 functions to solve sparse system of \linear equitions in the case of contact ana...
subroutine, public solve_lineq_contact_init(hecMESH, hecMAT, fstrMAT, is_sym)
This subroutine.
subroutine, public solve_lineq_contact(hecMESH, hecMAT, fstrMAT, istat, rf, conMAT)
This subroutine.
subroutine solve_lineq_serial_lag_hecmw(hecMESH, hecMAT, fstrMAT)
subroutine solve_lineq_serial_lag_hecmw_init(hecMAT, fstrMAT, is_sym)
This module provides interface of iteratie linear equation solver for contact problems using Lagrange...
subroutine, public solve_lineq_iter_contact_init(hecMESH, hecMAT, fstrMAT, is_sym)
subroutine, public solve_lineq_iter_contact(hecMESH, hecMAT, fstrMAT, istat, conMAT)
This module provides functions to solve sparse system of \linear equitions using intel MKL direct spa...
subroutine, public solve_lineq_mkl_contact(hecMESH, hecMAT, fstrMAT, istat, conMAT)
This subroutine executes the MKL solver.
subroutine, public solve_lineq_mkl_contact_init(hecMESH, is_sym)
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)