FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_precond_DIAG_11.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 
6 !C
7 !C***
8 !C*** module hecmw_precond_DIAG_11
9 !C***
10 !C
12  use hecmw_util
13 
14  private
15 
19 
20  integer(kind=kint) :: N
21  real(kind=kreal), pointer :: alu(:) => null()
22 
23  logical, save :: INITIALIZED = .false.
24 
25 contains
26 
27  subroutine hecmw_precond_diag_11_setup(hecMAT)
29  implicit none
30  type(hecmwst_matrix), intent(inout) :: hecmat
31  integer(kind=kint ) :: np
32  real (kind=kreal) :: sigma_diag
33  real(kind=kreal), pointer:: d(:)
34 
35  real (kind=kreal):: alutmp(1,1), pw(1)
36  integer(kind=kint ):: ii, i, j, k
37 
38  if (initialized) then
39  if (hecmat%Iarray(98) == 0 .and. hecmat%Iarray(97) == 0) return
41  endif
42 
43  n = hecmat%N
44  np = hecmat%NP
45  d => hecmat%D
46  sigma_diag = hecmw_mat_get_sigma_diag(hecmat)
47 
48  allocate(alu(np))
49  alu = 0.d0
50 
51  do ii= 1, n
52  alu(ii) = d(ii)
53  enddo
54 
55  if (hecmat%cmat%n_val.gt.0) then
56  do k= 1, hecmat%cmat%n_val
57  if (hecmat%cmat%pair(k)%i.ne.hecmat%cmat%pair(k)%j) cycle
58  ii = hecmat%cmat%pair(k)%i
59  alu(ii) = alu(ii) + hecmat%cmat%pair(k)%val(1, 1)
60  enddo
61 
62  !call hecmw_cmat_LU( hecMAT )
63 
64  endif
65 
66  !$omp parallel default(none),private(ii,ALUtmp,k,i,j,PW),shared(N,ALU,SIGMA_DIAG)
67  !$omp do
68  do ii= 1, n
69  alutmp(1,1)= alu(ii) * sigma_diag
70  alutmp(1,1)= 1.d0/alutmp(1,1)
71  alu(ii)= alutmp(1,1)
72  enddo
73  !$omp end do
74  !$omp end parallel
75 
76  initialized = .true.
77  hecmat%Iarray(98) = 0 ! symbolic setup done
78  hecmat%Iarray(97) = 0 ! numerical setup done
79 
80  end subroutine hecmw_precond_diag_11_setup
81 
83  implicit none
84  real(kind=kreal), intent(inout) :: ww(:)
85  integer(kind=kint) :: i
86 
87  !C
88  !C== Block SCALING
89 
90  !$omp parallel default(none),private(i),shared(N,WW,ALU)
91  !$omp do
92  do i= 1, n
93  ww(i)= alu(i)*ww(i)
94  enddo
95  !$omp end do
96  !$omp end parallel
97 
98  end subroutine hecmw_precond_diag_11_apply
99 
101  implicit none
102  if (associated(alu)) deallocate(alu)
103  nullify(alu)
104  initialized = .false.
105  end subroutine hecmw_precond_diag_11_clear
106 
107 end module hecmw_precond_diag_11
real(kind=kreal) function, public hecmw_mat_get_sigma_diag(hecMAT)
subroutine, public hecmw_precond_diag_11_setup(hecMAT)
subroutine, public hecmw_precond_diag_11_apply(WW)
subroutine, public hecmw_precond_diag_11_clear()
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal