Multiscale Finite Volume Research Articles

A combined Markov Chain Monte Carlo and Levenberg–Marquardt inversion method for heterogeneous subsurface reservoir modeling

In this study, we systematically investigate an inverse method for heterogeneous porous media to obtain porosity and permeability fields considering only production well data. The forward modeling of the input data, namely pressure and saturation fields, is based on the motion equations of a coupled Darcy flow system involving two phases of isothermal fluid flow. We discretize these equations using a multiscale finite volume simulation technique in the spatial domain, and the backward Euler method in time domain. In the inversion procedure, we combine a global optimization method, the Markov Chain Monte Carlo (MCMC) method, with a local optimization method, the Levenberg–Marquardt (LM) method. The MCMC was implemented as the Random Walk algorithm, and to generate the samples of the porosity and permeability fields, we employed Karhunen–Loève (KL) expansion of second-order stationary fields with Gaussian covariance. The coefficients of the KL expansion are estimated by minimizing the norm of the residual dependent on these fields. At the end of the MCMC iterations, we refine the KL coefficients associated with the porosity and permeability fields using the LM method. We verify in the numerical experiments that the accuracy of the porosity and permeability fields is improved by the LM refinement step.

A Multi-Scale Numerical Simulation Method Considering Anisotropic Relative Permeability

Most of the oil reservoirs in China are fluvial deposits with firm reservoir heterogeneity, where differences in fluid flow capacity in individual directions should not be ignored; however, the available commercial reservoir simulation software cannot consider the anisotropy of the relative permeability. To handle this challenge, this paper takes full advantage of the parallelism of the multi-scale finite volume (MsFV) method and establishes a multi-scale numerical simulation approach that incorporates the effects of reservoir anisotropy. The methodology is initiated by constructing an oil–water black-oil model considering the anisotropic relative permeability. Subsequently, the base model undergoes decoupling through a sequential solution, formulating the pressure and transport equations. Following this, a multi-scale grid system is configured, within which the pressure and transport equations are progressively developed in the fine-scale grid domain. Ultimately, the improved multi-scale finite volume (IMsFV) method is applied to mitigate low-frequency error in the coarse-scale grid, thereby enhancing computational efficiency. This paper introduces two primary innovations. The first is the development of a multi-scale solution method for the pressure equation incorporating anisotropic relative permeability. Validated using the Egg model, a comparative analysis with traditional numerical simulations demonstrates a significant improvement in computational speed without sacrificing accuracy. The second innovation involves applying the multi-scale framework to investigate the impact of anisotropy relative permeability on waterflooding performance, uncovering distinct mechanisms by which absolute and relative permeability anisotropy influence waterflooding outcomes. Therefore, the IMsFV method can be used as an effective tool for high-resolution simulation and precise residual oil prediction in anisotropic reservoirs.

A 3-D extension of the Multiscale Control Volume method for the simulation of the single-phase flow in anisotropic and heterogeneous porous media

The level of detail on modern geological models requires higher resolution grids that may render the simulation of multiphase flow in porous media intractable. Moreover, these models may comprise highly heterogeneous media with phenomena taking place in different scales. Scale transferring techniques allow for the solution of such problems in a lower resolution scale at reduced computational cost. Among these techniques, the Multiscale Finite Volume (MsFV) method constructs a set of numerical operators in order to map quantities from the fine-scale mesh to a coarser one and vice-versa while maintaining flux conservation on both scales. However, the MsFV formulation, as originally stated, is only consistent on k-orthogonal grids since it uses a linear Two-point Flux Approximation (TPFA) method and may struggle to generate consistent primal-dual coarse grids pairs on unstructured grids. The Multiscale Restriction Smoothed-Basis method (MsRSB) improves on the MsFV by introducing a new iterative procedure to find the multiscale operators and the concept of support regions which reduces the method's complexity when applied to unstructured fine and coarse grids. The original version was only consistent on k-orthogonal fine grids due to the TPFA discretization, but filtering methods have been developed to also enable consistent multipoint schemes on the fine scale. Meanwhile, the Multiscale Control Volume method (MsCV) replaces the TPFA by the Multipoint Flux Approximation with a Diamond stencil (MPFA-D) scheme on the fine-scale while further enhancing the generation of the geometric entities to allow unstructured grids on the fine and coarse scales for two-dimensional simulations. In this work we propose an extension to three-dimensional geometries of both the MsCV and the algorithm to obtain the multiscale geometric entities based on the concept of a background grid, a coarser grid used as a proxy for the primal coarse grid. We modify the original MPFA-D method to use the very robust Global Least Squares (GLS) interpolation technique to obtain the required auxiliary nodal unknowns. We also introduce an enhanced version of the 3-D MsCV with the incorporation of the enhanced MsRSB (E-MsRSB) to enforce M-matrix properties and improve convergence. Finally, we employ the MsCV operators in a two-stage smoother to show how it can be used as a good iterative procedure to recover the fine-scale solution. We show that the 3-D MsCV method produces good results employing unstructured grids on both scales to handle the simulation of the single-phase flow in anisotropic, and heterogeneous porous media and that the smoothing procedure is able to accurately retrieve the fine-scale solution within a few iterations.

