  We propose an unfitted finite element method for flow in fractured porous
media. The coupling across the fracture uses a Nitsche type mortaring, allowing
for an accurate representation of the jump in the normal component of the
gradient of the discrete solution across the fracture. The flow field in the
fracture is modelled simultaneously, using the average of traces of the bulk
variables on the fractured. In particular the Laplace-Beltrami operator for the
transport in the fracture is included using the average of the projection on
the tangential plane of the fracture of the trace of the bulk gradient. Optimal
order error estimates are proven under suitable regularity assumptions on the
domain geometry. The extension to the case of bifurcating fractures is
discussed. Finally the theory is illustrated by a series of numerical examples.
