This paper reports a new numerical scheme to simulate the radio-frequency (RF) induced RF sheath, which is suitable for a large 3D simulation. In the RF sheath boundary model, the tangential component of the electric field () is given by the gradient of a scalar electric field potential. We introduce two additional scalar potentials for the tangential components of the magnetic field, which effectively impose the normal electric displacement (D n ) on the plasma sheath boundary condition via in-homogeneous Neumann boundary condition and constrain the tangential electric field on the surface as curl-free (). In our approach, the non-linear sheath impedance is formulated as a natural extension of the large thickness (or asymptotic) sheath limit (), allowing for handling both asymptotic and non-linear regimes seamlessly. The new scheme is implemented using the Petra-M finite element method analysis framework and is verified with simulations in the literature. The significance of non-linearity is discussed in various plasma conditions. An application of this scheme to asymptotic RF sheath simulation on the WEST ICRF antenna side limiters is also discussed.