MATLAB Reservoir Simulation Toolbox

**Mataleb** · 25/08/2021, 12h27

Bonjour tout le monde, je suis débutant en Matlab, je n'arrive pas à exécuter cet exemple disponible sur le site: https://www.sintef.no/projectweb/mrs...tion-tutorial/
Merci d'avance.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
%Animated example
%We will now simulate a simple transport problem between the wells. We will plot the empty grid and the wells, and in each time step we will visualize the saturation in the cells where the inejcted phase is above a threshold, here 0.05.
% Disable gravity
%gravity off
% Set up uniform permeability and constant porosity
G = cartGrid([10, 10, 3]);
G = computeGeometry(G);
rock.perm = repmat(100*milli*darcy, [G.cells.num, 1]);
rock.poro = repmat(0.5            , [G.cells.num, 1]);
% A simple two phase system.
fluid = initSimpleFluid("mu" , [   1,  10]*centi*poise     , ...
                            "rho", [1014, 859]*kilogram/meter^3, ...
                            "n"  , [   2,   2]);
% Two wells, an injector at 1 bar and a producer at 0 bar.
W = verticalWell([], G, rock, 1, 1, [], ...
            'Type', 'bhp' , 'Val', 1*barsa(), ...
            'Radius', 0.1, 'InnerProduct', 'ip_tpf', ...
            'Comp_i', [0, 1]);
 
W = verticalWell(W, G, rock, 10, 10, [], ...
            "Type", "bhp" , "Val", 0*barsa(), ...
            "Radius", 0.1, "InnerProduct", "ip_tpf", ...
            "Comp_i", [1, 0]);
 
% Create a initialized state and set initial saturation to phase 1.
sol = initState(G, [], 0, [1, 0]);
 
% Find transmissibility.
T = computeTrans(G, rock);
dT = 10*day;
% Reference TPFA
psolve = @(state) incompTPFA(state, G, T, fluid, "wells", W);
%psolve = @(state) incompTPFA(state, G, T, fluid, "wellSol", W);
 
% Implicit transport solver
tsolve = @(state, dT) implicitTransport(state, G, dT, rock, ...
                                 %              fluid, "wells", W);
                                               fluid, "wellSol", W);
for i = 1:20
%    sol = tsolve(sol, dT);  %le code signale error
    sol = psolve(sol);
    clf;
    %plotGrid(G) %It grid all the volume
    plotGrid(G, "FaceAlpha", 0, "EdgeAlpha", .1)
    plotCellData(G, sol.s(:,2), sol.s(:,2)>0.05)
    %plotCellData(G, sol.pressure, j == round(G.cartDims(2)/2))
    plotWell(G, W);
    view(30,50);
    pause(.5)
end