This paper describes a numerical formulation for calculating wave propagation with high precision in a three-dimensional system. Yee's discretization scheme is used to formulate a frequency domain method that is compatible with the finite-difference time-domain (FDTD) procedure. When the S-matrix satisfies a unitarity (power flow conservation) condition, the method enables arbitrary Smatrix elements to be obtained within a numerical error of less than 10 −8 (2 × 10 −13 ) for double precision format.