!
!  coefh_mg_zero.f90
!
       subroutine coefh_mg_zero(p,index1r,index2r,index1,index2,& 
                         acoef1,bcoef1,acoef2,bcoef2,&
                         ap1,ae1,aw1,an1,as1,at1,ab1,b1,&
                         ap2,ae2,aw2,an2,as2,at2,ab2,b2,&
                         icvs1,jcvs1,kcvs1,icvs2,jcvs2,kcvs2,ncs1,ncs2,&
                         rz,iz,dz,cz,rzp,izp,czp,ierr)
!
       implicit none
       real :: rzp(*),rz(*)
       integer :: izp(*),iz(*)
       character (len=*) :: czp(*),cz(*)
       real (selected_int_kind (5)) :: dz(*)
       integer :: ierr
!
       integer :: p,icvs1,jcvs1,kcvs1,icvs2,jcvs2,kcvs2,ncs1,ncs2
       integer :: index1(ncs1,4),index2(ncs2,4),index1r(6),index2r(6)
       real :: acoef1(ncs1),bcoef1(ncs1),acoef2(ncs2),bcoef2(ncs2)
!       
       real :: ap1(icvs1,jcvs1,kcvs1),ae1(icvs1,jcvs1,kcvs1),aw1(icvs1,jcvs1,kcvs1), &
               an1(icvs1,jcvs1,kcvs1),as1(icvs1,jcvs1,kcvs1), &
               at1(icvs1,jcvs1,kcvs1),ab1(icvs1,jcvs1,kcvs1),b1(icvs1,jcvs1,kcvs1)
       real :: ap2(icvs2,jcvs2,kcvs2),ae2(icvs2,jcvs2,kcvs2),aw2(icvs2,jcvs2,kcvs2), &
               an2(icvs2,jcvs2,kcvs2),as2(icvs2,jcvs2,kcvs2), &
               at2(icvs2,jcvs2,kcvs2),ab2(icvs2,jcvs2,kcvs2),b2(icvs2,jcvs2,kcvs2)
!
       integer :: is1,ie1,js1,je1,ks1,ke1,is2,ie2,js2,je2,ks2,ke2
       integer :: n,i1,j1,k1,i2,j2,k2 
       logical :: lfirst1,lfirst2,least,lwest,lnorth,lsouth,ltop,lbottom
!
       lfirst1 = .false.;lfirst2 = .false. 
       is1=index1r(1);ie1=index1r(2)
       js1=index1r(3);je1=index1r(4)
       ks1=index1r(5);ke1=index1r(6)
       is2=index2r(1);ie2=index2r(2)
       js2=index2r(3);je2=index2r(4)
       ks2=index2r(5);ke2=index2r(6)
!
       if (p==0) then
!
! Side 1
         if ( ncs1>1 ) then  
           do n = 1,ncs1
!
             i1 = index1(n,1)
             j1 = index1(n,2)
             k1 = index1(n,3)
!
             if ( .not.lfirst1 ) then
               lwest=.false. ;least=.false.
               if (is1==ie1) then
                 if (is1==1) then
                   lwest = .true.
                 else
                   least = .true.
                 endif
               endif  
!
               lsouth=.false. ;lnorth=.false.
               if (js1==je1) then
                 if (js1==1) then
                   lsouth = .true.
                 else
                   lnorth = .true.
                 endif
               endif
!
               lbottom=.false. ;ltop=.false.
               if (ks1==ke1) then
                 if (ks1==1) then
                   lbottom = .true.
                 else
                   ltop = .true.
                 endif
               endif             
               lfirst1 = .true.
             endif
!
! Adjust coefficients
             if (lwest) then
               aw1(i1,j1,k1) = 0.0
             elseif (least) then
               ae1(i1,j1,k1) = 0.0                          
             elseif (lsouth) then
               as1(i1,j1,k1) = 0.0             
             elseif (lnorth) then
               an1(i1,j1,k1) = 0.0             
             elseif (lbottom) then
               ab1(i1,j1,k1) = 0.0             
             elseif (ltop) then
               at1(i1,j1,k1) = 0.0             
             else  
               write(*,*) 'coefh_mg_zero error'
               stop             
             endif
!           
           enddo
         endif
!
! Side 2          
         if ( ncs2>1 ) then       
           do n = 1,ncs2
!
             i2 = index2(n,1)
             j2 = index2(n,2)
             k2 = index2(n,3)
!
             if ( .not.lfirst2 ) then
               lwest=.false. ;least=.false.
               if (is2==ie2) then
                 if (is2==1) then
                   lwest = .true.
                 else
                   least = .true.
                 endif
               endif  
!
               lsouth=.false. ;lnorth=.false.
               if (js2==je2) then
                 if (js2==1) then
                   lsouth = .true.
                 else
                   lnorth = .true.
                 endif
               endif
