A comprehensive methodology for gridding, discretizing, coarsening, and simulating discrete-fracture-matrix models of naturally fractured reservoirs is described and applied. The model representation considered here can be used to define the grid and transmissibilities, at either the original fine scale or at coarser scales, for any connectivity-list-based finite-volume flow simulator. For our fine-scale mesh, we use a polyhedral gridding technique to construct a conforming matrix grid with adaptive refinement near fractures, which are represented as faces of grid cells. The algorithm uses a single input parameter to obtain a suitable compromise between fine-grid cell quality and the fidelity of the fracture representation. Discretization using a two-point flux approximation is accomplished with an existing procedure that treats fractures as lower-dimensional entities (i.e., resolution in the transverse direction is not required). The upscaling method is an aggregation-based technique in which coarse control volumes are aggregates of fine-scale cells, and coarse transmissibilities are computed using a general flow-based procedure. Numerical results are presented for waterflood, sour gas injection, and gas condensate primary production. Coarse-model accuracy is shown to generally decrease with increasing levels of coarsening, as would be expected. We demonstrate, however, that by using our methodology, two orders of magnitude of speedup can be achieved with models that introduce less than about 10% error (with error appropriately defined). This suggests that the overall framework may be very useful for the simulation of realistic discrete-fracture-matrix models.