Multiscale simulation of highly heterogeneous subsurface flow with physically adaptive 3D unstructured grids

The multiscale finite volume framework is extended with physical-based grid adaptation (MSFV-GA) for three-dimensional fully unstructured grids. The MSFV method produces large errors for highly heterogeneous and anisotropic problems on standard coarse grids generated based on geometrical features. Adding constraints based on the fine-scale physical properties to construct coarse grids can increase the accuracy and robustness of the MSFV solutions. A novel and automated coarse grid generation strategy adaptable to the reservoir heterogeneity for 3D unstructured grids is presented. Adaptive coarse grids are generated based on local variations of fine scale permeability field. The proposed procedure is very flexible respect to the local physical and geometrical limitations in coarse blocks. Multiscale coarse grid adaptation notably improves the basis functions calculations. The performance of the MSFV-GA in handling challenging problems involving high media contrasts or impermeable shale layers is assessed. The results are very promising and show that the MSFV-GA can accurately approximate complex physics problems. The presented 3D unstructured procedure can notably improve the applicability of the multiscale finite volume method and provide a promising framework for the next generation of reservoir flow simulators.

Algorithmic monotone multiscale finite volume methods for porous media flow

Multiscale finite volume methods are known to produce reduced systems with multipoint stencils which, in turn, could give non-monotone and out-of-bound solutions. We propose a novel solution to the monotonicity issue of multiscale methods. The proposed algorithmic monotone (AM- MsFV/MsRSB) framework is based on an algebraic modification to the original MsFV/MsRSB coarse-scale stencil. The AM-MsFV/MsRSB guarantees monotonic and within bound solutions without compromising accuracy for various coarsening ratios; hence, it effectively addresses the challenge of multiscale methods' sensitivity to coarse grid partitioning choices. Moreover, by preserving the near null space of the original operator, the AM-MsRSB showed promising performance when integrated in iterative formulations using both the control volume and the Galerkin-type restriction operators. We also propose a new approach to enhance the performance of MsRSB for MPFA discretized systems, particularly targeting the construction of the prolongation operator. Results show the potential of our approach in terms of accuracy of the computed basis functions and the overall convergence behavior of the multiscale solver while ensuring a monotone solution at all times.

Two‐level multiscale enrichment finite volume method based on the Newton iteration for the stationary incompressible magnetohydrodynamics flow

In this paper, we propose a stabilized finite volume method based on the Newton's type iteration for the steady incompressible magnetohydrodynamics (MHD) problem. In order to reduce the computational complexity, the lowest order mixed finite element pair ( ‐ ‐ ) is adopted to approximate the velocity, pressure and magnetic fields, and the multiscale enrichment method is introduced to overcome the restriction of discrete inf‐sup (LBB) condition. We firstly solve the incompressible MHD equations by the Newton iterations on a coarse mesh with the mesh size , stability and convergence results of numerical solutions are provided. Secondly, the combination of the two‐level method with the Newton iteration is used to approximate the considered problem, a coupled nonlinear problem is solved on the coarse mesh , then the Stokes and Maxwell equations are considered on a fine mesh with the mesh size , the uniform stability and optimal error estimates of two‐level Newton iterative method are provided. Theoretical findings show that the two‐level method has the same order as the one‐level method in ‐norm as long as the mesh sizes satisfy . However, the two‐level method involves much less work than the one‐level method. Finally, some numerical examples are presented to demonstrate the effectiveness of the considered numerical schemes.

A fixed point multi-scale finite volume method: Application to two-phase incompressible fluid flow through highly heterogeneous porous media

