FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_precond_22.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
9 !use hecmw_precond_BILU_22
14 implicit none
15
16 private
17
21
22contains
23
24 subroutine hecmw_precond_22_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)
33 case(3)
35 !CASE(10,11,12)
36 ! call hecmw_precond_BILU_22_setup(hecMAT)
37 case default
38 call hecmw_precond_nn_setup(hecmat, hecmesh, sym)
39 end select
40
41 end subroutine hecmw_precond_22_setup
42
43 subroutine hecmw_precond_22_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)
50 case(3)
52 !CASE(10,11,12)
53 ! call hecmw_precond_BILU_22_clear()
54 case default
55 call hecmw_precond_nn_clear(hecmat)
56 end select
57
58 end subroutine hecmw_precond_22_clear
59
60 subroutine hecmw_precond_22_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)
80 ! call hecmw_precond_BILU_22_apply(ZP)
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_22 (hecmesh, hecmat, z, r, zp, commtime)
96 enddo
97 end subroutine hecmw_precond_22_apply
98end module hecmw_precond_22
integer(kind=kint) function, public hecmw_mat_get_iterpremax(hecmat)
integer(kind=kint) function, public hecmw_mat_get_precond(hecmat)
subroutine, public hecmw_precond_22_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
subroutine, public hecmw_precond_22_setup(hecmat, hecmesh, sym)
subroutine, public hecmw_precond_22_clear(hecmat)
subroutine, public hecmw_precond_diag_22_clear()
subroutine, public hecmw_precond_diag_22_apply(ww)
subroutine, public hecmw_precond_diag_22_setup(hecmat)
subroutine, public hecmw_precond_nn_clear(hecmat)
subroutine, public hecmw_precond_nn_apply(hecmesh, hecmat, r, z, zp, time_precond, commtime)
subroutine, public hecmw_precond_nn_setup(hecmat, hecmesh, sym)
subroutine, public hecmw_precond_ssor_22_clear(hecmat)
subroutine, public hecmw_precond_ssor_22_setup(hecmat)
subroutine, public hecmw_precond_ssor_22_apply(zp)
subroutine, public hecmw_matresid_22(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()