This paper studies an integration of the forward products flow, as in the vehicle routing problem with crossdocking (VRPCD), and the reverse products flow, as in the VRP with reverse cross-docking (VRP-RCD), namely the VRP with forward and reverse cross-docking (VRP-FRCD). It is modelled as a mixed-integer linear programming model that covers a four-level supply chain network, including suppliers, retailers, outlets, and a cross-dock facility. The objective is to determine the appropriate number of vehicles used to complete the delivery process, together with the route sequence of every vehicle, such that the total fixed vehicle and distance-related costs are minimized. In order to solve the problem, a two-phase matheuristic is developed. The first phase focuses on finding as many route combinations as possible. This is implemented by an adaptive large neighborhood search (ALNS) algorithm. Subsequently, a set partitioning formulation is formulated and solved in the second phase to determine the best combination over all routes found in the first phase that minimizes the total cost. The performance of the matheuristic in solving the newly developed benchmark instances is compared against those of CPLEX and a pure ALNS algorithm. Experimental results suggest that the matheuristic is more powerful and beneficial than both ALNS and CPLEX to obtain high-quality solutions within an acceptable computational time.