Multiphase flow in highly heterogeneous porous media holds substantial importance in the petroleum industry since there are large oil volumes in vuggy and fractured reservoirs. In this study, Brinkman equation is employed to model two-phase flow in highly heterogeneous porous media. Brinkman model represents flow dynamics through the porous matrix and vugs using a single differential equation without the need for interfacial conditions. However, accurately characterizing the vugular and porous matrix regions remains a challenging task, especially in reservoir-scale models that typically represent geological formations with large areas. Reservoir simulations often involve models discretized into computational cells, each with dimensions on the order of several meters that are not fine enough do describe in detail the geometry of the vugs and fractures embedded in the porous matrix. This investigation explores two distinct approaches for generating equivalent properties of vugular media. The first approach uses results from the Brinkman model to derive equivalent properties for the single-continuum Darcy model at the same level of discretization. The second approach employs upscaling by numerical optimization to generate equivalent properties on a coarse scale. This method is applied to both single-continuum and dual-continuum Darcy models. Singlecontinuum models can effectively describe vugular structures characterized by dispersed vugs in the porous matrix without strong fluid channeling. However, for systems exhibiting substantial channeling, upscaling can only be achieved through dual-continuum models. In both proposed approaches, the determination of equivalent absolute permeability and equivalent absolute permeability curves is imperative for an efficient upscaling process.