FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_precond.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 !-------------------------------------------------------------------------------
5 
7  use hecmw_util
15  implicit none
16 
17  private
18  public :: hecmw_precond_setup
19  public :: hecmw_precond_clear
20  public :: hecmw_precond_apply
22  public :: hecmw_precond_get_timer
23 
24  real(kind=kreal) :: time_precond = 0.d0
25 
26 contains
27 
28  subroutine hecmw_precond_setup(hecMAT, hecMESH, sym)
29  implicit none
30  type (hecmwst_matrix), intent(inout) :: hecmat
31  type (hecmwst_local_mesh), intent(inout) :: hecmesh
32  integer(kind=kint) :: sym
33 
34  if (hecmw_mat_get_iterpremax( hecmat ).le.0) return
35 
36  select case(hecmat%NDOF)
37  case(3)
38  call hecmw_precond_33_setup(hecmat, hecmesh, sym)
39  case(4)
40  call hecmw_precond_44_setup(hecmat, hecmesh, sym)
41  case(6)
42  call hecmw_precond_66_setup(hecmat, hecmesh, sym)
43  case(1)
44  call hecmw_precond_11_setup(hecmat, hecmesh, sym)
45  case(2)
46  call hecmw_precond_22_setup(hecmat, hecmesh, sym)
47  case default
48  call hecmw_precond_nn_setup(hecmat, hecmesh, sym)
49  end select
50  end subroutine hecmw_precond_setup
51 
52  subroutine hecmw_precond_clear(hecMAT)
53  implicit none
54  type (hecmwst_matrix), intent(inout) :: hecmat
55 
56  if (hecmw_mat_get_iterpremax( hecmat ).le.0) return
57 
58  select case(hecmat%NDOF)
59  case(3)
60  call hecmw_precond_33_clear(hecmat)
61  case(4)
62  call hecmw_precond_44_clear(hecmat)
63  case(6)
64  call hecmw_precond_66_clear(hecmat)
65  case(1)
66  call hecmw_precond_11_clear(hecmat)
67  case(2)
68  call hecmw_precond_22_clear(hecmat)
69  case default
70  call hecmw_precond_nn_clear(hecmat)
71  end select
72 
73  end subroutine hecmw_precond_clear
74 
75  subroutine hecmw_precond_apply(hecMESH, hecMAT, R, Z, ZP, COMMtime)
76  implicit none
77  type (hecmwst_local_mesh), intent(inout) :: hecmesh
78  type (hecmwst_matrix), intent(inout) :: hecmat
79  real(kind=kreal), intent(inout) :: r(:)
80  real(kind=kreal), intent(inout) :: z(:), zp(:)
81  real(kind=kreal), intent(inout) :: commtime
82  integer(kind=kint ) :: i, n, np, nndof, npndof
83  n = hecmat%N
84  np = hecmat%NP
85  nndof = n * hecmat%NDOF
86  npndof = np * hecmat%NDOF
87 
88  if (hecmw_mat_get_iterpremax( hecmat ).le.0) then
89  do i= 1, nndof
90  z(i)= r(i)
91  enddo
92  return
93  endif
94 
95  !C {z}= [Minv]{r}
96  do i= 1, nndof
97  zp(i)= r(i)
98  enddo
99  do i= nndof+1, npndof
100  zp(i) = 0.d0
101  enddo
102  do i= 1, npndof
103  z(i)= 0.d0
104  enddo
105 
106  select case(hecmat%NDOF)
107  case(3)
108  call hecmw_precond_33_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
109  case(4)
110  call hecmw_precond_44_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
111  case(6)
112  call hecmw_precond_66_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
113  case(1)
114  call hecmw_precond_11_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
115  case(2)
116  call hecmw_precond_22_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
117  case default
118  call hecmw_precond_nn_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
119  end select
120 
121  end subroutine hecmw_precond_apply
122 
124  implicit none
125  time_precond = 0.d0
126  end subroutine hecmw_precond_clear_timer
128  implicit none
129  real(kind=kreal) :: hecmw_precond_get_timer
130  hecmw_precond_get_timer = time_precond
131  end function hecmw_precond_get_timer
132 end module hecmw_precond
integer(kind=kint) function, public hecmw_mat_get_iterpremax(hecMAT)
subroutine, public hecmw_precond_11_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_11_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
subroutine, public hecmw_precond_11_clear(hecMAT)
subroutine, public hecmw_precond_22_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_22_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
subroutine, public hecmw_precond_22_clear(hecMAT)
subroutine, public hecmw_precond_33_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
subroutine, public hecmw_precond_33_clear(hecMAT)
subroutine, public hecmw_precond_33_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_44_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_44_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
subroutine, public hecmw_precond_44_clear(hecMAT)
subroutine, public hecmw_precond_66_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_66_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
subroutine, public hecmw_precond_66_clear(hecMAT)
subroutine, public hecmw_precond_nn_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_nn_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
subroutine, public hecmw_precond_nn_clear(hecMAT)
subroutine, public hecmw_precond_clear_timer
real(kind=kreal) function, public hecmw_precond_get_timer()
subroutine, public hecmw_precond_clear(hecMAT)
subroutine, public hecmw_precond_setup(hecMAT, hecMESH, sym)
subroutine, public hecmw_precond_apply(hecMESH, hecMAT, R, Z, ZP, COMMtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal