FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_precond_33.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
17  implicit none
18 
19  private
20 
21  public :: hecmw_precond_33_setup
22  public :: hecmw_precond_33_clear
23  public :: hecmw_precond_33_apply
24 
25 contains
26 
27  subroutine hecmw_precond_33_setup(hecMAT, hecMESH, sym)
28  implicit none
29  type (hecmwst_matrix), intent(inout) :: hecmat
30  type (hecmwst_local_mesh), intent(in) :: hecmesh
31  integer(kind=kint), intent(in) :: sym
32 
33  select case(hecmw_mat_get_precond( hecmat ))
34  case(1,2)
35  call hecmw_precond_ssor_33_setup(hecmat)
36  case(3)
37  call hecmw_precond_diag_33_setup(hecmat)
38  case(5)
39  call hecmw_precond_ml_33_setup(hecmat, hecmesh, sym)
40  case(10,11,12)
41  call hecmw_precond_bilu_33_setup(hecmat)
42  case(20)
43  call hecmw_precond_33_sainv_setup(hecmat)
44  case(21)
45  call hecmw_precond_rif_33_setup(hecmat)
46  case default
47  call hecmw_precond_nn_setup(hecmat, hecmesh, sym)
48  end select
49 
50  end subroutine hecmw_precond_33_setup
51 
52  subroutine hecmw_precond_33_clear(hecMAT)
53  implicit none
54  type (hecmwst_matrix), intent(inout) :: hecmat
55 
56  select case(hecmw_mat_get_precond( hecmat ))
57  case(1,2)
58  call hecmw_precond_ssor_33_clear(hecmat)
59  case(3)
61  case(5)
63  case(10:12)
65  case(20)
67  case(21)
69  case default
70  call hecmw_precond_nn_clear(hecmat)
71  end select
72 
73  end subroutine hecmw_precond_33_clear
74 
75  subroutine hecmw_precond_33_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
76  implicit none
77  type (hecmwst_local_mesh), intent(in) :: hecmesh
78  type (hecmwst_matrix), intent(in) :: hecmat
79  real(kind=kreal), intent(in) :: r(:)
80  real(kind=kreal), intent(out) :: z(:), zp(:)
81  real(kind=kreal), intent(inout) :: time_precond
82  real(kind=kreal), intent(inout) :: commtime
83  integer(kind=kint ) :: i, iterpre, iterpremax
84  real(kind=kreal) :: start_time, end_time
85 
86  iterpremax = hecmw_mat_get_iterpremax( hecmat )
87  do iterpre= 1, iterpremax
88  start_time = hecmw_wtime()
89  select case(hecmw_mat_get_precond( hecmat ))
90  case(1,2)
92  case(3)
94  case(5)
96  case(10:12)
98  case(20)
100  case(21)
102  case default
103  call hecmw_precond_nn_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
104  return
105  end select
106  end_time = hecmw_wtime()
107  time_precond = time_precond + end_time - start_time
108 
109  !C-- additive Schwartz
110  do i= 1, hecmat%N * hecmat%NDOF
111  z(i)= z(i) + zp(i)
112  enddo
113  if (iterpre.eq.iterpremax) exit
114 
115  !C-- {ZP} = {R} - [A] {Z}
116  call hecmw_matresid_33 (hecmesh, hecmat, z, r, zp, commtime)
117  enddo
118  end subroutine hecmw_precond_33_apply
119 end module hecmw_precond_33
integer(kind=kint) function, public hecmw_mat_get_iterpremax(hecMAT)
integer(kind=kint) function, public hecmw_mat_get_precond(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_bilu_33_clear()
subroutine, public hecmw_precond_bilu_33_apply(WW)
subroutine, public hecmw_precond_bilu_33_setup(hecMAT)
subroutine, public hecmw_precond_diag_33_apply(WW)
subroutine, public hecmw_precond_diag_33_clear()
subroutine, public hecmw_precond_diag_33_setup(hecMAT)
subroutine, public hecmw_precond_ml_33_apply(WW)
subroutine, public hecmw_precond_ml_33_clear()
subroutine, public hecmw_precond_ml_33_setup(hecMAT, hecMESH, sym)
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_rif_33_apply(ZP)
subroutine, public hecmw_precond_rif_33_setup(hecMAT)
subroutine, public hecmw_precond_rif_33_clear()
subroutine, public hecmw_precond_33_sainv_setup(hecMAT)
subroutine, public hecmw_precond_33_sainv_apply(R, ZP)
subroutine, public hecmw_precond_33_sainv_clear()
subroutine, public hecmw_precond_ssor_33_apply(ZP)
subroutine, public hecmw_precond_ssor_33_setup(hecMAT)
subroutine, public hecmw_precond_ssor_33_clear(hecMAT)
subroutine, public hecmw_matresid_33(hecMESH, hecMAT, X, B, R, COMMtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
real(kind=kreal) function hecmw_wtime()