FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_precond_DIAG_44.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_44
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_44_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(4,4), pw(4)
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(16*np))
49  alu = 0.d0
50 
51  do ii= 1, n
52  alu(16*ii-15) = d(16*ii-15)
53  alu(16*ii-14) = d(16*ii-14)
54  alu(16*ii-13) = d(16*ii-13)
55  alu(16*ii-12) = d(16*ii-12)
56  alu(16*ii-11) = d(16*ii-11)
57  alu(16*ii-10) = d(16*ii-10)
58  alu(16*ii-9 ) = d(16*ii-9 )
59  alu(16*ii-8 ) = d(16*ii-8 )
60  alu(16*ii-7 ) = d(16*ii-7 )
61  alu(16*ii-6 ) = d(16*ii-6 )
62  alu(16*ii-5 ) = d(16*ii-5 )
63  alu(16*ii-4 ) = d(16*ii-4 )
64  alu(16*ii-3 ) = d(16*ii-3 )
65  alu(16*ii-2 ) = d(16*ii-2 )
66  alu(16*ii-1 ) = d(16*ii-1 )
67  alu(16*ii ) = d(16*ii )
68  enddo
69 
70  do ii= 1, n
71  alutmp(1,1)= alu(16*ii-15) * sigma_diag
72  alutmp(1,2)= alu(16*ii-14)
73  alutmp(1,3)= alu(16*ii-13)
74  alutmp(1,4)= alu(16*ii-12)
75 
76  alutmp(2,1)= alu(16*ii-11)
77  alutmp(2,2)= alu(16*ii-10) * sigma_diag
78  alutmp(2,3)= alu(16*ii- 9)
79  alutmp(2,4)= alu(16*ii- 8)
80 
81  alutmp(3,1)= alu(16*ii- 7)
82  alutmp(3,2)= alu(16*ii- 6)
83  alutmp(3,3)= alu(16*ii- 5) * sigma_diag
84  alutmp(3,4)= alu(16*ii- 4)
85 
86  alutmp(4,1)= alu(16*ii- 3)
87  alutmp(4,2)= alu(16*ii- 2)
88  alutmp(4,3)= alu(16*ii- 1)
89  alutmp(4,4)= alu(16*ii ) * sigma_diag
90 
91  do k= 1, 4
92  alutmp(k,k)= 1.d0/alutmp(k,k)
93  do i= k+1, 4
94  alutmp(i,k)= alutmp(i,k) * alutmp(k,k)
95  do j= k+1, 4
96  pw(j)= alutmp(i,j) - alutmp(i,k)*alutmp(k,j)
97  enddo
98  do j= k+1, 4
99  alutmp(i,j)= pw(j)
100  enddo
101  enddo
102  enddo
103 
104  alu(16*ii-15)= alutmp(1,1)
105  alu(16*ii-14)= alutmp(1,2)
106  alu(16*ii-13)= alutmp(1,3)
107  alu(16*ii-12)= alutmp(1,4)
108  alu(16*ii-11)= alutmp(2,1)
109  alu(16*ii-10)= alutmp(2,2)
110  alu(16*ii- 9)= alutmp(2,3)
111  alu(16*ii- 8)= alutmp(2,4)
112  alu(16*ii- 7)= alutmp(3,1)
113  alu(16*ii- 6)= alutmp(3,2)
114  alu(16*ii- 5)= alutmp(3,3)
115  alu(16*ii- 4)= alutmp(3,4)
116  alu(16*ii- 3)= alutmp(4,1)
117  alu(16*ii- 2)= alutmp(4,2)
118  alu(16*ii- 1)= alutmp(4,3)
119  alu(16*ii )= alutmp(4,4)
120  enddo
121 
122  initialized = .true.
123  hecmat%Iarray(98) = 0 ! symbolic setup done
124  hecmat%Iarray(97) = 0 ! numerical setup done
125  end subroutine hecmw_precond_diag_44_setup
126 
128  implicit none
129  real(kind=kreal), intent(inout) :: ww(:)
130  integer(kind=kint) :: i
131  real(kind=kreal) :: x1, x2, x3, x4
132 
133  !C
134  !C== Block SCALING
135  do i= 1, n
136  x1= ww(4*i-3)
137  x2= ww(4*i-2)
138  x3= ww(4*i-1)
139  x4= ww(4*i )
140  x2= x2 -alu(16*i-11)*x1
141  x3= x3 -alu(16*i- 7)*x1 -alu(16*i- 6)*x2
142  x4= x4 -alu(16*i- 3)*x1 -alu(16*i- 2)*x2 -alu(16*i- 1)*x3
143  x4= alu(16*i )* x4
144  x3= alu(16*i- 5)*(x3 -alu(16*i- 4)*x4)
145  x2= alu(16*i-10)*(x2 -alu(16*i- 8)*x4 -alu(16*i- 9)*x3)
146  x1= alu(16*i-15)*(x1 -alu(16*i-12)*x4 -alu(16*i-13)*x3 -alu(16*i-14)*x2)
147  ww(4*i-3)= x1
148  ww(4*i-2)= x2
149  ww(4*i-1)= x3
150  ww(4*i )= x4
151  enddo
152  end subroutine hecmw_precond_diag_44_apply
153 
155  implicit none
156  if (associated(alu)) deallocate(alu)
157  nullify(alu)
158  initialized = .false.
159  end subroutine hecmw_precond_diag_44_clear
160 
161 end module hecmw_precond_diag_44
real(kind=kreal) function, public hecmw_mat_get_sigma_diag(hecMAT)
subroutine, public hecmw_precond_diag_44_apply(WW)
subroutine, public hecmw_precond_diag_44_setup(hecMAT)
subroutine, public hecmw_precond_diag_44_clear()
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal