1!------------------------------------------------------------------------------- 2! Copyright (c) 2019 FrontISTR Commons 3! This software is released under the MIT License, see LICENSE.txt 4!------------------------------------------------------------------------------- 5 6module hecmw_precond_33 7 use hecmw_util 8 use hecmw_matrix_misc 9 use hecmw_precond_BILU_33 10 use hecmw_precond_DIAG_33 11 use hecmw_precond_SSOR_33 12 use hecmw_precond_ML_33 13 use hecmw_precond_SAINV_33 14 use hecmw_precond_RIF_33 15 use hecmw_precond_nn 16 use hecmw_solver_las_33 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 25contains 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) 60 call hecmw_precond_DIAG_33_clear() 61 case(5) 62 call hecmw_precond_ML_33_clear() 63 case(10:12) 64 call hecmw_precond_BILU_33_clear() 65 case(20) 66 call hecmw_precond_33_SAINV_clear() 67 case(21) 68 call hecmw_precond_RIF_33_clear() 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) 91 call hecmw_precond_SSOR_33_apply(ZP) 92 case(3) 93 call hecmw_precond_DIAG_33_apply(ZP) 94 case(5) 95 call hecmw_precond_ML_33_apply(ZP) 96 case(10:12) 97 call hecmw_precond_BILU_33_apply(ZP) 98 case(20) 99 call hecmw_precond_33_SAINV_apply(R,ZP) 100 case(21) 101 call hecmw_precond_RIF_33_apply(ZP) 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 119end module hecmw_precond_33 120