FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_precond_44.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
14  implicit none
15 
16  private
17 
18  public :: hecmw_precond_44_setup
19  public :: hecmw_precond_44_clear
20  public :: hecmw_precond_44_apply
21 
22 contains
23 
24  subroutine hecmw_precond_44_setup(hecMAT, hecMESH, sym)
25  implicit none
26  type (hecmwst_matrix), intent(inout) :: hecmat
27  type (hecmwst_local_mesh), intent(inout) :: hecmesh
28  integer(kind=kint), intent(in) :: sym
29 
30  select case(hecmw_mat_get_precond( hecmat ))
31  case(1,2)
32  call hecmw_precond_ssor_44_setup(hecmat)
33  case(3)
34  call hecmw_precond_diag_44_setup(hecmat)
35  case(10,11,12)
36  call hecmw_precond_bilu_44_setup(hecmat)
37  case default
38  call hecmw_precond_nn_setup(hecmat, hecmesh, sym)
39  end select
40 
41  end subroutine hecmw_precond_44_setup
42 
43  subroutine hecmw_precond_44_clear(hecMAT)
44  implicit none
45  type (hecmwst_matrix), intent(inout) :: hecmat
46 
47  select case(hecmw_mat_get_precond( hecmat ))
48  case(1,2)
49  call hecmw_precond_ssor_44_clear(hecmat)
50  case(3)
52  case(10,11,12)
54  case default
55  call hecmw_precond_nn_clear(hecmat)
56  end select
57 
58  end subroutine hecmw_precond_44_clear
59 
60  subroutine hecmw_precond_44_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
61  implicit none
62  type (hecmwst_local_mesh), intent(in) :: hecmesh
63  type (hecmwst_matrix), intent(in) :: hecmat
64  real(kind=kreal), intent(in) :: r(:)
65  real(kind=kreal), intent(out) :: z(:), zp(:)
66  real(kind=kreal), intent(inout) :: time_precond
67  real(kind=kreal), intent(inout) :: commtime
68  integer(kind=kint ) :: i, iterpre, iterpremax
69  real(kind=kreal) :: start_time, end_time
70 
71  iterpremax = hecmw_mat_get_iterpremax( hecmat )
72  do iterpre= 1, iterpremax
73  start_time = hecmw_wtime()
74  select case(hecmw_mat_get_precond( hecmat ))
75  case(1,2)
77  case(3)
79  case(10,11,12)
81  case default
82  call hecmw_precond_nn_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
83  return
84  end select
85  end_time = hecmw_wtime()
86  time_precond = time_precond + end_time - start_time
87 
88  !C-- additive Schwartz
89  do i= 1, hecmat%N * hecmat%NDOF
90  z(i)= z(i) + zp(i)
91  enddo
92  if (iterpre.eq.iterpremax) exit
93 
94  !C-- {ZP} = {R} - [A] {Z}
95  call hecmw_matresid_44 (hecmesh, hecmat, z, r, zp, commtime)
96  enddo
97  end subroutine hecmw_precond_44_apply
98 end module hecmw_precond_44
integer(kind=kint) function, public hecmw_mat_get_iterpremax(hecMAT)
integer(kind=kint) function, public hecmw_mat_get_precond(hecMAT)
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_bilu_44_apply(WW)
subroutine, public hecmw_precond_bilu_44_clear()
subroutine, public hecmw_precond_bilu_44_setup(hecMAT)
subroutine, public hecmw_precond_diag_44_apply(WW)
subroutine, public hecmw_precond_diag_44_setup(hecMAT)
subroutine, public hecmw_precond_diag_44_clear()
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_ssor_44_clear(hecMAT)
subroutine, public hecmw_precond_ssor_44_apply(ZP)
subroutine, public hecmw_precond_ssor_44_setup(hecMAT)
subroutine, public hecmw_matresid_44(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()