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