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