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