anchor.diff:            call mapc2p(xcell - 0.5d0*dx, ycell + 0.5d0*dy, zcell - 0.5d0*dz, xpcorn(4), ypcorn(4), zpcorn(4))
             ! normal information and area ratio are computed
             ! note normal will be in xz plane
  -          mag = ((ypcorn(3) - ypcorn(1))*(zpcorn(2) - zpcorn(4)) - (ypcorn(2) - ypcorn(4))*(zpcorn(3) - ypcorn(1)))**2
  +          mag = ((ypcorn(3) - ypcorn(1))*(zpcorn(2) - zpcorn(4)) - (ypcorn(2) - ypcorn(4))*(zpcorn(3) - zpcorn(1)))**2
             mag = mag + ((xpcorn(3) - xpcorn(1))*(ypcorn(2) - ypcorn(4)) - (xpcorn(2) - xpcorn(4))*(ypcorn(3) - ypcorn(1)))**2
             mag = dsqrt(mag)
  -          aux(6,i,j,k) = ((ypcorn(3) - ypcorn(1))*(zpcorn(2) - zpcorn(4)) - (ypcorn(2) - ypcorn(4))*(zpcorn(3) - ypcorn(1)))/mag
  +          aux(6,i,j,k) = ((ypcorn(3) - ypcorn(1))*(zpcorn(2) - zpcorn(4)) - (ypcorn(2) - ypcorn(4))*(zpcorn(3) - zpcorn(1)))/mag
             aux(7,i,j,k) = ((xpcorn(3) - xpcorn(1))*(ypcorn(2) - ypcorn(4)) - (xpcorn(2) - xpcorn(4))*(ypcorn(3) - ypcorn(1)))/mag
             aux(8,i,j,k) = 0.5d0*mag/(dy*dz)
   
  
positive.code:           call mapc2p(xcell + 0.5d0*dx, ycell - 0.5d0*dy, zcell - 0.5d0*dz, xpcorn(4), ypcorn(4), zpcorn(4))
            ! for this face, the normal information is needed in addition to the area ratio
            ! the normal will reside in the xz plane
            mag = ((ypcorn(3) - ypcorn(1))*(zpcorn(2) - zpcorn(4)) - (ypcorn(2) - ypcorn(4))*(zpcorn(3) - ypcorn(1)))**2
            mag = mag + ((xpcorn(3) - xpcorn(1))*(ypcorn(2) - ypcorn(4)) - (xpcorn(2) - xpcorn(4))*(ypcorn(3) - ypcorn(1)))**2
            mag = dsqrt(mag)
            aux(10,i,j,k) = ((ypcorn(3) - ypcorn(1))*(zpcorn(2) - zpcorn(4)) - (ypcorn(2) - ypcorn(4))*(zpcorn(3) - ypcorn(1)))/mag
            aux(11,i,j,k) = ((xpcorn(3) - xpcorn(1))*(ypcorn(2) - ypcorn(4)) - (xpcorn(2) - xpcorn(4))*(ypcorn(3) - ypcorn(1)))/mag
            aux(12,i,j,k) = 0.5d0*mag/(dx*dy)
  
  
negative.code:           call mapc2p(xcell + 0.5d0*dx, ycell - 0.5d0*dy, zcell + 0.5d0*dz, xpcorn(3), ypcorn(3), zpcorn(3))
            call mapc2p(xcell + 0.5d0*dx, ycell - 0.5d0*dy, zcell - 0.5d0*dz, xpcorn(4), ypcorn(4), zpcorn(4))
            ! only need area ratio from cross-product of diagonals, which will point in the y direction
            mag = (xpcorn(3) - xpcorn(1))*(zpcorn(2) - zpcorn(4)) - (xpcorn(2) - xpcorn(4))*(zpcorn(3) - zpcorn(1))
            aux(9,i,j,k) = 0.5d0*mag/(dx*dz)
  
            ! compute mapping info for lower face in z direction
  