An enhanced Multi-scale Finite Volume (MsFV) method is proposed to efficiently simulate two phase flow through highly heterogeneous porous media. Here, the accuracy of the MsFV method is significantly improved by enhancing its so-called basis functions, while its computational cost remains in the same order of the basic MsFV method. First, a fixed point is defined along each edge of local problems producing the basis functions and a proper basis function value at this point is estimated based on the local absolute permeability data. Then, the variable boundary condition is independently calculated for both sides of each fixed point. The proposed numerical method is validated by using the analytic solution of a heterogeneous problem. In addition, the convergence of the MsFV method is discussed. Considering highly heterogeneous permeability domains derived from 35 top layers of the tenth SPE comparative study problem, the accuracy of different localization schemes is compared for a set of imbibition problems with different global boundary conditions and mobility ratios. Numerical results have indicated that the overall accuracy of multi-scale velocity solutions is increased noticeably (up to 45%) by using the proposed localization scheme in the MsFV method.

A Hybrid Multiscale Finite Cloud Method and Finite Volume Method in Solving High Gradient Problem

Meshfree methods are always the alternatives to be considered besides the finite element method in dealing with complex engineering problems. Among these complex problems, there are problems where discontinuities or high gradients may only emerge at a small region of a domain and the remaining parts of the domain are smooth and continuous. In this case, adding nodes to the critical region is essential for the better accuracy of solutions. It would be effective if a meshfree method is applied to the critical region and a mesh-based method is employed to solve the remaining part of the domain since node refinement is easier for a meshfree method whereas a mesh-based method is simple to develop and capable of producing accurate solutions for smooth regions. In this study, the finite cloud method (FCM) is coupled with the finite volume method (FVM) such that node refinement only takes place at the region where the FCM is applied. A handshaking approach is proposed to enforce the continuity near the interface where both methods are applied. Numerical examples and computation methods are presented in this paper and the numerical results are examined.

Scalable Graphics Processing Unit–Based Multiscale Linear Solvers for Reservoir Simulation

Summary In this work, the scalability of two key multiscale solvers for the pressure equation arising from incompressible flow in heterogeneous porous media, namely, the multiscale finite volume (MSFV) solver, and the restriction-smoothed basis multiscale (MsRSB) solver, are investigated on the graphics processing unit (GPU) massively parallel architecture. The robustness and scalability of both solvers are compared against their corresponding carefully optimized implementation on the shared-memory multicore architecture in a structured problem setting. Although several components in MSFV and MsRSB algorithms are directly parallelizable, their scalability on the GPU architecture depends heavily on the underlying algorithmic details and data-structure design of every step, where one needs to ensure favorable control and data flow on the GPU, while extracting enough parallel work for a massively parallel environment. In addition, the type of algorithm chosen for each step greatly influences the overall robustness of the solver. Thus, we extend the work on the parallel multiscale methods of Manea et al. (2016) to map the MSFV and MsRSB special kernels to the massively parallel GPU architecture. The scalability of our optimized parallel MSFV and MsRSB GPU implementations are demonstrated using highly heterogeneous structured 3D problems derived from the SPE10 Benchmark (Christie and Blunt 2001). Those problems range in size from millions to tens of millions of cells. For both solvers, the multicore implementations are benchmarked on a shared-memory multicore architecture consisting of two packages of Intel® Cascade Lake Xeon Gold 6246 central processing unit (CPU), whereas the GPU implementations are benchmarked on a massively parallel architecture consisting of NVIDIA Volta V100 GPUs. We compare the multicore implementations to the GPU implementations for both the setup and solution stages. Finally, we compare the parallel MsRSB scalability to the scalability of MSFV on the multicore (Manea et al. 2016) and GPU architectures. To the best of our knowledge, this is the first parallel implementation and demonstration of these versatile multiscale solvers on the GPU architecture. NOTE: This paper is also published as part of the 2021 SPE Reservoir Simulation Conference Special Issue.

An algebraic multiscale solver for the simulation of two-phase flow in heterogeneous and anisotropic porous media using general unstructured grids (AMS-U)

