FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_precond_DIAG_66.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_66
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 contains
24 
25  subroutine hecmw_precond_diag_66_setup(hecMAT)
27  implicit none
28  type(hecmwst_matrix), intent(in) :: hecmat
29  integer(kind=kint ) :: np
30  real (kind=kreal) :: sigma_diag
31  real(kind=kreal), pointer:: d(:)
32 
33  real (kind=kreal):: alutmp(6,6), pw(6)
34  integer(kind=kint ):: ii, i, j, k
35 
36  n = hecmat%N
37  np = hecmat%NP
38  d => hecmat%D
39  sigma_diag = hecmw_mat_get_sigma_diag(hecmat)
40 
41  allocate(alu(36*np))
42  alu = 0.d0
43 
44  do ii= 1, n
45  alu(36*ii-35) = d(36*ii-35)
46  alu(36*ii-34) = d(36*ii-34)
47  alu(36*ii-33) = d(36*ii-33)
48  alu(36*ii-32) = d(36*ii-32)
49  alu(36*ii-31) = d(36*ii-31)
50  alu(36*ii-30) = d(36*ii-30)
51  alu(36*ii-29) = d(36*ii-29)
52  alu(36*ii-28) = d(36*ii-28)
53  alu(36*ii-27) = d(36*ii-27)
54  alu(36*ii-26) = d(36*ii-26)
55  alu(36*ii-25) = d(36*ii-25)
56  alu(36*ii-24) = d(36*ii-24)
57  alu(36*ii-23) = d(36*ii-23)
58  alu(36*ii-22) = d(36*ii-22)
59  alu(36*ii-21) = d(36*ii-21)
60  alu(36*ii-20) = d(36*ii-20)
61  alu(36*ii-19) = d(36*ii-19)
62  alu(36*ii-18) = d(36*ii-18)
63  alu(36*ii-17) = d(36*ii-17)
64  alu(36*ii-16) = d(36*ii-16)
65  alu(36*ii-15) = d(36*ii-15)
66  alu(36*ii-14) = d(36*ii-14)
67  alu(36*ii-13) = d(36*ii-13)
68  alu(36*ii-12) = d(36*ii-12)
69  alu(36*ii-11) = d(36*ii-11)
70  alu(36*ii-10) = d(36*ii-10)
71  alu(36*ii-9 ) = d(36*ii-9 )
72  alu(36*ii-8 ) = d(36*ii-8 )
73  alu(36*ii-7 ) = d(36*ii-7 )
74  alu(36*ii-6 ) = d(36*ii-6 )
75  alu(36*ii-5 ) = d(36*ii-5 )
76  alu(36*ii-4 ) = d(36*ii-4 )
77  alu(36*ii-3 ) = d(36*ii-3 )
78  alu(36*ii-2 ) = d(36*ii-2 )
79  alu(36*ii-1 ) = d(36*ii-1 )
80  alu(36*ii ) = d(36*ii )
81  enddo
82 
83  ! if (hecMAT%cmat%n_val.gt.0) then
84  ! do k= 1, hecMAT%cmat%n_val
85  ! if (hecMAT%cmat%pair(k)%i.ne.hecMAT%cmat%pair(k)%j) cycle
86  ! ii = hecMAT%cmat%pair(k)%i
87  ! ALU(9*ii-8) = ALU(9*ii-8) + hecMAT%cmat%pair(k)%val(1, 1)
88  ! ALU(9*ii-7) = ALU(9*ii-7) + hecMAT%cmat%pair(k)%val(1, 2)
89  ! ALU(9*ii-6) = ALU(9*ii-6) + hecMAT%cmat%pair(k)%val(1, 3)
90  ! ALU(9*ii-5) = ALU(9*ii-5) + hecMAT%cmat%pair(k)%val(2, 1)
91  ! ALU(9*ii-4) = ALU(9*ii-4) + hecMAT%cmat%pair(k)%val(2, 2)
92  ! ALU(9*ii-3) = ALU(9*ii-3) + hecMAT%cmat%pair(k)%val(2, 3)
93  ! ALU(9*ii-2) = ALU(9*ii-2) + hecMAT%cmat%pair(k)%val(3, 1)
94  ! ALU(9*ii-1) = ALU(9*ii-1) + hecMAT%cmat%pair(k)%val(3, 2)
95  ! ALU(9*ii ) = ALU(9*ii ) + hecMAT%cmat%pair(k)%val(3, 3)
96  ! enddo
97  !
98  ! !call hecmw_cmat_LU( hecMAT )
99  !
100  ! endif
101 
102  do ii= 1, n
103  alutmp(1,1)= alu(36*ii-35) * sigma_diag
104  alutmp(1,2)= alu(36*ii-34)
105  alutmp(1,3)= alu(36*ii-33)
106  alutmp(1,4)= alu(36*ii-32)
107  alutmp(1,5)= alu(36*ii-31)
108  alutmp(1,6)= alu(36*ii-30)
109 
110  alutmp(2,1)= alu(36*ii-29)
111  alutmp(2,2)= alu(36*ii-28) * sigma_diag
112  alutmp(2,3)= alu(36*ii-27)
113  alutmp(2,4)= alu(36*ii-26)
114  alutmp(2,5)= alu(36*ii-25)
115  alutmp(2,6)= alu(36*ii-24)
116 
117  alutmp(3,1)= alu(36*ii-23)
118  alutmp(3,2)= alu(36*ii-22)
119  alutmp(3,3)= alu(36*ii-21) * sigma_diag
120  alutmp(3,4)= alu(36*ii-20)
121  alutmp(3,5)= alu(36*ii-19)
122  alutmp(3,6)= alu(36*ii-18)
123 
124  alutmp(4,1)= alu(36*ii-17)
125  alutmp(4,2)= alu(36*ii-16)
126  alutmp(4,3)= alu(36*ii-15)
127  alutmp(4,4)= alu(36*ii-14) * sigma_diag
128  alutmp(4,5)= alu(36*ii-13)
129  alutmp(4,6)= alu(36*ii-12)
130 
131  alutmp(5,1)= alu(36*ii-11)
132  alutmp(5,2)= alu(36*ii-10)
133  alutmp(5,3)= alu(36*ii-9 )
134  alutmp(5,4)= alu(36*ii-8 )
135  alutmp(5,5)= alu(36*ii-7 ) * sigma_diag
136  alutmp(5,6)= alu(36*ii-6 )
137 
138  alutmp(6,1)= alu(36*ii-5 )
139  alutmp(6,2)= alu(36*ii-4 )
140  alutmp(6,3)= alu(36*ii-3 )
141  alutmp(6,4)= alu(36*ii-2 )
142  alutmp(6,5)= alu(36*ii-1 )
143  alutmp(6,6)= alu(36*ii ) * sigma_diag
144 
145  do k= 1, 6
146  alutmp(k,k)= 1.d0/alutmp(k,k)
147  do i= k+1, 6
148  alutmp(i,k)= alutmp(i,k) * alutmp(k,k)
149  do j= k+1, 6
150  pw(j)= alutmp(i,j) - alutmp(i,k)*alutmp(k,j)
151  enddo
152  do j= k+1, 6
153  alutmp(i,j)= pw(j)
154  enddo
155  enddo
156  enddo
157 
158  alu(36*ii-35)= alutmp(1,1)
159  alu(36*ii-34)= alutmp(1,2)
160  alu(36*ii-33)= alutmp(1,3)
161  alu(36*ii-32)= alutmp(1,4)
162  alu(36*ii-31)= alutmp(1,5)
163  alu(36*ii-30)= alutmp(1,6)
164  alu(36*ii-29)= alutmp(2,1)
165  alu(36*ii-28)= alutmp(2,2)
166  alu(36*ii-27)= alutmp(2,3)
167  alu(36*ii-26)= alutmp(2,4)
168  alu(36*ii-25)= alutmp(2,5)
169  alu(36*ii-24)= alutmp(2,6)
170  alu(36*ii-23)= alutmp(3,1)
171  alu(36*ii-22)= alutmp(3,2)
172  alu(36*ii-21)= alutmp(3,3)
173  alu(36*ii-20)= alutmp(3,4)
174  alu(36*ii-19)= alutmp(3,5)
175  alu(36*ii-18)= alutmp(3,6)
176  alu(36*ii-17)= alutmp(4,1)
177  alu(36*ii-16)= alutmp(4,2)
178  alu(36*ii-15)= alutmp(4,3)
179  alu(36*ii-14)= alutmp(4,4)
180  alu(36*ii-13)= alutmp(4,5)
181  alu(36*ii-12)= alutmp(4,6)
182  alu(36*ii-11)= alutmp(5,1)
183  alu(36*ii-10)= alutmp(5,2)
184  alu(36*ii-9 )= alutmp(5,3)
185  alu(36*ii-8 )= alutmp(5,4)
186  alu(36*ii-7 )= alutmp(5,5)
187  alu(36*ii-6 )= alutmp(5,6)
188  alu(36*ii-5 )= alutmp(6,1)
189  alu(36*ii-4 )= alutmp(6,2)
190  alu(36*ii-3 )= alutmp(6,3)
191  alu(36*ii-2 )= alutmp(6,4)
192  alu(36*ii-1 )= alutmp(6,5)
193  alu(36*ii )= alutmp(6,6)
194  enddo
195  end subroutine hecmw_precond_diag_66_setup
196 
198  implicit none
199  real(kind=kreal), intent(inout) :: ww(:)
200  integer(kind=kint) :: i
201  real(kind=kreal) :: x1, x2, x3, x4, x5, x6
202 
203  !C
204  !C== Block SCALING
205  do i= 1, n
206  x1= ww(6*i-5)
207  x2= ww(6*i-4)
208  x3= ww(6*i-3)
209  x4= ww(6*i-2)
210  x5= ww(6*i-1)
211  x6= ww(6*i )
212  x2= x2 -alu(36*i-29)*x1
213  x3= x3 -alu(36*i-23)*x1 -alu(36*i-22)*x2
214  x4= x4 -alu(36*i-17)*x1 -alu(36*i-16)*x2 -alu(36*i-15)*x3
215  x5= x5 -alu(36*i-11)*x1 -alu(36*i-10)*x2 -alu(36*i-9 )*x3 -alu(36*i-8)*x4
216  x6= x6 -alu(36*i-5 )*x1 -alu(36*i-4 )*x2 -alu(36*i-3 )*x3 -alu(36*i-2)*x4 -alu(36*i-1)*x5
217  x6= alu(36*i )* x6
218  x5= alu(36*i-7 )*( x5 -alu(36*i-6 )*x6 )
219  x4= alu(36*i-14)*( x4 -alu(36*i-12)*x6 -alu(36*i-13)*x5)
220  x3= alu(36*i-21)*( x3 -alu(36*i-18)*x6 -alu(36*i-19)*x5 -alu(36*i-20)*x4)
221  x2= alu(36*i-28)*( x2 -alu(36*i-24)*x6 -alu(36*i-25)*x5 -alu(36*i-26)*x4 -alu(36*i-27)*x3)
222  x1= alu(36*i-35)*( x1 -alu(36*i-30)*x6 -alu(36*i-31)*x5 -alu(36*i-32)*x4 -alu(36*i-33)*x3 -alu(36*i-34)*x2)
223  ww(6*i-5)= x1
224  ww(6*i-4)= x2
225  ww(6*i-3)= x3
226  ww(6*i-2)= x4
227  ww(6*i-1)= x5
228  ww(6*i )= x6
229  enddo
230  end subroutine hecmw_precond_diag_66_apply
231 
233  implicit none
234  if (associated(alu)) deallocate(alu)
235  nullify(alu)
236  end subroutine hecmw_precond_diag_66_clear
237 
238 end module hecmw_precond_diag_66
real(kind=kreal) function, public hecmw_mat_get_sigma_diag(hecMAT)
subroutine, public hecmw_precond_diag_66_clear()
subroutine, public hecmw_precond_diag_66_setup(hecMAT)
subroutine, public hecmw_precond_diag_66_apply(WW)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal