We present a finite-element solution to the 3D electromagnetic forward-modeling problem in the frequency domain. The method is based on decomposing the electric field into vector and scalar potentials in the Helmholtz equation and in the equation of conservation of charge. Edge element and nodal element basis functions were used, respectively, for the vector and scalar potentials. This decomposition was performed with the intention of satisfying the continuity of the tangential component of the electric field and the normal component of the current density across the interelement boundaries, therefore finding an efficient solution to the problem. The computational domain was subdivided into unstructured tetrahedral elements. The system of equations was discretized using the Galerkin variant of the weighted residuals method, with the approximated vector and scalar potentials as the unknowns of a sparse linear system. A generalized minimum residual solver with an incomplete LU preconditioner was used to iteratively solve the system. The solution method was validated using five examples. In the first and second examples, the fields generated by small dipoles on the surface of a homogeneous half-space were compared against their corresponding analytic solutions. The third example provided a comparison with the results from an integral equation method for a long grounded wire source on a model with a conductive block buried in a less conductive half-space. The fourth example concerned verifying the method for a large conductivity contrast where a magnetic dipole transmitter-receiver pair moves over a graphite cube immersed in brine. Solutions from the numerical approach were in good agreement with the data from physical scale modeling of this scenario. The last example verified the solution for a resistive disk model buried in marine conductive sediments. For all examples, convergence of the solution that used potentials were significantly quicker than that using the electric field.