We introduce a new multiscale finite volume framework for the simulation of multi-phase flow in heterogenous and anisotropic porous media on quite general unstructured grids that enable geophysical grid defined properties to be used directly on a high definition grid.A novel background grid strategy is presented, where an auxiliary mesh aids the creation of multiscale primal and dual coarse grids. Due to the unstructured grid connectivity, the dual grids do not retain the natural structured grid decomposition. This may induce errors in the computed basis functions, and contribute to loss of mass conservation. A generalization of the Algebraic Multiscale Solver (AMS) is proposed to prevent the basis function induced leakage outside the support region enabling coupling with a Control Volume Distributed MultiPoint Flux Approximation (CVD-MPFA), ensuring consistency on non k-orthogonal unstructured grids.We present three different problems in which the accuracy of the framework is validated by comparing the multiscale solver with the direct simulation on the fine-scale. The results obtained show that the method can produce well resolved solutions for two-phase flow in highly heterogeneous and anisotropic porous media using general unstructured grids on all scales. As a consequence, the method captures the most important flow features that result on high-resolution geophysical models while providing suitable discretization for complex geological formations found in current petroleum reservoir problems.The novelty of this scheme involves three components which are combined to enable consistent multiscale simulation on unstructured grids: i) a new approach to primal and dual coarse grid generation, ii) a novel technique to prevent basis function induced leakage and iii) the coupling of the Algebraic Multiscale Solver (AMS) with a Multipoint Flux approximation with a Diamond Stencil (MPFA-D).

Multiscale finite volume method with adaptive unstructured grids for flow simulation in heterogeneous fractured porous media

Multiscale finite volume method with adaptive unstructured grids for flow simulation in heterogeneous fractured porous media

Upscaling numérico baseado em uma metodologia multiescala para a simulação de escoamentos em reservatórios de petróleo heterogêneos e anisotrópicos / Numerical upscaling based on a multiscale methodology for simulating flows in heterogeneous and anisotropic oil reservoirs

This work aims to develop a scale transfer methodology in highly heterogeneous and anisotropic porous media, that makes it possible to solve elliptical problems with accuracy and low computational cost. A scale transfer method is developed that uses the Multiscale Finite Volume Method (MSFVM) to perform a flow based numerical upscaling of absolute permeability. The transmissibility matrix was built in the coarse mesh over the projection and restriction operators obtaining permeability in the coarse scale. This allows to solve only local problems in the coarse mesh volumes and to use pressure multiscale solution in the fine mesh. The permeability field can be obtained by different strategies, and enabling simulations directly on the coarse mesh, controlling the accuracy of the solutions and also with the possibility of updating the permeability field over the recalculated locally conservative multiscale solution.

Simulation of two-phase incompressible fluid flow in highly heterogeneous porous media by considering localization assumption in multiscale finite volume method

Two-phase incompressible fluid flow through highly heterogeneous porous media is simulated by using the Multiscale Finite Volume (MsFV) method. Effects of the localization assumption on the accuracy of the MsFV are investigated by comparing the results associated with different boundary conditions of local problems producing the basis functions. The total number of six boundary conditions of two general types, including Dirichlet and Dirichlet-Neumann types, are compared. For the former, the linear, variable (reduced), and step-type boundary conditions are considered and a modified variable boundary condition is proposed. For the latter, a basic and a step-type Neumann-Dirichlet boundary condition are suggested. To estimate the errors in the MsFV solutions for continuous problems, a heterogeneous two-dimensional problem with continuous permeability field is designed and solved analytically. Synthetic two-scale permeability fields as well as highly heterogeneous random fields are used to assess the accuracy of the MsFV solutions with different localization schemes, in comparison with the fine-scale reference solution. Numerical results indicate that the modified variable boundary condition, with a proper value of its weighting factor, can generally produce the most accurate results, when compared with the other localization schemes.

A benchmark study of the multiscale and homogenization methods for fully implicit multiphase flow simulations

