Need some guidance to solve this problem for velocity and pressure variables

3 次查看（过去 30 天）
Kumaresh Kumaresh
Kumaresh Kumaresh 2022-8-29
评论： Kumaresh Kumaresh 2022-8-29
Ran in:
Hello all,
In my problem, I need to solve the following equations to calculate the variables velocity (from Eq .(1 and 2)) and pressure (from Eq .(3)).
∇.G=S (1)
∇×G-G×∇(ln(F1+F2*G))=0 (2)
∇^2 (P^2)=-(F1+F2*G)S - G∙∇(F1+F2*G) (3)
where, G, G = mass velocity vector and its magnitude, P = Pressure and S = source term.
Here is my code below and queries to follow:
1) While defining continuity Eq. (1), I'm facing an error,
% ERROR: Array indices must be positive integers or logical values.
2) Without Eq. (1), the velocity and pressure values obtained are not good as expected.
I'm not sure about my approach as how to find the appropriate velocity and pressure.
Someone kindly examine the code and judge my approach and assist me further.
Thank you
clear; clc; close all;
L=0.25; W =5;
Nx=10; Ny=200;
x = linspace(0,L,Nx);
y = linspace(0,W,Ny);
[X, Y] = meshgrid(x,y);
dx = x/Nx; dy = y/Ny;
rhoo = 0.8;
por = 0.3; mu =0.8E-5; dia = 0.003;
P = 101325.*ones(Ny,Nx);
V = zeros(Ny,Nx);
Gvec = zeros(Ny,Nx);
TimeStep =1;
Vx = gradient(X, TimeStep);
Vy = gradient(Y, TimeStep);
Vmag = sqrt(Vx.^2 + Vy.^2);
% Boundary conditions
V(:,1) = 0;
V(Ny,:) = 0;
V(1,:) = V(Ny-1,:);
V(:,Nx) = V(:,Nx-1);
V(Ny,1) = 0;
V(1,1) = 0;
S = [0 0 0.1458 0.2058 0.2161 0.2171 0 0 0 0];
f1 = (300*((1-por).^2)*P*mu)/(rhoo.*(dia.^2)*(por.^3));
f2 = (3.5.*(1-por).*P)/(rhoo.*dia.*(por.^3));
F1 = f1.*ones(Ny,Nx);
F2 = f2.*ones(Ny,Nx);
Gvec =rhoo.*V;
Gmag =rhoo.*Vmag;
%divergence(Gvec) = 0 % Eq. (1)
% ERROR: Array indices must be positive integers or logical values.
d1 = log (F1 + (F2.*Gmag));
d2 = gradient(d1);
A = curl(Gvec,X);
B = curl(d2,Gvec);
V = A - B % Eq. (2)
V = 200×10
0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139 0.0139
d3 = F1 + (F2.*Gmag);
d4 = gradient(d3);
P = gradient(gradient((P.^2))) + divergence(d3,Gvec) + d3.*S % Eq. (3)
P = 200×10
1.0e+08 * 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000 0 0 1.0178 1.4367 1.5086 1.5155 0 0 0.0000 0.0000
  4 个评论
显示 2更早的评论
Torsten
Torsten 2022-8-29
编辑：Torsten 2022-8-29
May I know why you are guiding to use black box tools like ANSYS or COSMOL ?
The laminar Navier Stokes equations need a special treatment because of the velocity-pressure coupling. They are hard to solve for a beginner. Or how much time do you plan to invest in reading before you start programming ?
Is that impossible to solve by MATLAB ?
If it can be done in OPENFOAM, ANSYS and COMSOL, it can also be done in MATLAB. But you will have to program it from scratch if you don't find code somewhere else. So again the question: how much time do you plan to invest ? I would estimate half a year if you have to do all by yourself.
Kumaresh Kumaresh
Kumaresh Kumaresh 2022-8-29
Before working on this problem in MATLAB, I tried implementing it in OpenFOAM but failed to achieve the expected outcome. As you mentioned, the above equations are coupled sets of equations, which need special treatment. I don't find code anywhere or any default solvers related to it. This is a work by Ergun known by his name as the Ergun equation to examine the flow distribution in a porous medium.
Yes, I'm a beginner in MATLAB and this is part of my thesis work. I don't have too much time to invest in it. However, I will do my best in all possible ways to make it possible.
Kindly guide me.
Thank you for your valuable comments.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Fluid Dynamics 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by