The regional-scale erosion of hundreds to thousands of meters of stratigraphy remain one of the most pivotal yet debated arguments within Flood Geology. Despite the many previous studies of erosion of anticlines as one illustration of regional-scale erosion, no rigorous mathematical model has been proposed for systematic and repeatable modeling of fold surfaces for erosion estimates in data scarce locations. As such, this study proposes a boundary-value problem approach for modeling symmetric and non-verging fold systems. Applied to folds in the Mount St. Helens region, the model performed well in describing characteristics of the half-wavelength of the fold system. The modeled surface resulted in a calculation of 6.16 km of vertical relief eroded from the current topography, a value that could be increased to 10.1 km when transferring the modeled surface to the outermost observed fold surface. Site-specific geology suggests an additional 1 to 4 km of stratigraphy may have rested atop this modeled surface. This application of boundary-value problems represents a promising technique to systematically reconstruct fold systems for erosion estimates. The approach requires minimal inputs that are easily acquired from geologic maps although this limits its application to approximately symmetric and non-verging fold systems. Even so, this technique represents a first step towards developing an easily deployable yet rigorous approach to model fold systems for repeatable and consistent erosion estimates.