CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   Main CFD Forum (https://www.cfd-online.com/Forums/main/)
-   -   Using the Simple Algorithm for 2D Staggered Grid in Matlab (https://www.cfd-online.com/Forums/main/82954-using-simple-algorithm-2d-staggered-grid-matlab.html)

jal121387 December 9, 2010 15:29

Using the Simple Algorithm for 2D Staggered Grid in Matlab
 
I am trying to write code for a project in Matlab. It is a pipe flow problem using the simple algorithm for a 2D staggered grid in matlab. I am quite new to the programming world and have sort of been tossed into it quickly in grad school. I need some help coding, so if anyone can offer it please let me know.

banmakadachha April 1, 2014 14:00

have you figured this out yet
 
If you have figured this out, could you please send me the code at banmakadachha@gmail.com

tas38 April 1, 2014 15:56

You may be able to modify the staggered lid driven cavity code available in the MATLAB file exchange to suit your needs:

http://www.mathworks.com/matlabcentr...en-cavity-flow

engrmansoor2534 April 9, 2014 15:41

Matlab code for SIMPLEC Algortihm
 
I am trying to write Matlab code for one dimensional convergent nozzle having 5 pressure pressure nodes and 4 velocity nodes using SIMPLEC algorithm.I have done the necessary changes required in SIMPLE i.e changed the value of d term for both momentum equation and pressure correction equations. my problem is that my code is not converging it is giving NaN error.I am using an under relaxation of 0.8. SIMPLE code is running properly and converges after 18 iterations. I have kept the criteria of convergence the residual for the momentum equations.some help required....
thanks

tas38 April 9, 2014 16:55

Can you please post both the SIMPLE and SIMPLEC versions of the codes? As they are one-dimensional, it should be relatively easy to find the error.

engrmansoor2534 April 10, 2014 02:16

Matlab code for SIMPLEC Algortihm
 
here are the values for d term for SIMPLEC algorithm and the corresponding residuals.
D_2=((A_B+A_C)/(2))/((AP(2)-(AP(3)+AP(1))));
D_1=((A_B+A_A)/(2))/((AP(1)-AP(2)));
and similar for D_3, D_4
where AP(1) and AP(2) and so on are the coefficients of velocities
As in SIMPLE d is only equal to
D_2=((A_B+A_C)/(2*AP(2)));

I think the error is in TDMA solution of pressure correction equation
its code
A_W_B=(density*D_1*((A_B+A_A)/2));
A_E_B=(density*D_2*((A_B+A_C)/2));

F_W_B=(density*U1_STAR*((A_B+A_A)/2));
F_E_B=(density*U2_STAR*((A_B+A_C)/2));

A_P_B=A_W_B+A_E_B;
b_B=F_W_B-F_E_B;


A_W_C=(density*D_2*((A_B+A_C)/2));
A_E_C=(density*D_3*((A_D+A_C)/2));

F_W_C=(density*U2_STAR*((A_B+A_C)/2));
F_E_C=(density*U3_STAR*((A_D+A_C)/2));

A_P_C=A_W_C+A_E_C;
b_C=F_W_C-F_E_C;

A_W_D=(density*D_3*((A_D+A_C)/2));
A_E_D=(density*D_4*((A_D+A_E)/2));

F_W_D=(density*U3_STAR*((A_C+A_D)/2));
F_E_D=(density*U4_STAR*((A_D+A_E)/2));

A_P_D=A_W_D+A_E_D;
b_D=F_W_D-F_E_D;

P_DASH_A=0;
P_DASH_E=0;

A=zeros(3);
A(1,1)=A_P_B;
A(1,2)=-A_E_B;
A(2,1)=-A_W_C;
A(2,2)=A_P_C;
A(2,3)=-A_E_C;
A(3,2)=-A_W_D;
A(3,3)=A_P_D;
b(1)=b_B;
b(2)=b_C;
b(3)=b_D;
m=3;
X=zeros(m,1);
A(1,2)= A(1,2) / A(1,1);
b(1)= b(1) / A(1,1);
for i=2:m-1

press= A(i,i) - (A(i,i-1) * A(i-1,i));

A(i,i+1)= A(i,i+1) / press;

b(i)= ( b(i) - A(i,i-1) * b(i-1) ) / press;
end

i=m;

X(m)=(b(i) - (A(i,i-1) * b(i-1))) / (A(i,i) - (A(i,i-1) * A(i-1,i)));

for a=m-1:-1:1

X(a)= b(a)-(A(a,a+1) * X(a+1));

end

and the residual I have set for convergence is
r(1)=(AP(1)*FU(1))-SU(1);
r(2)=(AP(2)*FU(2))-((F_w_2*FU(1))+SU(2));
r(3)=(AP(3)*FU(3))-((F_w_3*FU(2))+SU(3));
r(4)=(AP(4)*FU(4))-((F_w_4*FU(3))+SU(4));
res= abs(r(1)+r(2)+r(3)+r(4));
if (res<=0.00001)
break
end
kindly check if there are any errors in my code
thanks

ahr February 18, 2019 13:56

HI,

Im writting the same problem in fortran 90; I got the same (under-relaxed) results until the first iteration for 5 nodes, it means:

u1=1.78856
u2=2.29959
u3=3.21942
u4=5.36571

and also I got:

pA=9.08536
pB=8.8115
pC=8.3396
pD=7.4664
pE=0

the issue is that once I run the next iteration, I dont get the mentioned equation in the book, at node 1: 1.20425(U1)=1.98592, I get: 1.14848(U1)=1.412218

I did it by hand, I dont know what I am doing wrong.

duri February 19, 2019 04:38

Attach a debugger check along with your hand calculation.


All times are GMT -4. The time now is 00:02.