Facing problem while implementing dual porosity

[Shadab Hasan]

For my project work, I have been trying to make a CBM reservoir having dual porosity but I am unable to make it even after using ChatGPT if anyone can help it will be very helpful for me as I am new to mrst there might be more error in my code if possible debugging it might be useful for me

%% Jharia Coalfield Reservoir Simulation using MRST

mrstModule add ad-blackoil ad-props dual-porosity

%% Grid Setup

dims = [20, 20, 1];

physDims = [304.8, 304.8, 3.35]*meter;

G = cartGrid(dims, physDims);

G = computeGeometry(G);

%% Rock Properties (Dual Porosity for Jharia)

% From Table 2 in Asif et al. (2022):

% Matrix porosity = 0.745% and Cleat (fracture) porosity = 0.18%

matrixPorosity = 0.00745; % 0.745%

fracturePorosity = 0.0018; % 0.18%

% Permeabilities:

matrixPerm = 1.5 * milli*darcy; % Matrix permeability ~1.5 mD

fracturePerm = 6.09 * milli*darcy; % Cleat permeability ~6.09 mD

rock = dualPorosityRock(G, matrixPorosity, fracturePorosity, matrixPerm, fracturePerm);

%% Fluid Properties (CO2-ECBM System)

% Define a water-gas system; the gas phase represents CH4 and CO2.

% Fluid densities are set as typical values.

fluid = initSimpleADIFluid('phases', 'WG', ... % Water - Gas system

'n', [2 2], ... % Corey exponents

'rho', [1000 1.98]*kilogram/meter^3);

% Extended Langmuir Adsorption parameters (using values from Table 2):

% For CH4: VL = 1428.57 scf/ton, PL = 961.43 psi

% For CO2: VL = 3142.85 scf/ton, PL = 1282.4 psi

% (Ensure units are consistent with the simulation framework)

langmuirCH4 = struct('VL', 1428.57, 'PL', 961.43*psia);

langmuirCO2 = struct('VL', 3142.85, 'PL', 1282.4*psia);

fluid.adsorption = @(p, z, flag) extendedLangmuirAdsorption(p, z, langmuirCH4, langmuirCO2, flag);

%% Simulation Model

% Adjust sorption time and compressibility parameters if needed.

model = DualPorosityModel(G, rock, fluid, ...

'diffusion', true, ... % Enable Fickian diffusion

'sorptionTime', 10*day, ... % Sorption time scale (adjust as needed)

'rockCompressibility', 1.213e-3/psia, ... % Using pore compressibility from Table 2

'matrixShrinkage', 5e-6/psia); % Retain or update based on further calibration

%% Initial State

% The study reports a reservoir pressure of ~135 psi in Jharia.

pressure = 135*psia;

% Use an initial saturation guess (e.g., 50% water, 50% gas)

state = initResSol(G, pressure, [0.5 0.5]);

%% Wells Setup

% Well placement mimics the paper’s geometry: a central injector and a production well at a corner.

% CO2 Injector: set a rate based on cumulative injection of ~2.25 BSCF over 3650 days.

% Convert ~2.25e9 scf over 3650 days ≈ 17450 m^3/day (using 1 m^3 ≈ 35.3147 scf)

injRate = 17450*meter^3/day;

% Place the injector in the center cell:

centerCell = round(G.cells.num/2);

W = addWell([], G, rock, centerCell, ...

'type', 'rate', 'val', injRate, ...

'comp_i', [0 1], ... % CO2 injection (assign all gas phase to CO2)

'sign', 1);

% CH4 Producer: set bottom-hole pressure to ~75 psi (as per abandoned pressure in Table 6).

prodBHP = 75*psia;

% Place the producer at one of the corner cells (e.g., cell 1)

W = addWell(W, G, rock, 1, ...

'type', 'bhp', 'val', prodBHP, ...

'comp_i', [1 0], ... % Production of CH4 (assign gas phase to CH4)

'sign', -1);

%% Simulation Parameters

dt = 30*day; % Time step of 30 days

totalTime = 3650*day; % Simulation for 10 years

%% Run Simulation

[wellSols, states, schedule] = simulateScheduleAD(state, model, simpleSchedule(dt, 'stepCount', totalTime/dt));

%% Post-processing

% Plot production curves (e.g., gas production over time)

plotWellSols(wellSols, cumsum(schedule.step.val));