Mechanical nonlinearities dominate the motion of nanoresonators already at relatively small oscillation amplitudes. Although single and coupled two-degree-of-freedom models have been used to account for experimentally observed nonlinear effects, it is shown that these models quickly deviate from experimental findings when multiple modes influence the nonlinear response. Here, we present a nonlinear reduced-order modeling methodology based on finite-element method simulations for capturing the global nonlinear dynamics of nanomechanical resonators. Our physics-based approach obtains the quadratic and cubic nonlinearities of resonators over a wide frequency range that spans 70 MHz. To qualitatively validate our approach, we perform experiments on a graphene nanodrum driven optothermally and show that the model can replicate diverse ranges of nonlinear phenomena, including multistability, parametric resonance, and different internal resonances without considering any empirical nonlinear fitting parameters. By providing a direct link between microscopic geometry, material parameters, and nonlinear dynamic response, we clarify the physical significance of nonlinear parameters that are obtained from fitting the dynamics of nanomechanical systems, and provide a route for designing devices with desired nonlinear behavior.