Graphene nanoplatelets (GNP), dispersed in epoxy resins, have been reported to show remarkable improvements in fracture toughness. This work seeks to investigate the nanoscale mechanisms near the crack tip in the presence of graphene nanoplatelets that promote the improvement in fracture toughness. The size of the fracture process zone formed in epoxy near a crack, typically on the order of tens of microns, is beyond the range of pure molecular dynamics (MD) modeling. Hence, MD simulations, concurrently coupled with finite element method (FEM) to bridge length and time scales, are used to simulate the crack initiation behavior under fracture loading. A novel concurrent multiscale coupling methodology is used, that partitions the full simulation domain into MD and FEM subdomains, and iteratively computes solutions for a given load step. The nonlinear behavior of the polymer in the MD domain is approximated in the FEM domain using a Ramberg–Osgood constitutive model, parameterized based on MD representative volume elements (RVE). K-dominant displacement fields, derived from linear elastic fracture mechanics (LEFM) are used as boundary conditions, away from the crack tip and surrounding the fracture process zone (FPZ). Stresses and strains in the MD domain (near-crack-tip region) are compared at a given level of far-field loading for baseline specimens and GNP-embedded specimens. The overarching goal is to be able to predict fracture toughness at crack initiation at room temperature in a brittle epoxy polymer. The effect of GNP size and position on crack tip shielding and hence on enhanced fracture properties will be investigated and the quantified using MD simulations.