FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_precond_DIAG_22.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_22
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
25contains
26
27 subroutine hecmw_precond_diag_22_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(2,2), pw(2)
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(4*np))
49 alu = 0.d0
50
51 do ii= 1, n
52 alu(4*ii-3) = d(4*ii-3)
53 alu(4*ii-2) = d(4*ii-2)
54 alu(4*ii-1) = d(4*ii-1)
55 alu(4*ii-0) = d(4*ii-0)
56 enddo
57
58 if (hecmat%cmat%n_val.gt.0) then
59 do k= 1, hecmat%cmat%n_val
60 if (hecmat%cmat%pair(k)%i.ne.hecmat%cmat%pair(k)%j) cycle
61 ii = hecmat%cmat%pair(k)%i
62 alu(4*ii-3) = alu(4*ii-3) + hecmat%cmat%pair(k)%val(1, 1)
63 alu(4*ii-2) = alu(4*ii-2) + hecmat%cmat%pair(k)%val(1, 2)
64 alu(4*ii-1) = alu(4*ii-1) + hecmat%cmat%pair(k)%val(2, 1)
65 alu(4*ii-0) = alu(4*ii-0) + hecmat%cmat%pair(k)%val(2, 2)
66 enddo
67
68 !call hecmw_cmat_LU( hecMAT )
69
70 endif
71
72 !$omp parallel default(none),private(ii,ALUtmp,k,i,j,PW),shared(N,ALU,SIGMA_DIAG)
73 !$omp do
74 do ii= 1, n
75 alutmp(1,1)= alu(4*ii-3) * sigma_diag
76 alutmp(1,2)= alu(4*ii-2)
77 alutmp(2,1)= alu(4*ii-1)
78 alutmp(2,2)= alu(4*ii-0) * sigma_diag
79
80 do k= 1, 2
81 alutmp(k,k)= 1.d0/alutmp(k,k)
82 do i= k+1, 2
83 alutmp(i,k)= alutmp(i,k) * alutmp(k,k)
84 do j= k+1, 2
85 pw(j)= alutmp(i,j) - alutmp(i,k)*alutmp(k,j)
86 enddo
87 do j= k+1, 2
88 alutmp(i,j)= pw(j)
89 enddo
90 enddo
91 enddo
92 alu(4*ii-3)= alutmp(1,1)
93 alu(4*ii-2)= alutmp(1,2)
94 alu(4*ii-1)= alutmp(2,1)
95 alu(4*ii-0)= alutmp(2,2)
96 enddo
97 !$omp end do
98 !$omp end parallel
99
100 initialized = .true.
101 hecmat%Iarray(98) = 0 ! symbolic setup done
102 hecmat%Iarray(97) = 0 ! numerical setup done
103
104 end subroutine hecmw_precond_diag_22_setup
105
107 implicit none
108 real(kind=kreal), intent(inout) :: ww(:)
109 integer(kind=kint) :: i
110 real(kind=kreal) :: x1, x2
111
112 !C
113 !C== Block SCALING
114
115 !$omp parallel default(none),private(i,X1,X2),shared(N,WW,ALU)
116 !$omp do
117 do i= 1, n
118 x1= ww(2*i-1)
119 x2= ww(2*i-0)
120 x2= x2 - alu(4*i-1)*x1
121 x2= alu(4*i )* x2
122 x1= alu(4*i-3)*( x1 - alu(4*i-2)*x2 )
123 ww(2*i-1)= x1
124 ww(2*i-0)= x2
125 enddo
126 !$omp end do
127 !$omp end parallel
128
129 end subroutine hecmw_precond_diag_22_apply
130
132 implicit none
133 if (associated(alu)) deallocate(alu)
134 nullify(alu)
135 initialized = .false.
136 end subroutine hecmw_precond_diag_22_clear
137
138end module hecmw_precond_diag_22
real(kind=kreal) function, public hecmw_mat_get_sigma_diag(hecmat)
subroutine, public hecmw_precond_diag_22_clear()
subroutine, public hecmw_precond_diag_22_apply(ww)
subroutine, public hecmw_precond_diag_22_setup(hecmat)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal