We present a second kind integral equation (SKIE) formulation

for calculating the electromagnetic modes of optical waveguides,

where the unknowns are only on material interfaces. The resulting numerical

algorithm can handle optical waveguides with a large number of inclusions of 

arbitrary irregular cross section. It is capable of finding the bound, leaky, and 

complex modes for optical fibers and waveguides including photonic crystal 

fibers (PCF), dielectric fibers and waveguides. Most importantly, the formulation 

is well conditioned even in the case of nonsmooth geometries. Our method is highly 

accurate and thus can be used to calculate the propagation loss of the electromagnetic 

modes accurately, which provides the photonics industry a reliable tool for the design 

of more compact and efficient photonic devices. We illustrate and validate the 

performance of our method through extensive numerical studies and by comparison 

with semi-analytical results and previously published results.