Semi-analytical methods, such as rigorous coupled wave analysis, have been pivotal for numerical analysis of photonic structures. In comparison to other methods, they offer much faster computation, especially for structures with constant cross-sectional shapes (such as metasurface units). However, when the cross-sectional shape varies even mildly (such as a taper), existing semi-analytical methods suffer from high computational cost. We show that the existing methods can be viewed as a zeroth-order approximation with respect to the structure's cross-sectional variation. We instead derive a high-order perturbative expansion with respect to the cross-sectional variation. Based on this expansion, we propose a new semi-analytical method that is fast to compute even in presence of large cross-sectional shape variation. Furthermore, we design an algorithm that automatically discretizes the structure in a way that achieves a user specified accuracy level while at the same time reducing the computational cost.