Abstract: This research is concerned with the electromagnetics of plasmonic layered media. First, the surface wave modes, i.e. surface plasmon polaritons, of the source-free stack of layers have been solved. For this purpose, a generalized expanded matrix formulation has been introduced. This formulation has the advantage of being easy in implementation as a computer program. Moreover, the proposed formulation shows explicitly the dependence of the characteristic equation on the physical structure parameters, which facilitates the optimization process of surface modes' propagation characteristics. These characteristics include the phase constant, attenuation constant, and field distribution of each mode. This technique has been applied for pure dielectric structures as well as a variety of plasmonic structures, such as: Insulator/Metal (IM), Insulator/Metal/Insulator (IMI), and Metal/Insulator/Metal (MIM), and a stack of five layers made up of plasmonic and insulator layers. These structures are investigated throughout the thesis. The dependence of the surface plasmon polaritons of these plasmonic media on the physical structure's parameters is explained. Moreover, the dependency on frequency, i.e. dispersion, is also investigated. Second, a new generalized formulation for the spectral domain Green's functions of unit current sources embedded within a stack of plasmonic layers is presented in this thesis. The current source can be either electric or magnetic. It can be oriented either parallel to the layers or perpendicular to them. All electromagnetic fields due to these generalized sources are calculated. The proposed technique utilizes inward and outward recurrence processes in order to solve for the desired Green's functions. The new approach requires the calculation of a limited number of spectral coefficients, from which all types of spectral Green's functions can be calculated. Unit current sources are inserted within the dielectric and plasmonic stacks under investigation. The Green's functions of these sources are calculated and their poles are linked to the surface wave modes previously obtained. The impact of varying the physical parameters of the layer structure on the spectral Green's functions is presented and explained. Finally, the spectral domain Green's functions of the generalized sources are inverted to the spatial domain. For this purpose, the two-level Discrete Complex Image Method (DCIM) is used. This technique is based on expanding the spectral function into a number of exponentials whose spatial counterparts are also in the form of another series of exponentials. Performing the DCIM in two levels makes it possible to deal numerically with the spectral function's asymptotic behavior and to avoid the singular behavior of the function around the poles and branch-cuts. The spatial Green's functions of electric and magnetic unit current sources inserted within the dielectric and plasmonic structures of interest are obtained and discussed. All results presented within the thesis are compared with the literature, whenever possible, and very good agreement is always observed.