FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
hecmw_precond_RIF_nn.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 module hecmw_precond_rif_nn
7  use hecmw_util
8  use m_hecmw_comm_f
9  use hecmw_matrix_misc
10 
11  private
12 
13  public:: hecmw_precond_rif_nn_setup
14  public:: hecmw_precond_rif_nn_apply
15  public:: hecmw_precond_rif_nn_clear
16 
17  integer(4),parameter :: krealp = 8
18 
19  integer(kind=kint) :: NPFIU, NPFIL
20  integer(kind=kint) :: N
21  integer(kind=kint), pointer :: inumFI1L(:) => null()
22  integer(kind=kint), pointer :: inumFI1U(:) => null()
23  integer(kind=kint), pointer :: FI1L(:) => null()
24  integer(kind=kint), pointer :: FI1U(:) => null()
25 
26  integer(kind=kint), pointer :: indexL(:) => null()
27  integer(kind=kint), pointer :: indexU(:) => null()
28  integer(kind=kint), pointer :: itemL(:) => null()
29  integer(kind=kint), pointer :: itemU(:) => null()
30  real(kind=kreal), pointer :: d(:) => null()
31  real(kind=kreal), pointer :: al(:) => null()
32  real(kind=kreal), pointer :: au(:) => null()
33 
34  real(kind=krealp), pointer :: sainvl(:) => null()
35  real(kind=krealp), pointer :: sainvd(:) => null()
36  real(kind=krealp), pointer :: rifu(:) => null()
37  real(kind=krealp), pointer :: rifl(:) => null()
38  real(kind=krealp), pointer :: rifd(:) => null()
39 
40 contains
41 
42  !C***
43  !C*** hecmw_precond_nn_sainv_setup
44  !C***
45  subroutine hecmw_precond_rif_nn_setup(hecMAT)
46  implicit none
47  type(hecmwst_matrix) :: hecmat
48 
49  integer(kind=kint ) :: precond, ndof, ndof2
50 
51  real(kind=krealp) :: filter
52 
53  n = hecmat%N
54  ndof = hecmat%NDOF
55  ndof2 = ndof*ndof
56  precond = hecmw_mat_get_precond(hecmat)
57 
58  d => hecmat%D
59  au=> hecmat%AU
60  al=> hecmat%AL
61  indexl => hecmat%indexL
62  indexu => hecmat%indexU
63  iteml => hecmat%itemL
64  itemu => hecmat%itemU
65 
66  if (precond.eq.21) call form_ilu1_rif_nn(hecmat)
67 
68  allocate (sainvd(ndof2*hecmat%NP))
69  allocate (sainvl(ndof2*npfiu))
70  allocate (rifd(ndof2*hecmat%NP))
71  allocate (rifu(ndof2*npfiu))
72  sainvd = 0.0d0
73  sainvl = 0.0d0
74  rifd = 0.0d0
75  rifu = 0.0d0
76 
77  filter= hecmat%Rarray(5)
78 
79  write(*,"(a,F15.8)")"### RIF FILTER :",filter
80 
81  call hecmw_rif_nn(hecmat)
82 
83  allocate (rifl(ndof2*npfiu))
84  rifl = 0.0d0
85 
86  call hecmw_rif_make_u_nn(hecmat)
87 
88  end subroutine hecmw_precond_rif_nn_setup
89 
90  subroutine hecmw_precond_rif_nn_apply(ZP, NDOF)
91  implicit none
92  real(kind=kreal), intent(inout) :: zp(:)
93  integer(kind=kint) :: in, i, j, isl, iel,ndof,ndof2,idof,jdof
94  real(kind=kreal) :: sw(ndof),x(ndof)
95  ndof2=ndof*ndof
96  !C-- FORWARD
97  do i= 1, n
98  do idof = 1, ndof
99  sw(idof) = zp(ndof*(i-1)+idof)
100  end do
101  isl= inumfi1l(i-1)+1
102  iel= inumfi1l(i)
103  do j= isl, iel
104  in= fi1l(j)
105  do idof = 1, ndof
106  x(idof) = zp(ndof*(in-1)+idof)
107  end do
108  do idof = 1, ndof
109  do jdof = 1, ndof
110  sw(idof) = sw(idof) - rifl(ndof2*(j-1)+ndof*(idof-1)+jdof)*x(jdof)
111  end do
112  end do
113  enddo
114  x = sw
115  do idof = 2,ndof
116  do jdof = 1, idof-1
117  x(idof) = x(idof) - rifd(ndof2*(i-1)+ndof*(idof-1)+jdof )*x(jdof)
118  end do
119  end do
120  zp(ndof*(i-1)+1:ndof*(i-1)+ndof) = x(1:ndof)
121  enddo
122 
123  do i=1, n
124  do idof=1,ndof
125  zp(ndof*(i-1)+idof)= zp(ndof*(i-1)+idof)*rifd(ndof2*(i-1)+(idof-1)*ndof+idof)
126  end do
127  enddo
128  !C-- BACKWARD
129  do i= n, 1, -1
130  do idof = 1, ndof
131  sw(idof) = zp(ndof*(i-1)+idof)
132  end do
133  isl= inumfi1u(i-1) + 1
134  iel= inumfi1u(i)
135  do j= iel, isl, -1
136  in= fi1u(j)
137  do idof = 1, ndof
138  x(idof) = zp(ndof*(in-1)+idof)
139  end do
140  do idof = 1, ndof
141  do jdof = 1, ndof
142  sw(idof) = sw(idof) - rifu(ndof2*(j-1)+ndof*(idof-1)+jdof)*x(jdof)
143  end do
144  end do
145  enddo
146  x=sw
147  do idof = ndof, 1, -1
148  do jdof = ndof, idof+1, -1
149  x(idof) = x(idof) - rifd(ndof*ndof*(i-1)+ndof*(jdof-1)+idof)*x(jdof)
150  end do
151  end do
152  do idof = 1, ndof
153  zp(ndof*(i-1)+idof) = x(idof)
154  end do
155  enddo
156  end subroutine hecmw_precond_rif_nn_apply
157 
158 
159  !C***
160  !C*** hecmw_rif_nn
161  !C***
162  subroutine hecmw_rif_nn(hecMAT)
163  implicit none
164  type (hecmwst_matrix) :: hecmat
165  integer(kind=kint) :: i, j, k, js, je, in, itr, np, ndof, ndof2, idof, jdof, iitr
166  real(kind=krealp) :: dd, dtmp(hecmat%NDOF), x(hecmat%NDOF)
167  real(kind=krealp) :: filter
168  real(kind=krealp), allocatable :: zz(:), vv(:)
169 
170  filter= hecmat%Rarray(5)
171 
172  np = hecmat%NP
173  ndof = hecmat%NDOF
174  ndof2= ndof*ndof
175  allocate(vv(ndof*np))
176  allocate(zz(ndof*np))
177 
178  do itr=1,np
179  do iitr=1,ndof
180  zz(:) = 0.0d0
181  vv(:) = 0.0d0
182  !{v}=[A]{zi}
183  do idof = 1,ndof
184  zz(ndof*(itr-1)+idof)= sainvd(ndof2*(itr-1)+ndof*(idof-1)+iitr)
185  end do
186  zz(ndof*(itr-1)+iitr)= 1.0d0
187 
188  js= inumfi1l(itr-1) + 1
189  je= inumfi1l(itr )
190  do j= js, je
191  in = fi1l(j)
192  do idof = 1, ndof
193  zz(ndof*(in-1)+idof)=sainvl(ndof2*(j-1)+ndof*(iitr-1)+idof)
194  end do
195  enddo
196 
197  do i= 1, itr
198  do idof = 1,ndof
199  x(idof)=zz(ndof*(i-1)+idof)
200  end do
201  do idof = 1, ndof
202  do jdof = 1, ndof
203  vv(ndof*(i-1)+idof) = vv(ndof*(i-1)+idof) + d(ndof2*(i-1)+ndof*(idof-1)+jdof)*x(jdof)
204  end do
205  end do
206 
207  js= indexl(i-1) + 1
208  je= indexl(i )
209  do j=js,je
210  in = iteml(j)
211  do idof = 1, ndof
212  do jdof = 1, ndof
213  vv(ndof*(in-1)+idof) = vv(ndof*(in-1)+idof) + al(ndof2*(j-1)+ndof*(jdof-1)+idof)*x(jdof)
214  end do
215  end do
216  enddo
217  js= indexu(i-1) + 1
218  je= indexu(i )
219  do j= js, je
220  in = itemu(j)
221  do idof = 1, ndof
222  do jdof = 1, ndof
223  vv(ndof*(in-1)+idof) = vv(ndof*(in-1)+idof) + au(ndof2*(j-1)+ndof*(jdof-1)+idof)*x(jdof)
224  end do
225  end do
226  enddo
227  enddo
228 
229  !{d}={v^t}{z_j}
230  do idof = 1, ndof
231  dtmp(idof)= sainvd(ndof2*(itr-1)+ndof*(idof-1)+idof)
232  end do
233 
234  do i= itr,np
235  do idof = 1,ndof
236  sainvd(ndof2*(i-1)+ndof*(idof-1)+idof) = vv(ndof*(i-1)+idof)
237  do jdof = 1, idof-1
238  sainvd(ndof2*(i-1)+ndof*(idof-1)+idof) = sainvd(ndof2*(i-1)+ndof*(idof-1)+idof) + &
239  & sainvd(ndof2*(i-1)+ndof*(jdof-1)+idof)*vv(ndof*(i-1)+jdof)
240  end do
241  end do
242  js= inumfi1l(i-1) + 1
243  je= inumfi1l(i )
244  do j= js, je
245  in = fi1l(j)
246  do idof = 1,ndof
247  x(idof)=vv(ndof*(in-1)+idof)
248  end do
249  do idof = 1, ndof
250  do jdof = 1, ndof
251  sainvd(ndof2*(i-1)+ndof*(idof-1)+idof) = sainvd(ndof2*(i-1)+ndof*(idof-1)+idof) + &
252  & sainvl(ndof2*(j-1)+ndof*(idof-1)+jdof)*x(jdof)
253  end do
254  end do
255  enddo
256  enddo
257 
258  !Update D
259  dd = 1.0d0/sainvd(ndof2*(itr-1)+ndof*(iitr-1)+iitr)
260  do idof=1,iitr-1
261  sainvd(ndof2*(itr-1)+ndof*(idof-1)+idof) = dtmp(idof)
262  end do
263  sainvd(ndof2*(itr-1)+ndof*(iitr-1)+iitr)=dd
264  do idof = iitr+1, ndof
265  sainvd(ndof2*(itr-1)+ndof*(idof-1)+idof) = sainvd(ndof2*(itr-1)+ndof*(idof-1)+idof)*dd
266  end do
267 
268  do i =itr+1,np
269  do idof = 1, ndof
270  sainvd(ndof2*(i-1)+ndof*(idof-1)+idof) = sainvd(ndof2*(i-1)+ndof*(idof-1)+idof)*dd
271  end do
272  enddo
273 
274  do idof = iitr, ndof
275  rifd(ndof2*(itr-1)+ndof*(idof-1)+iitr) = sainvd(ndof2*(itr-1)+ndof*(idof-1)+idof) !RIF UPDATE
276  end do
277 
278  js= inumfi1u(itr-1) + 1
279  je= inumfi1u(itr )
280  do j= js, je
281  do idof = 1, ndof
282  rifu(ndof2*(j-1)+ndof*(iitr-1)+idof) = sainvd(ndof2*(fi1u(j)-1)+ndof*(idof-1)+idof)
283  end do
284  enddo
285 
286  !Update Z
287  do k=iitr+1,ndof
288  dd = sainvd(ndof2*(itr-1)+ndof*(k-1)+k)
289  if(abs(dd) > filter)then
290  do jdof = 1, iitr
291  sainvd(ndof2*(itr-1)+ndof*(jdof-1)+k)= sainvd(ndof2*(itr-1)+ndof*(jdof-1)+k) - dd*zz(ndof*(itr-1)+jdof)
292  end do
293  js= inumfi1l(itr-1) + 1
294  je= inumfi1l(itr )
295  do j= js, je
296  in = fi1l(j)
297  do idof = 1, ndof
298  sainvl(ndof2*(j-1)+ndof*(k-1)+idof) = sainvl(ndof2*(j-1)+ndof*(k-1)+idof)-dd*zz(ndof*(in-1)+idof)
299  end do
300  enddo
301  endif
302  end do
303 
304  do i= itr +1,np
305  js= inumfi1l(i-1) + 1
306  je= inumfi1l(i )
307  do idof = 1, ndof
308  dd = sainvd(ndof2*(i-1)+ndof*(idof-1)+idof)
309  if(abs(dd) > filter)then
310  do j= js, je
311  in = fi1l(j)
312  if (in > itr) exit
313  do jdof=1,ndof
314  sainvl(ndof2*(j-1)+ndof*(idof-1)+jdof)=sainvl(ndof2*(j-1)+ndof*(idof-1)+jdof)-dd*zz(ndof*(in-1)+jdof)
315  end do
316  enddo
317  endif
318  end do
319  enddo
320  end do
321  enddo
322  deallocate(vv)
323  deallocate(zz)
324 
325  end subroutine hecmw_rif_nn
326 
327  subroutine hecmw_rif_make_u_nn(hecMAT)
328  implicit none
329  type (hecmwst_matrix) :: hecmat
330  integer(kind=kint) i,j,k,n,m,o,idof,jdof,ndof,ndof2
331  integer(kind=kint) is,ie,js,je
332  ndof=hecmat%NDOF
333  ndof2=ndof*ndof
334  n = 1
335  do i= 1, hecmat%NP
336  is=inumfi1l(i-1) + 1
337  ie=inumfi1l(i )
338  flag1:do k= is, ie
339  m = fi1l(k)
340  js=inumfi1u(m-1) + 1
341  je=inumfi1u(m )
342  do j= js,je
343  o = fi1u(j)
344  if (o == i)then
345  do idof = 1, ndof
346  do jdof = 1, ndof
347  rifl(ndof2*(n-1)+ndof*(idof-1)+jdof)=rifu(ndof2*(j-1)+ndof*(jdof-1)+idof)
348  end do
349  end do
350  n = n + 1
351  cycle flag1
352  endif
353  enddo
354  enddo flag1
355  enddo
356  end subroutine hecmw_rif_make_u_nn
357 
358  !C***
359  !C*** FORM_ILU1_nn
360  !C*** form ILU(1) matrix
361  subroutine form_ilu0_rif_nn(hecMAT)
362  implicit none
363  type(hecmwst_matrix) :: hecmat
364 
365  allocate (inumfi1l(0:hecmat%NP), inumfi1u(0:hecmat%NP))
366  allocate (fi1l(hecmat%NPL), fi1u(hecmat%NPU))
367 
368  inumfi1l = 0
369  inumfi1u = 0
370  fi1l = 0
371  fi1u = 0
372 
373  inumfi1l(:) = hecmat%indexL(:)
374  inumfi1u(:) = hecmat%indexU(:)
375  fi1l(:) = hecmat%itemL(:)
376  fi1u(:) = hecmat%itemU(:)
377 
378  npfiu = hecmat%NPU
379  npfil = hecmat%NPL
380 
381  end subroutine form_ilu0_rif_nn
382 
383  !C***
384  !C*** FORM_ILU1_nn
385  !C*** form ILU(1) matrix
386  subroutine form_ilu1_rif_nn(hecMAT)
387  implicit none
388  type(hecmwst_matrix) :: hecmat
389 
390  integer(kind=kint),allocatable :: iwsl(:), iwsu(:), iw1(:), iw2(:)
391  integer(kind=kint) :: nplf1,npuf1
392  integer(kind=kint) :: i,jj,kk,l,isk,iek,isj,iej
393  integer(kind=kint) :: icou,icou0,icouu,icouu1,icouu2,icouu3,icoul,icoul1,icoul2,icoul3
394  integer(kind=kint) :: j,k,isl,isu
395  !C
396  !C +--------------+
397  !C | find fill-in |
398  !C +--------------+
399  !C===
400 
401  !C
402  !C-- count fill-in
403  allocate (iw1(hecmat%NP) , iw2(hecmat%NP))
404  allocate (inumfi1l(0:hecmat%NP), inumfi1u(0:hecmat%NP))
405 
406  inumfi1l= 0
407  inumfi1u= 0
408 
409  nplf1= 0
410  npuf1= 0
411  do i= 2, hecmat%NP
412  icou= 0
413  iw1= 0
414  iw1(i)= 1
415  do l= indexl(i-1)+1, indexl(i)
416  iw1(iteml(l))= 1
417  enddo
418  do l= indexu(i-1)+1, indexu(i)
419  iw1(itemu(l))= 1
420  enddo
421 
422  isk= indexl(i-1) + 1
423  iek= indexl(i)
424  do k= isk, iek
425  kk= iteml(k)
426  isj= indexu(kk-1) + 1
427  iej= indexu(kk )
428  do j= isj, iej
429  jj= itemu(j)
430  if (iw1(jj).eq.0 .and. jj.lt.i) then
431  inumfi1l(i)= inumfi1l(i)+1
432  iw1(jj)= 1
433  endif
434  if (iw1(jj).eq.0 .and. jj.gt.i) then
435  inumfi1u(i)= inumfi1u(i)+1
436  iw1(jj)= 1
437  endif
438  enddo
439  enddo
440  nplf1= nplf1 + inumfi1l(i)
441  npuf1= npuf1 + inumfi1u(i)
442  enddo
443 
444  !C
445  !C-- specify fill-in
446  allocate (iwsl(0:hecmat%NP), iwsu(0:hecmat%NP))
447  allocate (fi1l(hecmat%NPL+nplf1), fi1u(hecmat%NPU+npuf1))
448 
449  npfiu = hecmat%NPU+npuf1
450  npfil = hecmat%NPL+nplf1
451 
452  fi1l= 0
453  fi1u= 0
454 
455  iwsl= 0
456  iwsu= 0
457  do i= 1, hecmat%NP
458  iwsl(i)= indexl(i)-indexl(i-1) + inumfi1l(i) + iwsl(i-1)
459  iwsu(i)= indexu(i)-indexu(i-1) + inumfi1u(i) + iwsu(i-1)
460  enddo
461 
462  do i= 2, hecmat%NP
463  icoul= 0
464  icouu= 0
465  inumfi1l(i)= inumfi1l(i-1) + inumfi1l(i)
466  inumfi1u(i)= inumfi1u(i-1) + inumfi1u(i)
467  icou= 0
468  iw1= 0
469  iw1(i)= 1
470  do l= indexl(i-1)+1, indexl(i)
471  iw1(iteml(l))= 1
472  enddo
473  do l= indexu(i-1)+1, indexu(i)
474  iw1(itemu(l))= 1
475  enddo
476 
477  isk= indexl(i-1) + 1
478  iek= indexl(i)
479  do k= isk, iek
480  kk= iteml(k)
481  isj= indexu(kk-1) + 1
482  iej= indexu(kk )
483  do j= isj, iej
484  jj= itemu(j)
485  if (iw1(jj).eq.0 .and. jj.lt.i) then
486  icoul = icoul + 1
487  fi1l(icoul+iwsl(i-1)+indexl(i)-indexl(i-1))= jj
488  iw1(jj) = 1
489  endif
490  if (iw1(jj).eq.0 .and. jj.gt.i) then
491  icouu = icouu + 1
492  fi1u(icouu+iwsu(i-1)+indexu(i)-indexu(i-1))= jj
493  iw1(jj) = 1
494  endif
495  enddo
496  enddo
497  enddo
498 
499  isl = 0
500  isu = 0
501  do i= 1, hecmat%NP
502  icoul1= indexl(i) - indexl(i-1)
503  icoul2= inumfi1l(i) - inumfi1l(i-1)
504  icoul3= icoul1 + icoul2
505  icouu1= indexu(i) - indexu(i-1)
506  icouu2= inumfi1u(i) - inumfi1u(i-1)
507  icouu3= icouu1 + icouu2
508  !C
509  !C-- LOWER part
510  icou0= 0
511  do k= indexl(i-1)+1, indexl(i)
512  icou0 = icou0 + 1
513  iw1(icou0)= iteml(k)
514  enddo
515 
516  do k= inumfi1l(i-1)+1, inumfi1l(i)
517  icou0 = icou0 + 1
518  iw1(icou0)= fi1l(icou0+iwsl(i-1))
519  enddo
520 
521  do k= 1, icoul3
522  iw2(k)= k
523  enddo
524  call rif_sort_nn (iw1, iw2, icoul3, hecmat%NP)
525 
526  do k= 1, icoul3
527  fi1l(k+isl)= iw1(k)
528  enddo
529  !C
530  !C-- UPPER part
531  icou0= 0
532  do k= indexu(i-1)+1, indexu(i)
533  icou0 = icou0 + 1
534  iw1(icou0)= itemu(k)
535  enddo
536 
537  do k= inumfi1u(i-1)+1, inumfi1u(i)
538  icou0 = icou0 + 1
539  iw1(icou0)= fi1u(icou0+iwsu(i-1))
540  enddo
541 
542  do k= 1, icouu3
543  iw2(k)= k
544  enddo
545  call rif_sort_nn (iw1, iw2, icouu3, hecmat%NP)
546 
547  do k= 1, icouu3
548  fi1u(k+isu)= iw1(k)
549  enddo
550 
551  isl= isl + icoul3
552  isu= isu + icouu3
553  enddo
554 
555  !C===
556  do i= 1, hecmat%NP
557  inumfi1l(i)= iwsl(i)
558  inumfi1u(i)= iwsu(i)
559  enddo
560 
561  deallocate (iw1, iw2)
562  deallocate (iwsl, iwsu)
563  !C===
564  end subroutine form_ilu1_rif_nn
565 
566  !C
567  !C***
568  !C*** fill_in_S33_SORT
569  !C***
570  !C
571  subroutine rif_sort_nn(STEM, INUM, N, NP)
572  use hecmw_util
573  implicit none
574  integer(kind=kint) :: n, np
575  integer(kind=kint), dimension(NP) :: stem
576  integer(kind=kint), dimension(NP) :: inum
577  integer(kind=kint), dimension(:), allocatable :: istack
578  integer(kind=kint) :: m,nstack,jstack,l,ir,ip,i,j,k,ss,ii,temp,it
579 
580  allocate (istack(-np:+np))
581 
582  m = 100
583  nstack= np
584 
585  jstack= 0
586  l = 1
587  ir = n
588 
589  ip= 0
590  1 continue
591  ip= ip + 1
592 
593  if (ir-l.lt.m) then
594  do j= l+1, ir
595  ss= stem(j)
596  ii= inum(j)
597 
598  do i= j-1,1,-1
599  if (stem(i).le.ss) goto 2
600  stem(i+1)= stem(i)
601  inum(i+1)= inum(i)
602  end do
603  i= 0
604 
605  2 continue
606  stem(i+1)= ss
607  inum(i+1)= ii
608  end do
609 
610  if (jstack.eq.0) then
611  deallocate (istack)
612  return
613  endif
614 
615  ir = istack(jstack)
616  l = istack(jstack-1)
617  jstack= jstack - 2
618  else
619 
620  k= (l+ir) / 2
621  temp = stem(k)
622  stem(k) = stem(l+1)
623  stem(l+1)= temp
624 
625  it = inum(k)
626  inum(k) = inum(l+1)
627  inum(l+1)= it
628 
629  if (stem(l+1).gt.stem(ir)) then
630  temp = stem(l+1)
631  stem(l+1)= stem(ir)
632  stem(ir )= temp
633  it = inum(l+1)
634  inum(l+1)= inum(ir)
635  inum(ir )= it
636  endif
637 
638  if (stem(l).gt.stem(ir)) then
639  temp = stem(l)
640  stem(l )= stem(ir)
641  stem(ir)= temp
642  it = inum(l)
643  inum(l )= inum(ir)
644  inum(ir)= it
645  endif
646 
647  if (stem(l+1).gt.stem(l)) then
648  temp = stem(l+1)
649  stem(l+1)= stem(l)
650  stem(l )= temp
651  it = inum(l+1)
652  inum(l+1)= inum(l)
653  inum(l )= it
654  endif
655 
656  i= l + 1
657  j= ir
658 
659  ss= stem(l)
660  ii= inum(l)
661 
662  3 continue
663  i= i + 1
664  if (stem(i).lt.ss) goto 3
665 
666  4 continue
667  j= j - 1
668  if (stem(j).gt.ss) goto 4
669 
670  if (j.lt.i) goto 5
671 
672  temp = stem(i)
673  stem(i)= stem(j)
674  stem(j)= temp
675 
676  it = inum(i)
677  inum(i)= inum(j)
678  inum(j)= it
679 
680  goto 3
681 
682  5 continue
683 
684  stem(l)= stem(j)
685  stem(j)= ss
686  inum(l)= inum(j)
687  inum(j)= ii
688 
689  jstack= jstack + 2
690 
691  if (jstack.gt.nstack) then
692  write (*,*) 'NSTACK overflow'
693  stop
694  endif
695 
696  if (ir-i+1.ge.j-1) then
697  istack(jstack )= ir
698  istack(jstack-1)= i
699  ir= j-1
700  else
701  istack(jstack )= j-1
702  istack(jstack-1)= l
703  l= i
704  endif
705 
706  endif
707 
708  goto 1
709 
710  end subroutine rif_sort_nn
711 
712  subroutine hecmw_precond_rif_nn_clear()
713  implicit none
714 
715  if (associated(sainvd)) deallocate(sainvd)
716  if (associated(sainvl)) deallocate(sainvl)
717  if (associated(rifd)) deallocate(rifd)
718  if (associated(rifu)) deallocate(rifu)
719  if (associated(rifl)) deallocate(rifl)
720  if (associated(inumfi1l)) deallocate(inumfi1l)
721  if (associated(inumfi1u)) deallocate(inumfi1u)
722  if (associated(fi1l)) deallocate(fi1l)
723  if (associated(fi1u)) deallocate(fi1u)
724  nullify(inumfi1l)
725  nullify(inumfi1u)
726  nullify(fi1l)
727  nullify(fi1u)
728  nullify(d)
729  nullify(al)
730  nullify(au)
731  nullify(indexl)
732  nullify(indexu)
733  nullify(iteml)
734  nullify(itemu)
735 
736  end subroutine hecmw_precond_rif_nn_clear
737 end module hecmw_precond_rif_nn
hecmw_precond_rif_nn::hecmw_precond_rif_nn_apply
subroutine, public hecmw_precond_rif_nn_apply(ZP, NDOF)
Definition: hecmw_precond_RIF_nn.f90:91
hecmw_matrix_misc::hecmw_mat_get_precond
integer(kind=kint) function, public hecmw_mat_get_precond(hecMAT)
Definition: hecmw_matrix_misc.f90:321
hecmw_precond_rif_nn
Definition: hecmw_precond_RIF_nn.f90:6
hecmw_util
I/O and Utility.
Definition: hecmw_util_f.F90:7
m_hecmw_comm_f
Definition: hecmw_comm_f.F90:6
hecmw_util::kreal
integer(kind=4), parameter kreal
Definition: hecmw_util_f.F90:16
hecmw_precond_rif_nn::hecmw_precond_rif_nn_clear
subroutine, public hecmw_precond_rif_nn_clear()
Definition: hecmw_precond_RIF_nn.f90:713
hecmw_matrix_misc
Definition: hecmw_matrix_misc.f90:6
hecmw_precond_rif_nn::hecmw_precond_rif_nn_setup
subroutine, public hecmw_precond_rif_nn_setup(hecMAT)
Definition: hecmw_precond_RIF_nn.f90:46
hecmw_util::hecmwst_matrix
Definition: hecmw_util_f.F90:444