Simulations of fluid and porous media couplings, multi-particles & multi-scale problems are important and challenging because of different governing equations, muilti-scales and multi-connected domain, and complicated interface conditions such as BJ and BJS relations. Most of numerical methods in the literature are based on direct approaches that leads large system of equations, often coupled and ill-conditioned, with lower order accuracy near the interface or boundary. In this talk, we propose a finite difference approach with unfitted meshes. By introducing several augmented variables along the interface, we can decouple the original problem as several Poisson/Helmholtz equations with intermediate jump conditions in the solution and the normal derivatives. One obvious advantage is that a fast Poisson/Helmholtz solver can be utilized. The augmented variables should be chosen such that the Beavers-Joseph-Saffman (BJS) and other interface conditions are satisfied. Another significant strategy is to enforce the divergence condition at the interface from the fluid side. We have shown that the original and transformed systems are equivalent. Because the interface conditions are enforced in strong form, we have observed second order convergence for both of the velocity and the pressure for our constructed non-trivial analytic solutions with circular interfaces. The proposed new method has also been utilized to simulate different flow/porous media setting with complicated interfaces which leads to some interesting simulations results such as effect of corners, orientation effect etc.