FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
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
41contains
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
738end 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_setup(hecmat)
subroutine, public hecmw_precond_rif_nn_apply(zp, ndof)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal