  The paper develops a hybrid method for solving a system of
advection--diffusion equations in a bulk domain coupled to advection--diffusion
equations on an embedded surface. A monotone nonlinear finite volume method for
equations posed in the bulk is combined with a trace finite element method for
equations posed on the surface. In our approach, the surface is not fitted by
the mesh and is allowed to cut through the background mesh in an arbitrary way.
Moreover, a triangulation of the surface into regular shaped elements is not
required. The background mesh is an octree grid with cubic cells. As an example
of an application, we consider the modeling of contaminant transport in
fractured porous media. One standard model leads to a coupled system of
advection--diffusion equations in a bulk (matrix) and along a surface
(fracture). A series of numerical experiments with both steady and unsteady
problems and different embedded geometries illustrate the numerical properties
of the hybrid approach. The method demonstrates great flexibility in handling
curvilinear or branching lower dimensional embedded structures.
