Modeling and numerical simulations of fractured, vuggy, porus media is a challenging problem which occurs frequently in reservoir engineering. The problem is especially relevant in flow simulations of karst reservoirs where vugs and caves are embedded in a porous rock and are connected via fracture networks at multiple scales. In this paper we propose a unified approach to this problem by using the Stokes-Brinkman equations at the fine scale. These equations are capable of representing porous media such as rock as well as free flow regions (fractures, vugs, caves) in a single system of equations. We then consider upscaling these equations to a coarser scale. The cell problems, needed to compute coarse-scale permeability of Representative Element of Volume (REV) are discussed. A mixed finite element method is then used to solve the Stokes-Brinkman equation at the fine scale for a number of flow problems, representative for different types of vuggy reservoirs. Upscaling is also performed by numerical solutions of Stokes-Brinkman cell problems in selected REVs. Both isolated vugs in porous matrix as well as vugs connected by fracture networks are analyzed by comparing fine-scale and coarse-scale flow fields. Several different types of fracture networks, representative of short- and long-range fractures are studied numerically. It is also shown that the Stokes-Brinkman equations can naturally be used to model additional physical effects pertaining to vugular media such as partial fracture with fill-in by some material and/or fluids with suspended solid particles.