Accurate simulation of multiphase flow in subsurface formations is challenging, as the formations span large length scales (km) with high-resolution heterogeneous properties. To deal with this challenge, different multiscale methods have been developed. Such methods construct coarse-scale systems, based on a given high-resolution fine-scale system. Furthermore, they are amenable to parallel computing and allow for a-posteriori error control. The multiscale methods differ from each other in the way the transition between the different scales is made. Multiscale (finite element and finite volume) methods compute local basis functions to map the solutions (e.g. pressure) between coarse and fine scales. Instead, homogenization methods solve local periodic problems to determine effective models and parameters (e.g. permeability) at a coarser scale. It is yet unknown how these two methods compare with each other, especially when applied to complex geological formations, with no clear scale separation in the property fields. This paper develops the first comparison benchmark study of these two methods and extends their applicability to fully implicit simulations using the algebraic dynamic multilevel (ADM) method. At each time step, on the given fine-scale mesh and based on an error analysis, the fully implicit system is solved on a dynamic multilevel grid. The entries of this system are obtained by using multiscale local basis functions (ADM-MS), and, respectively, by homogenization over local domains (ADM-HO). Both sets of local basis functions (ADM-MS) and local effective parameters (ADM-HO) are computed at the beginning of the simulation, with no further updates during the multiphase flow simulation. The two methods are extended and implemented in the same open-source DARSim2 simulator (https://gitlab.com/darsim2simulator), to provide fair quality comparisons. The results reveal insightful understanding of the two approaches, and qualitatively benchmark their performance. It is re-emphasized that the test cases considered here include permeability fields with no clear scale separation. The development of this paper sheds new lights on advanced multiscale methods for simulation of coupled processes in porous media.

A Multiscale Control Volume framework using the Multiscale Restriction Smooth Basis and a non-orthodox Multi-Point Flux Approximation for the simulation of two-phase flows on truly unstructured grids

In these days, the advances in the fields of geology have led to geocellular models containing up to hundreds millions cells. Meanwhile, the standard petroleum reservoir simulations are limited to models with only a fraction of this amount. Thus, the multiples simulation of multiphase flow in heterogeneous and anisotropic media on the geological scale are simply impractical. To overcome this problem Multiscale Finite Volume (MsFV) uses operators to project the fine-scale system of equations onto a coarse-scale space originated from a lower-resolution grid. The resulting coarse system is solved, and by using the multiscale operator projected back onto the higher-resolution grid. Notwithstanding, the MsFV is not capable of capturing geologic properties on unstructured grids, as in its classical form, it relies on a Two Point Flux Approximation (TPFA) for fluxes, which is only consistent for k-orthogonal grids. Furthermore, the algorithm used by the MsFV to generate the geometric entities, such as, coarse-scale and dual meshes and its centers are not valid for general unstructured grids. The Multiscale Restricted Smoothed Basis (MsRSB) method expands this framework by creating the algorithms capable of generating these structures and by proposing an iterative multiscale operator capable of handling simulation on unstructured coarse grids. However, it fails to produce consistent solutions on fine-scale unstructured grids and for arbitrary permeability tensors, as it uses the TPFA method. In this article, we couple a Multi-Point Flux Approximation with a Diamond stencil to the MsRSB, to extend its use to general unstructured grids. We call our framework the “Multiscale Control-Volume” (MsCV) method. Moreover, we bring the iterative Dirichlet pressure correction scheme and the upscaling flux correction to the unstructured multiscale context. On the first examples, we study the impact of these corrections on the simulation of one-phase flow using the SPE Model 2 benchmark. In addition, we show the behavior of the MsCV framework on the simulation of two-phase flows on general unstructured grids and highly heterogeneous and isotropic permeability fields.

Conservative multirate multiscale simulation of multiphase flow in heterogeneous porous media

Accurate and efficient simulation of multiphase flow in heterogeneous porous media motivates the development of space-time multiscale strategies for the coupled nonlinear flow (pressure) and saturation transport equations. The flow equation entails heterogeneous high-resolution (fine-scale) coefficients and is global (elliptic or parabolic). The time-dependent saturation profile, on the other hand, may exhibit sharp local gradients or discontinuities (fronts) where the solution accuracy is highly sensitive to the time-step size. Therefore, accurate flow solvers need to address the multiscale spatial scales, while advanced transport solvers need to also tackle multiple time scales. This paper presents the first multirate multiscale method for space-time conservative multiscale simulation of sequentially coupled flow and transport equations. The method computes the pressure equation at the coarse spatial scale with a multiscale finite volume technique, while the transport equation is solved by taking variable time-step sizes at different locations of the domain. At each coarse time step, the developed local time-stepping technique employs an adaptive recursive time step refinement to capture the fronts accurately. The applicability (accuracy and efficiency) of the method is investigated for a wide range of two-phase flow simulations in heterogeneous porous media. For the studied cases, the proposed method is found to provide 3 to 4 times faster simulations. Therefore, it provides a promising strategy to minimize the tradeoff between accuracy and efficiency for field-scale applications.

Unstructured grid adaptation for multiscale finite volume method

The multiscale finite volume (MSFV) method produces non-monotone solutions in the presence of highly heterogeneous and channelized permeability fields, or in media with impermeable barriers. The accuracy of the MSFV method can be significantly improved by generating adapted coarse grids according to the fine-scale permeability field. This paper presents efficient algorithms for generating adaptive unstructured primal and dual coarse grids based on permeability features to efficiently improve the MSFV results. The primal coarse grid is generated based on a multilevel tabu search algorithm and its boundaries are adapted to reduced non-physical coarse-scale transmissibilities. Adaptive dual coarse grid is generated based on Dijkstra’s routing algorithm. The performance of the proposed algorithms is assessed for challenging test cases with highly heterogeneous and channelized permeability fields as well as impermeable shale layers. Numerical results show that permeability-adapted coarse grids significantly improve the accuracy of the MSFV method for simulation of multiphase flow in highly heterogeneous porous media.

Multiscale Finite Element Simulation of Thermal Properties and Mechanical Strength of Reduced Graphene Oxide Reinforced Aluminium Matrix Composite

The effective properties of metal matrix composites (MMCs) depend on matrix material and reinforcement property specifications as well as bonding at interphase. The use of numerical methods such as finite element (FE) and mean field homogenization (MFH) can assist in predicting MMC properties thus reducing time and cost of optimizing composite properties through experiments. In the present work, a multiscale representative volume element (RVE) of the microstructure of reduced graphene oxide (rGO) reinforced Aluminium (Al) matrix composite (rGO/Al) is created in MSC DigiMat and analysed using Abaqus software. The effect of porosity and rGO reinforcement on thermal conductivity and strength of the rGO/Al composites is studied. The variation in thermal conductivity between FE-RVE and experimental data is a maximum of 2.2% and a minimum of 0.07% for rGO reinforcement of 1 wt.% and 3 wt.% respectively. The results show good agreement between FE-RVE simulation, MFH and experimental data. This approach can provide an efficient technique for selecting matrix and reinforcement phase properties for MMC fabrication. Keywords: Al/rGO composite, Multiscale finite element-representative volume, Thermal and mechanical properties

Nonlinear acceleration of sequential fully implicit (SFI) method for coupled flow and transport in porous media

The sequential fully implicit (SFI) method was introduced along with the development of the multiscale finite volume (MSFV) framework, and has received considerable attention in recent years. Each time step for SFI consists of an outer loop to solve the coupled system, in which there is one inner Newton loop to implicitly solve the pressure equation and another loop to implicitly solve the transport equations. Limited research has been conducted that deals with the outer coupling level to investigate the convergence performance. In this paper we extend the basic SFI method with several nonlinear acceleration techniques for improving the outer-loop convergence. Specifically, we consider numerical relaxation, quasi-Newton (QN) and Anderson acceleration (AA) methods. The acceleration techniques are adapted and studied for the first time within the context of SFI for coupled flow and transport in porous media. We reveal that the iterative form of SFI is equivalent to a nonlinear block Gauss–Seidel (BGS) process.The effectiveness of the acceleration techniques is demonstrated using several challenging examples. The results show that the basic SFI method is quite inefficient, suffering from slow convergence or even convergence failure. In order to better understand the behaviors of SFI, we carry out detailed analysis on the coupling mechanisms between the sub-problems. Compared with the basic SFI method, superior convergence performance is achieved by the acceleration techniques, which can resolve the convergence difficulties associated with various types of coupling effects. We show across a wide range of flow conditions that the acceleration techniques can stabilize the iterative process, and largely reduce the outer iteration count.

Multiscale finite volume method for finite-volume-based simulation of poroelasticity

We propose a multiscale finite volume method (MSFV) for simulation of coupled flow-deformation in heterogeneous porous media under elastic deformation (i.e., poroelastic model). The fine-scale fully resolved system of equations is obtained based on a conservative finite-volume method in which the displacement and pore pressure unknowns are located in a staggered configuration. The coupling is treated through a fully-coupled fully-implicit formulation. On this fully-coupled finite-volume system, coarse-scale grids for flow and deformation are imposed. Local basis functions for scalar pore pressure and vectorial displacement unknowns are then solved over their respective local domains at the beginning of the simulation, and reused for the rest of the time-dependent simulations. These local basis functions are then clustered to form the prolongation operator. As for the finite-volume nature of the proposed multiscale system, finite-volume restriction operators for poroelastic systems are utilised. Once the coarse-scale system is solved, its solution is prolonged back to the original fine-scale resolution, providing approximate fine-scale solution. The finite-volume multiscale formulation provides conservative stress and mass flux both at fine and coarse scale. Several numerical test cases are provided first to validate the fine-scale finite-volume discrete fully-implicit simulation, and then to investigate the accuracy of the proposed multiscale formulation. Moreover we compare our fully implicit MSFV method with hybrid multiscale Finite Element-Finite Volume (h-MSFE-FV). Our multiscale method allows for quantification of the elastic geomechanical behaviour with using only a fraction of the fine-scale grid cells, even for highly heterogeneous time-dependent models. As such, it casts a promising approach for field-scale quantification of the mechanical deformation and stress field due to injection and production in a subsurface formation.