FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_precond_nn.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
24
25contains
26
27 subroutine hecmw_precond_nn_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)
36 case(3)
38 case(5)
39 call hecmw_precond_ml_nn_setup(hecmat, hecmesh, sym)
40 case(7)
41 call hecmw_solve_direct_mumps(hecmesh, hecmat)
42 case(10,11,12)
44 case(20)
46 case(21)
48 case default
49 write (*,'(/a )')'#### HEC-MW-SOLVER-E-1001'
50 write (*,'( a/)')' inconsistent solver/preconditioning'
52 end select
53
54
55 end subroutine hecmw_precond_nn_setup
56
57 subroutine hecmw_precond_nn_clear(hecMAT)
58 implicit none
59 type (hecmwst_matrix), intent(inout) :: hecmat
60
61 select case(hecmw_mat_get_precond( hecmat ))
62 case(1,2)
64 case(3)
66 case(5)
68 case(10:12)
70 case(20)
72 case(21)
74 case default
75 end select
76
77 end subroutine hecmw_precond_nn_clear
78
79 !C
80 !C***
81 !C*** hecmw_precond_nn_apply
82 !C***
83 !C
84 subroutine hecmw_precond_nn_apply(hecMESH, hecMAT, R, Z, ZP, time_precond, COMMtime)
85 implicit none
86 type (hecmwst_local_mesh), intent(in) :: hecmesh
87 type (hecmwst_matrix), intent(in) :: hecmat
88 real(kind=kreal), intent(in) :: r(:)
89 real(kind=kreal), intent(out) :: z(:), zp(:)
90 real(kind=kreal), intent(inout) :: time_precond
91 real(kind=kreal), intent(inout) :: commtime
92 integer(kind=kint ) :: i, iterpre, iterpremax
93 real(kind=kreal) :: start_time, end_time
94
95 iterpremax = hecmw_mat_get_iterpremax( hecmat )
96 do iterpre= 1, iterpremax
97 start_time = hecmw_wtime()
98 select case(hecmw_mat_get_precond( hecmat ))
99 case(1,2)
100 call hecmw_precond_ssor_nn_apply(zp,hecmat%NDOF)
101 case(3)
102 call hecmw_precond_diag_nn_apply(zp,hecmat%NDOF)
103 case(5)
105 case(10:12)
106 call hecmw_precond_bilu_nn_apply(zp,hecmat%NDOF)
107 case(20)
109 case(21)
110 call hecmw_precond_rif_nn_apply(zp,hecmat%NDOF)
111 case default
112 end select
113 end_time = hecmw_wtime()
114 time_precond = time_precond + end_time - start_time
115
116 !C-- additive Schwartz
117 do i= 1, hecmat%N * hecmat%NDOF
118 z(i)= z(i) + zp(i)
119 enddo
120 if (iterpre.eq.iterpremax) exit
121
122 !C-- {ZP} = {R} - [A] {Z}
123 call hecmw_matresid_nn (hecmesh, hecmat, z, r, zp, commtime)
124 enddo
125 end subroutine hecmw_precond_nn_apply
126end module hecmw_precond_nn
integer(kind=kint) function, public hecmw_mat_get_iterpremax(hecmat)
integer(kind=kint) function, public hecmw_mat_get_precond(hecmat)
subroutine, public hecmw_precond_bilu_nn_setup(hecmat)
subroutine, public hecmw_precond_bilu_nn_apply(ww, ndof)
subroutine, public hecmw_precond_bilu_nn_clear()
subroutine, public hecmw_precond_diag_nn_clear()
subroutine, public hecmw_precond_diag_nn_setup(hecmat)
subroutine, public hecmw_precond_diag_nn_apply(ww, ndof)
subroutine, public hecmw_precond_ml_nn_setup(hecmat, hecmesh, sym)
subroutine, public hecmw_precond_ml_nn_clear()
subroutine, public hecmw_precond_ml_nn_apply(ww)
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_rif_nn_clear()
subroutine, public hecmw_precond_rif_nn_setup(hecmat)
subroutine, public hecmw_precond_rif_nn_apply(zp, ndof)
subroutine, public hecmw_precond_nn_sainv_apply(r, zp)
subroutine, public hecmw_precond_nn_sainv_clear()
subroutine, public hecmw_precond_nn_sainv_setup(hecmat)
subroutine, public hecmw_precond_ssor_nn_setup(hecmat)
subroutine, public hecmw_precond_ssor_nn_clear(hecmat)
subroutine, public hecmw_precond_ssor_nn_apply(zp, ndof)
This module provides linear equation solver interface for MUMPS.
subroutine, public hecmw_solve_direct_mumps(hecmesh, hecmat)
subroutine, public hecmw_matresid_nn(hecmesh, hecmat, x, b, r, commtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=kint) function hecmw_comm_get_comm()
integer(kind=4), parameter kreal
subroutine hecmw_abort(comm)
real(kind=kreal) function hecmw_wtime()