Using EDFM with Large Grids

[Mahmoud Heshmat]

Hello,

I hope you're all doing well.

I am trying to use EDFM in a large grid (physical dimensions = (1800, 880, 1100) meters, with arbitrary grid resolution. 

The problem is the EDFMgrid(..) function gives an error that I find hard to explain. This issue only shows in large grids.

I tried changing tolerance but it still doesn't work.

I tried debugging it using help from GPT but still couldn't solve the issue. 

what I found is that the calculated volume of the grid generated doesn't match the target volume of the fracture within tolerance and the code keeps decreasing the tolerance but it fails in the end. 

I would appreciate your help.

best,

Mahmoud Abdellatif

the code:

%% LOAD NECESSARY MODULES
mrstModule add hfm;             % hybrid fracture module
mrstModule add mrst-gui;        % plotting routines
mrstModule add ad-blackoil;     % AD blackoil solver
mrstModule add ad-core;         % AD core module
mrstModule add ad-props;        % AD properties

%%
xRef = -1000; yRef = -200; refDepth = 1900;
[Lx, Ly, Lz] = deal(1800, 800, 1100);
[nx, ny, nz] = deal(1800/50, 800/50, 1100/50);

physdim = [Lx, Ly, Lz];
G = cartGrid([nx, ny, nz], physdim);

G.nodes.coords(:, 1) = G.nodes.coords(:, 1) + xRef;
G.nodes.coords(:, 2) = G.nodes.coords(:, 2) + yRef;
G.nodes.coords(:, 3) = G.nodes.coords(:, 3) + refDepth;
G = computeGeometry(G);

%%

figure;
plotGrid(G, 'facecolor', 'none', 'edgealpha', 0.4);

view(3);

%% medium properties
perm = 1e-6*darcy;
poro = 0.01;  
G.rock = makeRock(G, perm, poro);  
% G.rock = addThermalRockProps(G.rock, 'lambdaR', 3.05, 'rhoR', 2700, 'CpR', 850);

%% medium properties
fluid = initSimpleADIFluid('mu', 0.001, 'rho', 1000*kilogram/meter^3, 'phases', 'W');%, 'c', 4.4e-10, 'pRef', 0);
% fluid = addThermalFluidProps(fluid, 'Cp'     , 4.2e3, ...
%                                     'lambdaF', 0.6, ...
%                                     'useEOS' , false );

%%

fracplanes(1).points = [125, -125, 2000;
                        125, 525, 2000;
                        -625, 525, 2900;
                        -625, -125, 2900];

fracplanes(1).aperture = 0.01;
fracplanes(1).poro = 1.0;
fracplanes(1).perm = 1e1*darcy;


fracplanes(2).points = [25, -225, 1925;
                            25, 625, 1925;
                            25, 625, 1925+1000;
                            25, -225, 1925+1000];
fracplanes(2).aperture = 1/25;
fracplanes(2).poro=0.8;
fracplanes(2).perm=1e1*darcy;
%%
figure;
plotfracongrid(G,fracplanes,'label',false,...
    'fracfacealpha',1,'fracfacecolor','g');
view(30,45);
axis equal tight


%%
tol = 1e3; % Tolerance to check for plane intersection
[G,fracplanes] = EDFMgrid(G, fracplanes, 'Tolerance', tol, 'verbose', false);

The error:

Index in position 1 is invalid. Array indices must be positive integers or logical
values.

Error in checkplaneAABBintersect (line 28)
endpoints=nodes(nodepair,:);
^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Error in EDFMgrid (line 106)
        parfor j=1:nummatcells % check for plane intersection. PARFOR HERE.
        ^^^^^^^^^^^^^^^^^^^^^^^^
Error in test (line 109)
[G,fracplanes] = EDFMgrid(G, fracplanes, 'Tolerance', tol, 'verbose', false);
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^