FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
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
41contains
42
43 !C***
44 !C*** hecmw_precond_nn_sainv_setup
45 !C***
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
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
728end module hecmw_precond_sainv_nn
integer(kind=kint) function, public hecmw_mat_get_precond(hecmat)
subroutine, public hecmw_precond_nn_sainv_apply(r, zp)
subroutine, public hecmw_precond_nn_sainv_clear()
subroutine, public hecmw_precond_nn_sainv_setup(hecmat)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal