This paper presents an accurate and efficient modeling approach for high-fidelity 3D stress prediction and failure initiation and propagation modeling of tapered composite flexbeams with multiple ply drops subjected to multiaxial loading. To achieve the maximum level of modeling versatility and solution capability, four new modeling capabilities are developed: (1) mesh generation of a tapered laminate, (2) accurate extraction of the stress components at potential delamination interfaces, (3) nonlocal stress characterization for delamination initiation prediction at a ply drop, and (4) application of a discrete damage informed continuum damage model for failure prediction under multiaxial loading. To support the toolkit development and demonstration, six representative flexbeams with 129 plies and 30 ply drops are designed, fabricated, and tested under combined tension and bending loads with Digital Image Correlation (DIC) and strain gauges to monitor damage evolution. Data on the load-displacement curves and the final fracture patterns are collected for the validation of the high-fidelity toolkit. In addition to the comparison of the load-displacement curves, the simulated damage patterns at different loading stages are used to reveal the failure location and failure mode that resulted in the final rupture.