In the 3D bounded domain \(\Omega \) with smooth boundary \(\partial \Omega \) , we consider stability problem of solutions to the MHD equations for \(U = (u, B)\) denoting the velocity and magnetic field, respectively. Introducing the harmonic vector field \(X_{\tiny{\mbox{har}}}(\Omega) = \{h \in C^\infty (\bar{\Omega }); \text { rot }h=0, \text { div }h=0, h\cdot \nu |_{\partial \Omega }=0\}\) , we first show that \(U_*\equiv (0, B_*)\) with \(B_*\in X_{\tiny{\mbox{har}}}(\Omega)\) gives an equilibrium state, where \(\nu \) denotes the unit outer normal to \(\partial \Omega \) . Next, we prove that if \(B_*\) is small in \(H^{1, \infty }(\Omega )\) and if the initial disturbance \(U_0= (u_0, B_0)\) is sufficiently small in \(L^3(\Omega )\) with \(B_0 - B_*\) perpendicular to \(X_{\tiny{\mbox{har}}}(\Omega)\) , then \(U_*\) is exponentially stable, which means that there exists a global solution \(U(t)= (u(t), B(t))\) with \(U(0) = U_0\) such that \(\Vert U(t) - U_*\Vert _{L^p(\Omega )} = O(e^{-\alpha t})\) for all \(1\le p \le \infty \) with some \(\alpha >0\) as \(t\rightarrow \infty \) .