18 integer(4),
parameter :: krealp = 8
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()
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()
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()
50 integer(kind=kint ) :: precond
52 real(kind=krealp) :: filter
60 indexl => hecmat%indexL
61 indexu => hecmat%indexU
65 if (precond.eq.21)
call form_ilu1_rif_33(hecmat)
67 allocate (sainvd(9*hecmat%NP))
68 allocate (sainvl(9*npfiu))
69 allocate (rifd(9*hecmat%NP))
70 allocate (rifu(9*npfiu))
76 filter= hecmat%Rarray(5)
78 write(*,
"(a,F15.8)")
"### RIF FILTER :",filter
80 call hecmw_rif_33(hecmat)
82 allocate (rifl(9*npfiu))
85 call hecmw_rif_make_u_33(hecmat)
91 real(kind=
kreal),
intent(inout) :: zp(:)
92 integer(kind=kint) :: in, i, j, isl, iel
93 real(kind=
kreal) :: sw1, sw2, sw3, x1, x2, x3
108 sw1= sw1 - rifl(9*j-8)*x1 - rifl(9*j-7)*x2 - rifl(9*j-6)*x3
109 sw2= sw2 - rifl(9*j-5)*x1 - rifl(9*j-4)*x2 - rifl(9*j-3)*x3
110 sw3= sw3 - rifl(9*j-2)*x1 - rifl(9*j-1)*x2 - rifl(9*j )*x3
115 x2= x2 - rifd(9*i-5)*x1
116 x3= x3 - rifd(9*i-2)*x1 - rifd(9*i-1)*x2
123 zp(3*i-2)= zp(3*i-2)*rifd(9*i-8)
124 zp(3*i-1)= zp(3*i-1)*rifd(9*i-4)
125 zp(3*i )= zp(3*i )*rifd(9*i )
134 isl= inumfi1u(i-1) + 1
141 sw1= sw1 - rifu(9*j-8)*x1 - rifu(9*j-7)*x2 - rifu(9*j-6)*x3
142 sw2= sw2 - rifu(9*j-5)*x1 - rifu(9*j-4)*x2 - rifu(9*j-3)*x3
143 sw3= sw3 - rifu(9*j-2)*x1 - rifu(9*j-1)*x2 - rifu(9*j )*x3
149 x2= x2 - rifd(9*i-1)*x3
150 x1= x1 - rifd(9*i-2)*x3 - rifd(9*i-5)*x2
161 subroutine hecmw_rif_33(hecMAT)
164 integer(kind=kint) :: i, j, js, je, in, itr, np
165 real(kind=krealp) :: x1, x2, x3, dd, dd1, dd2, dd3, dtemp(3)
166 real(kind=krealp) :: filter
167 real(kind=krealp),
allocatable :: zz(:), vv(:)
169 filter= hecmat%Rarray(5)
184 zz(3*itr-2)= sainvd(9*itr-8)
185 zz(3*itr-1)= sainvd(9*itr-5)
186 zz(3*itr )= sainvd(9*itr-2)
190 js= inumfi1l(itr-1) + 1
194 zz(3*in-2)= sainvl(9*j-8)
195 zz(3*in-1)= sainvl(9*j-7)
196 zz(3*in )= sainvl(9*j-6)
203 vv(3*i-2) = vv(3*i-2) + d(9*i-8)*x1 + d(9*i-7)*x2 + d(9*i-6)*x3
204 vv(3*i-1) = vv(3*i-1) + d(9*i-5)*x1 + d(9*i-4)*x2 + d(9*i-3)*x3
205 vv(3*i ) = vv(3*i ) + d(9*i-2)*x1 + d(9*i-1)*x2 + d(9*i )*x3
211 vv(3*in-2)= vv(3*in-2) + al(9*j-8)*x1 + al(9*j-5)*x2 + al(9*j-2)*x3
212 vv(3*in-1)= vv(3*in-1) + al(9*j-7)*x1 + al(9*j-4)*x2 + al(9*j-1)*x3
213 vv(3*in )= vv(3*in ) + al(9*j-6)*x1 + al(9*j-3)*x2 + al(9*j )*x3
220 vv(3*in-2)= vv(3*in-2) + au(9*j-8)*x1 + au(9*j-5)*x2 + au(9*j-2)*x3
221 vv(3*in-1)= vv(3*in-1) + au(9*j-7)*x1 + au(9*j-4)*x2 + au(9*j-1)*x3
222 vv(3*in )= vv(3*in ) + au(9*j-6)*x1 + au(9*j-3)*x2 + au(9*j )*x3
228 sainvd(9*i-8) = vv(3*i-2)
229 sainvd(9*i-4) = vv(3*i-2)*sainvd(9*i-7) + vv(3*i-1)
230 sainvd(9*i ) = vv(3*i-2)*sainvd(9*i-6) + vv(3*i-1)*sainvd(9*i-3) + vv(3*i)
231 js= inumfi1l(i-1) + 1
238 sainvd(9*i-8)= sainvd(9*i-8) + x1*sainvl(9*j-8) + x2*sainvl(9*j-7) + x3*sainvl(9*j-6)
239 sainvd(9*i-4)= sainvd(9*i-4) + x1*sainvl(9*j-5) + x2*sainvl(9*j-4) + x3*sainvl(9*j-3)
240 sainvd(9*i )= sainvd(9*i ) + x1*sainvl(9*j-2) + x2*sainvl(9*j-1) + x3*sainvl(9*j )
245 dd = 1.0d0/sainvd(9*itr-8)
247 sainvd(9*itr-4) =sainvd(9*itr-4)*dd
248 sainvd(9*itr ) =sainvd(9*itr )*dd
251 sainvd(9*i-8) = sainvd(9*i-8)*dd
252 sainvd(9*i-4) = sainvd(9*i-4)*dd
253 sainvd(9*i ) = sainvd(9*i )*dd
257 rifd(9*itr-5) = sainvd(9*itr-4)
258 rifd(9*itr-2) = sainvd(9*itr )
260 js= inumfi1u(itr-1) + 1
263 rifu(9*j-8) = sainvd(9*fi1u(j)-8)
264 rifu(9*j-7) = sainvd(9*fi1u(j)-4)
265 rifu(9*j-6) = sainvd(9*fi1u(j) )
271 if(abs(dd2) > filter)
then
272 sainvd(9*itr-7)= sainvd(9*itr-7) - dd2*zz(3*itr-2)
273 js= inumfi1l(itr-1) + 1
277 sainvl(9*j-5) = sainvl(9*j-5)-dd2*zz(3*in-2)
278 sainvl(9*j-4) = sainvl(9*j-4)-dd2*zz(3*in-1)
279 sainvl(9*j-3) = sainvl(9*j-3)-dd2*zz(3*in )
284 if(abs(dd3) > filter)
then
285 sainvd(9*itr-6)= sainvd(9*itr-6) - dd3*zz(3*itr-2)
286 js= inumfi1l(itr-1) + 1
290 sainvl(9*j-2) = sainvl(9*j-2)-dd3*zz(3*in-2)
291 sainvl(9*j-1) = sainvl(9*j-1)-dd3*zz(3*in-1)
292 sainvl(9*j ) = sainvl(9*j )-dd3*zz(3*in )
297 js= inumfi1l(i-1) + 1
300 if(abs(dd1) > filter)
then
304 sainvl(9*j-8) = sainvl(9*j-8)-dd1*zz(3*in-2)
305 sainvl(9*j-7) = sainvl(9*j-7)-dd1*zz(3*in-1)
306 sainvl(9*j-6) = sainvl(9*j-6)-dd1*zz(3*in )
310 if(abs(dd2) > filter)
then
314 sainvl(9*j-5) = sainvl(9*j-5)-dd2*zz(3*in-2)
315 sainvl(9*j-4) = sainvl(9*j-4)-dd2*zz(3*in-1)
316 sainvl(9*j-3) = sainvl(9*j-3)-dd2*zz(3*in )
320 if(abs(dd3) > filter)
then
324 sainvl(9*j-2) = sainvl(9*j-2)-dd3*zz(3*in-2)
325 sainvl(9*j-1) = sainvl(9*j-1)-dd3*zz(3*in-1)
326 sainvl(9*j ) = sainvl(9*j )-dd3*zz(3*in )
338 zz(3*itr-2)= sainvd(9*itr-7)
339 zz(3*itr-1)= sainvd(9*itr-4)
340 zz(3*itr )= sainvd(9*itr-1)
344 js= inumfi1l(itr-1) + 1
348 zz(3*in-2)= sainvl(9*j-5)
349 zz(3*in-1)= sainvl(9*j-4)
350 zz(3*in )= sainvl(9*j-3)
357 vv(3*i-2) = vv(3*i-2) + d(9*i-8)*x1 + d(9*i-7)*x2 + d(9*i-6)*x3
358 vv(3*i-1) = vv(3*i-1) + d(9*i-5)*x1 + d(9*i-4)*x2 + d(9*i-3)*x3
359 vv(3*i ) = vv(3*i ) + d(9*i-2)*x1 + d(9*i-1)*x2 + d(9*i )*x3
365 vv(3*in-2)= vv(3*in-2) + al(9*j-8)*x1 + al(9*j-5)*x2 + al(9*j-2)*x3
366 vv(3*in-1)= vv(3*in-1) + al(9*j-7)*x1 + al(9*j-4)*x2 + al(9*j-1)*x3
367 vv(3*in )= vv(3*in ) + al(9*j-6)*x1 + al(9*j-3)*x2 + al(9*j )*x3
374 vv(3*in-2)= vv(3*in-2) + au(9*j-8)*x1 + au(9*j-5)*x2 + au(9*j-2)*x3
375 vv(3*in-1)= vv(3*in-1) + au(9*j-7)*x1 + au(9*j-4)*x2 + au(9*j-1)*x3
376 vv(3*in )= vv(3*in ) + au(9*j-6)*x1 + au(9*j-3)*x2 + au(9*j )*x3
381 dtemp(1) = sainvd(9*itr-8)
384 sainvd(9*i-8) = vv(3*i-2)
385 sainvd(9*i-4) = vv(3*i-2)*sainvd(9*i-7) + vv(3*i-1)
386 sainvd(9*i ) = vv(3*i-2)*sainvd(9*i-6) + vv(3*i-1)*sainvd(9*i-3) + vv(3*i)
387 js= inumfi1l(i-1) + 1
394 sainvd(9*i-8)= sainvd(9*i-8) + x1*sainvl(9*j-8) + x2*sainvl(9*j-7) + x3*sainvl(9*j-6)
395 sainvd(9*i-4)= sainvd(9*i-4) + x1*sainvl(9*j-5) + x2*sainvl(9*j-4) + x3*sainvl(9*j-3)
396 sainvd(9*i )= sainvd(9*i ) + x1*sainvl(9*j-2) + x2*sainvl(9*j-1) + x3*sainvl(9*j )
401 dd = 1.0d0/sainvd(9*itr-4)
403 sainvd(9*itr-8) = dtemp(1)
404 sainvd(9*itr ) =sainvd(9*itr )*dd
407 sainvd(9*i-8) = sainvd(9*i-8)*dd
408 sainvd(9*i-4) = sainvd(9*i-4)*dd
409 sainvd(9*i ) = sainvd(9*i )*dd
413 rifd(9*itr-1) = sainvd(9*itr )
415 js= inumfi1u(itr-1) + 1
418 rifu(9*j-5) = sainvd(9*fi1u(j)-8)
419 rifu(9*j-4) = sainvd(9*fi1u(j)-4)
420 rifu(9*j-3) = sainvd(9*fi1u(j) )
426 if(abs(dd3) > filter)
then
427 sainvd(9*itr-6)= sainvd(9*itr-6) - dd3*zz(3*itr-2)
428 sainvd(9*itr-3)= sainvd(9*itr-3) - dd3*zz(3*itr-1)
430 js= inumfi1l(itr-1) + 1
434 sainvl(9*j-2) = sainvl(9*j-2)-dd3*zz(3*in-2)
435 sainvl(9*j-1) = sainvl(9*j-1)-dd3*zz(3*in-1)
436 sainvl(9*j ) = sainvl(9*j )-dd3*zz(3*in )
441 js= inumfi1l(i-1) + 1
444 if(abs(dd1) > filter)
then
448 sainvl(9*j-8) = sainvl(9*j-8)-dd1*zz(3*in-2)
449 sainvl(9*j-7) = sainvl(9*j-7)-dd1*zz(3*in-1)
450 sainvl(9*j-6) = sainvl(9*j-6)-dd1*zz(3*in )
454 if(abs(dd2) > filter)
then
458 sainvl(9*j-5) = sainvl(9*j-5)-dd2*zz(3*in-2)
459 sainvl(9*j-4) = sainvl(9*j-4)-dd2*zz(3*in-1)
460 sainvl(9*j-3) = sainvl(9*j-3)-dd2*zz(3*in )
464 if(abs(dd3) > filter)
then
468 sainvl(9*j-2) = sainvl(9*j-2)-dd3*zz(3*in-2)
469 sainvl(9*j-1) = sainvl(9*j-1)-dd3*zz(3*in-1)
470 sainvl(9*j ) = sainvl(9*j )-dd3*zz(3*in )
483 zz(3*itr-2)= sainvd(9*itr-6)
484 zz(3*itr-1)= sainvd(9*itr-3)
485 zz(3*itr )= sainvd(9*itr )
489 js= inumfi1l(itr-1) + 1
493 zz(3*in-2)= sainvl(9*j-2)
494 zz(3*in-1)= sainvl(9*j-1)
495 zz(3*in )= sainvl(9*j )
502 vv(3*i-2) = vv(3*i-2) + d(9*i-8)*x1 + d(9*i-7)*x2 + d(9*i-6)*x3
503 vv(3*i-1) = vv(3*i-1) + d(9*i-5)*x1 + d(9*i-4)*x2 + d(9*i-3)*x3
504 vv(3*i ) = vv(3*i ) + d(9*i-2)*x1 + d(9*i-1)*x2 + d(9*i )*x3
510 vv(3*in-2)= vv(3*in-2) + al(9*j-8)*x1 + al(9*j-5)*x2 + al(9*j-2)*x3
511 vv(3*in-1)= vv(3*in-1) + al(9*j-7)*x1 + al(9*j-4)*x2 + al(9*j-1)*x3
512 vv(3*in )= vv(3*in ) + al(9*j-6)*x1 + al(9*j-3)*x2 + al(9*j )*x3
519 vv(3*in-2)= vv(3*in-2) + au(9*j-8)*x1 + au(9*j-5)*x2 + au(9*j-2)*x3
520 vv(3*in-1)= vv(3*in-1) + au(9*j-7)*x1 + au(9*j-4)*x2 + au(9*j-1)*x3
521 vv(3*in )= vv(3*in ) + au(9*j-6)*x1 + au(9*j-3)*x2 + au(9*j )*x3
527 dtemp(1) = sainvd(9*itr-8)
528 dtemp(2) = sainvd(9*itr-4)
531 sainvd(9*i-8) = vv(3*i-2)
532 sainvd(9*i-4) = vv(3*i-2)*sainvd(9*i-7) + vv(3*i-1)
533 sainvd(9*i ) = vv(3*i-2)*sainvd(9*i-6) + vv(3*i-1)*sainvd(9*i-3) + vv(3*i)
534 js= inumfi1l(i-1) + 1
541 sainvd(9*i-8)= sainvd(9*i-8) + x1*sainvl(9*j-8) + x2*sainvl(9*j-7) + x3*sainvl(9*j-6)
542 sainvd(9*i-4)= sainvd(9*i-4) + x1*sainvl(9*j-5) + x2*sainvl(9*j-4) + x3*sainvl(9*j-3)
543 sainvd(9*i )= sainvd(9*i ) + x1*sainvl(9*j-2) + x2*sainvl(9*j-1) + x3*sainvl(9*j )
548 dd = 1.0d0/sainvd(9*itr )
550 sainvd(9*itr-8) = dtemp(1)
551 sainvd(9*itr-4) = dtemp(2)
554 sainvd(9*i-8) = sainvd(9*i-8)*dd
555 sainvd(9*i-4) = sainvd(9*i-4)*dd
556 sainvd(9*i ) = sainvd(9*i )*dd
561 js= inumfi1u(itr-1) + 1
564 rifu(9*j-2) = sainvd(9*fi1u(j)-8)
565 rifu(9*j-1) = sainvd(9*fi1u(j)-4)
566 rifu(9*j ) = sainvd(9*fi1u(j) )
572 js= inumfi1l(i-1) + 1
575 if(abs(dd1) > filter)
then
579 sainvl(9*j-8) = sainvl(9*j-8)-dd1*zz(3*in-2)
580 sainvl(9*j-7) = sainvl(9*j-7)-dd1*zz(3*in-1)
581 sainvl(9*j-6) = sainvl(9*j-6)-dd1*zz(3*in )
585 if(abs(dd2) > filter)
then
589 sainvl(9*j-5) = sainvl(9*j-5)-dd2*zz(3*in-2)
590 sainvl(9*j-4) = sainvl(9*j-4)-dd2*zz(3*in-1)
591 sainvl(9*j-3) = sainvl(9*j-3)-dd2*zz(3*in )
595 if(abs(dd3) > filter)
then
599 sainvl(9*j-2) = sainvl(9*j-2)-dd3*zz(3*in-2)
600 sainvl(9*j-1) = sainvl(9*j-1)-dd3*zz(3*in-1)
601 sainvl(9*j ) = sainvl(9*j )-dd3*zz(3*in )
609 end subroutine hecmw_rif_33
611 subroutine hecmw_rif_make_u_33(hecMAT)
614 integer(kind=kint) i,j,k,n,m,o
615 integer(kind=kint) is,ie,js,je
628 rifl(9*n-8)=rifu(9*j-8)
629 rifl(9*n-7)=rifu(9*j-5)
630 rifl(9*n-6)=rifu(9*j-2)
631 rifl(9*n-5)=rifu(9*j-7)
632 rifl(9*n-4)=rifu(9*j-4)
633 rifl(9*n-3)=rifu(9*j-1)
634 rifl(9*n-2)=rifu(9*j-6)
635 rifl(9*n-1)=rifu(9*j-3)
636 rifl(9*n )=rifu(9*j )
643 end subroutine hecmw_rif_make_u_33
648 subroutine form_ilu0_rif_33(hecMAT)
652 allocate (inumfi1l(0:hecmat%NP), inumfi1u(0:hecmat%NP))
653 allocate (fi1l(hecmat%NPL), fi1u(hecmat%NPU))
660 inumfi1l(:) = hecmat%indexL(:)
661 inumfi1u(:) = hecmat%indexU(:)
662 fi1l(:) = hecmat%itemL(:)
663 fi1u(:) = hecmat%itemU(:)
668 end subroutine form_ilu0_rif_33
673 subroutine form_ilu1_rif_33(hecMAT)
677 integer(kind=kint),
allocatable :: iwsl(:), iwsu(:), iw1(:), iw2(:)
678 integer(kind=kint) :: nplf1,npuf1
679 integer(kind=kint) :: i,jj,kk,l,isk,iek,isj,iej
680 integer(kind=kint) :: icou,icou0,icouu,icouu1,icouu2,icouu3,icoul,icoul1,icoul2,icoul3
681 integer(kind=kint) :: j,k,isl,isu
690 allocate (iw1(hecmat%NP) , iw2(hecmat%NP))
691 allocate (inumfi1l(0:hecmat%NP), inumfi1u(0:hecmat%NP))
702 do l= indexl(i-1)+1, indexl(i)
705 do l= indexu(i-1)+1, indexu(i)
713 isj= indexu(kk-1) + 1
717 if (iw1(jj).eq.0 .and. jj.lt.i)
then
718 inumfi1l(i)= inumfi1l(i)+1
721 if (iw1(jj).eq.0 .and. jj.gt.i)
then
722 inumfi1u(i)= inumfi1u(i)+1
727 nplf1= nplf1 + inumfi1l(i)
728 npuf1= npuf1 + inumfi1u(i)
733 allocate (iwsl(0:hecmat%NP), iwsu(0:hecmat%NP))
734 allocate (fi1l(hecmat%NPL+nplf1), fi1u(hecmat%NPU+npuf1))
736 npfiu = hecmat%NPU+npuf1
737 npfil = hecmat%NPL+nplf1
745 iwsl(i)= indexl(i)-indexl(i-1) + inumfi1l(i) + iwsl(i-1)
746 iwsu(i)= indexu(i)-indexu(i-1) + inumfi1u(i) + iwsu(i-1)
752 inumfi1l(i)= inumfi1l(i-1) + inumfi1l(i)
753 inumfi1u(i)= inumfi1u(i-1) + inumfi1u(i)
757 do l= indexl(i-1)+1, indexl(i)
760 do l= indexu(i-1)+1, indexu(i)
768 isj= indexu(kk-1) + 1
772 if (iw1(jj).eq.0 .and. jj.lt.i)
then
774 fi1l(icoul+iwsl(i-1)+indexl(i)-indexl(i-1))= jj
777 if (iw1(jj).eq.0 .and. jj.gt.i)
then
779 fi1u(icouu+iwsu(i-1)+indexu(i)-indexu(i-1))= jj
789 icoul1= indexl(i) - indexl(i-1)
790 icoul2= inumfi1l(i) - inumfi1l(i-1)
791 icoul3= icoul1 + icoul2
792 icouu1= indexu(i) - indexu(i-1)
793 icouu2= inumfi1u(i) - inumfi1u(i-1)
794 icouu3= icouu1 + icouu2
798 do k= indexl(i-1)+1, indexl(i)
803 do k= inumfi1l(i-1)+1, inumfi1l(i)
805 iw1(icou0)= fi1l(icou0+iwsl(i-1))
811 call rif_sort_33 (iw1, iw2, icoul3, hecmat%NP)
819 do k= indexu(i-1)+1, indexu(i)
824 do k= inumfi1u(i-1)+1, inumfi1u(i)
826 iw1(icou0)= fi1u(icou0+iwsu(i-1))
832 call rif_sort_33 (iw1, iw2, icouu3, hecmat%NP)
848 deallocate (iw1, iw2)
849 deallocate (iwsl, iwsu)
851 end subroutine form_ilu1_rif_33
858 subroutine rif_sort_33(STEM, INUM, N, NP)
861 integer(kind=kint) :: n, np
862 integer(kind=kint),
dimension(NP) :: stem
863 integer(kind=kint),
dimension(NP) :: inum
864 integer(kind=kint),
dimension(:),
allocatable :: istack
865 integer(kind=kint) :: m,nstack,jstack,l,ir,ip,i,j,k,ss,ii,temp,it
867 allocate (istack(-np:+np))
886 if (stem(i).le.ss)
goto 2
897 if (jstack.eq.0)
then
916 if (stem(l+1).gt.stem(ir))
then
925 if (stem(l).gt.stem(ir))
then
934 if (stem(l+1).gt.stem(l))
then
951 if (stem(i).lt.ss)
goto 3
955 if (stem(j).gt.ss)
goto 4
978 if (jstack.gt.nstack)
then
979 write (*,*)
'NSTACK overflow'
983 if (ir-i+1.ge.j-1)
then
997 end subroutine rif_sort_33
1002 if (
associated(sainvd))
deallocate(sainvd)
1003 if (
associated(sainvl))
deallocate(sainvl)
1004 if (
associated(rifd))
deallocate(rifd)
1005 if (
associated(rifu))
deallocate(rifu)
1006 if (
associated(rifl))
deallocate(rifl)
1007 if (
associated(inumfi1l))
deallocate(inumfi1l)
1008 if (
associated(inumfi1u))
deallocate(inumfi1u)
1009 if (
associated(fi1l))
deallocate(fi1l)
1010 if (
associated(fi1u))
deallocate(fi1u)
integer(kind=kint) function, public hecmw_mat_get_precond(hecMAT)
subroutine, public hecmw_precond_rif_33_apply(ZP)
subroutine, public hecmw_precond_rif_33_setup(hecMAT)
subroutine, public hecmw_precond_rif_33_clear()
integer(kind=4), parameter kreal