!
               lbottom=.false. ;ltop=.false.
               if (ks2==ke2) then
                 if (ks2==1) then
                   lbottom = .true.
                 else
                   ltop = .true.
                 endif
               endif             
               lfirst2 = .true.
             endif
!
! Adjust coefficients
             if (lwest) then
               aw2(i2,j2,k2) = 0.0
             elseif (least) then
               ae2(i2,j2,k2) = 0.0                          
             elseif (lsouth) then
               as2(i2,j2,k2) = 0.0             
             elseif (lnorth) then
               an2(i2,j2,k2) = 0.0             
             elseif (lbottom) then
               ab2(i2,j2,k2) = 0.0             
             elseif (ltop) then
               at2(i2,j2,k2) = 0.0             
             else
               write(*,*) 'coefh_mg_zero error'
               stop             
             endif
!
           enddo           
         endif
!
       elseif (p>0) then
!
! Side 1
         if ( ncs1>1 ) then  
           do n = 1,ncs1
!
             i1 = index1(n,1)
             j1 = index1(n,2)
             k1 = index1(n,3)
!
             if ( .not.lfirst1 ) then
               lwest=.false. ;least=.false.
               if (is1==ie1) then
                 if (is1==1) then
                   lwest = .true.
                 else
                   least = .true.
                 endif
               endif  
!
               lsouth=.false. ;lnorth=.false.
               if (js1==je1) then
                 if (js1==1) then
                   lsouth = .true.
                 else
                   lnorth = .true.
                 endif
               endif
!
               lbottom=.false. ;ltop=.false.
               if (ks1==ke1) then
                 if (ks1==1) then
                   lbottom = .true.
                 else
                   ltop = .true.
                 endif
               endif             
               lfirst1 = .true.
             endif
!
! Adjust coefficients
             if (lwest) then
               acoef1(n) = aw1(i1,j1,k1)             
               aw1(i1,j1,k1) = 0.0
             elseif (least) then
               acoef1(n) = ae1(i1,j1,k1)              
               ae1(i1,j1,k1) = 0.0                          
             elseif (lsouth) then
               acoef1(n) = as1(i1,j1,k1)              
               as1(i1,j1,k1) = 0.0             
             elseif (lnorth) then
               acoef1(n) = an1(i1,j1,k1)              
               an1(i1,j1,k1) = 0.0             
             elseif (lbottom) then
               acoef1(n) = ab1(i1,j1,k1)              
               ab1(i1,j1,k1) = 0.0             
             elseif (ltop) then
               acoef1(n) = at1(i1,j1,k1)              
               at1(i1,j1,k1) = 0.0             
             else
               write(*,*) 'coefh_mg_zero error'
               stop             
             endif
!           
           enddo
         endif
!
! Side 2          
         if ( ncs2>1 ) then       
           do n = 1,ncs2
!
             i2 = index2(n,1)
             j2 = index2(n,2)
             k2 = index2(n,3)
!
             if ( .not.lfirst2 ) then
               lwest=.false. ;least=.false.
               if (is2==ie2) then
                 if (is2==1) then
                   lwest = .true.
                 else
                   least = .true.
                 endif
               endif  
!
               lsouth=.false. ;lnorth=.false.
               if (js2==je2) then
                 if (js2==1) then
                   lsouth = .true.
                 else
                   lnorth = .true.
                 endif
               endif
!
               lbottom=.false. ;ltop=.false.
               if (ks2==ke2) then
                 if (ks2==1) then
                   lbottom = .true.
                 else
                   ltop = .true.
                 endif
               endif             
               lfirst2 = .true.
             endif
!
! Adjust coefficients
             if (lwest) then
               acoef2(n) = aw2(i2,j2,k2)
               aw2(i2,j2,k2) = 0.0
             elseif (least) then
               acoef2(n) = ae2(i2,j2,k2)             
               ae2(i2,j2,k2) = 0.0                          
             elseif (lsouth) then
               acoef2(n) = as2(i2,j2,k2)             
               as2(i2,j2,k2) = 0.0             
             elseif (lnorth) then
               acoef2(n) = an2(i2,j2,k2)             
               an2(i2,j2,k2) = 0.0                            
             elseif (lbottom) then
               acoef2(n) = ab2(i2,j2,k2)             
               ab2(i2,j2,k2) = 0.0             
             elseif (ltop) then
               acoef2(n) = at2(i2,j2,k2)             
               at2(i2,j2,k2) = 0.0             
             else
               write(*,*) 'coefh_mg_zero error'
               stop             
             endif
!
           enddo           
         endif       
!       
       else
         write(*,*) 'coefh_mg_zero error'
         stop        
       endif
!
       end subroutine coefh_mg_zero                       