Fix bug in nodal solver that uses harmonic averageing (AMReX-Codes#1902)

This is a minor bug that does not affect the final result.  The convergence
rate might be slightly affected.

Author: WeiqunZhang

Src/LinearSolvers/MLMG/AMReX_MLNodeLap_2D_K.H

@@ -305,9 +305,9 @@ Real mlndlap_adotx_ha (int i, int j, int k, Array4<Real const> const& x,
                +   x(i-1,j+1,k)*(facx*sx(i-1,j  ,k)+facy*sy(i-1,j  ,k))
                +   x(i+1,j+1,k)*(facx*sx(i  ,j  ,k)+facy*sy(i  ,j  ,k))
                +   x(i-1,j,k)*(Real(2.0)*facx*(sx(i-1,j-1,k)+sx(i-1,j,k))
-                                   -     facy*(sy(i-1,j-1,k)+sx(i-1,j,k)))
+                                   -     facy*(sy(i-1,j-1,k)+sy(i-1,j,k)))
                +   x(i+1,j,k)*(Real(2.0)*facx*(sx(i  ,j-1,k)+sx(i  ,j,k))
-                                   -     facy*(sy(i  ,j-1,k)+sx(i  ,j,k)))
+                                   -     facy*(sy(i  ,j-1,k)+sy(i  ,j,k)))
                +   x(i,j-1,k)*(   -facx*(sx(i-1,j-1,k)+sx(i,j-1,k))
                         +Real(2.0)*facy*(sy(i-1,j-1,k)+sy(i,j-1,k)))
                +   x(i,j+1,k)*(   -facx*(sx(i-1,j  ,k)+sx(i,j  ,k))
@@ -550,9 +550,9 @@ void mlndlap_gauss_seidel_ha (Box const& bx, Array4<Real> const& sol,
                     + sol(i-1,j+1,k)*(facx*sx(i-1,j  ,k)+facy*sy(i-1,j  ,k))
                     + sol(i+1,j+1,k)*(facx*sx(i  ,j  ,k)+facy*sy(i  ,j  ,k))
                     + sol(i-1,j,k)*(Real(2.0)*facx*(sx(i-1,j-1,k)+sx(i-1,j,k))
-                                        -     facy*(sy(i-1,j-1,k)+sx(i-1,j,k)))
+                                        -     facy*(sy(i-1,j-1,k)+sy(i-1,j,k)))
                     + sol(i+1,j,k)*(Real(2.0)*facx*(sx(i  ,j-1,k)+sx(i  ,j,k))
-                                        -     facy*(sy(i  ,j-1,k)+sx(i  ,j,k)))
+                                        -     facy*(sy(i  ,j-1,k)+sy(i  ,j,k)))
                     + sol(i,j-1,k)*(   -facx*(sx(i-1,j-1,k)+sx(i,j-1,k))
                              +Real(2.0)*facy*(sy(i-1,j-1,k)+sy(i,j-1,k)))
                     + sol(i,j+1,k)*(   -facx*(sx(i-1,j  ,k)+sx(i,j  ,k))