
[Sponsors] 
May 5, 2020, 10:29 
Flat Plate Boundary Layer

#1 
New Member
Colin
Join Date: Apr 2020
Posts: 7
Rep Power: 2 
I want to solve Continuity and xmomentum Equations for a boundary layer over a flat plate. It is a 0 pressure gradient flow.
My equations are diverging, and I can't seem to find why they are. I am using an underrelaxation technique to obtain convergence. Could anyone help to spot my mistake? Using a central/forward/backward difference scheme based on where I am at on my grid. MATLAB Code is below: %% Parameters rho = 1.225; % kg/m^3  Air sea level standard day U = 10; % m/s  freestream mu = 1.789e5; % kg/(m*s)  Air sea level standard day nu = mu/rho; % m^2/s L = 0.25; % m H = 0.05; % m % nx and ny defines matrix grid, dx and dy for step size nx = 51; ny = 11; dx = .005; % 50 grid points in xdir  .005m spacing dy = .005; % 10 grid points in ydir  .005m spacing % Re based on length Re_L = U*L/nu; % Filling mesh with 0 to begin for v and 10 for all u v = zeros(ny,nx); % All initial v are 0 u = ones(ny,nx)*U; % Filling the entire mesh with all u = 10 % Boundary Conditions % Goes (j,i) each segment, fills mesh with known values u(:,1) = U; % Inlet v(:,1) = 0; % Inlet u(1, = 0; % No Slip wall v(1, = 0; % No Slip wall u(ny, = U; % Free stream (top of mesh) %% Solving % Guess values of alpha alphau = .05; %guess, will need to figure out best value for this alphav = .015; % Differences diff_u=[]; diff_v=[]; diff_u_max=[]; diff_v_max=[]; diff_both=[]; % Tolerance and convergence for main while loop tol = .001; z = 0; % Starting the index for iterations while z < 100 % This is number of iterations it performs for i = 2:1:nx1 for j = 2:1:ny1 if i == nx % use backwards v_old = v(j,i); %v(j,i) = ((u(j,i1)  u(j,i))/dx)*dy + (v(j1,i)); v(j,i) = ((u(j,i1)  u(j,i))/dx)*dy + (v(j1,i)); change_v = v_old  v(j,i); v(j,i)= v_old + (alphav * change_v); u_old = u(j,i); u(j,i) = u(j,i1) + (((v(j,i) + v(j1,i)) / dy) * dx); change_u = u_old  u(j,i+1); u(j,i+1) = u_old + (alphau * change_u); diff_u = [diff_u change_u]; diff_v = [diff_v change_v]; elseif j == ny v_old = v(j,i); v(j,i) = v(j1,i); u_old = u(j,i+1); u(j,i+1) = U; else % This is central differences u_old = u(j,i+1); u(j,i+1) = (2*dx/u(j,i)) * (((mu/rho)*((u(j1,i) + 2*u(j,i) + u(j+1,i))/(dy^2)))  (v(j,i)*((u(j1,i)  u(j+1,i))/(2*dy)))) + u(j,i1); change_u = u_old  u(j,i+1); u(j,i+1) = u_old + (alphau * change_u); v_old = v(j+1,i); v(j+1,i) = v(j1,i+1)  dy/2/dx*(u(j,i+1)u(j,i)+u(j1,i+1)u(j1,i)); change_v = v_old  v(j+1,i); v(j+1,i)= v_old + (alphav * change_v); diff_u = [diff_u change_u]; diff_v = [diff_v change_v]; end end end diff_U = max(abs(diff_u)); diff_V = max(abs(diff_v)); diff_u_max = [diff_u_max diff_U]; diff_v_max = [diff_v_max diff_V]; diff = [diff_U diff_V]; max_diff = max(diff); if abs(max_diff) < tol break end z = z + 1; end % Plotting to show convergence throughout the iterations figure hold on plot(diff_u_max) plot(diff_v_max) hold off xlabel('Iterations'); ylabel('Chnage in u and v'); title('Convergence of u and v Throughout Iterations'); legend('u difference', 'v difference'); 

May 5, 2020, 11:08 

#2 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 5,665
Rep Power: 60 
The Prandtl equations are parabolic, you need to integrate along x in a timelike fashion


May 5, 2020, 11:21 

#3 
New Member
Colin
Join Date: Apr 2020
Posts: 7
Rep Power: 2 

May 5, 2020, 11:28 

#4  
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 5,665
Rep Power: 60 
Quote:
if you have u(i,j) standing for u(x,y), you compute u(i+1,j) basing on the values at i. Thus, x is the timelike direction. Be aware of the numerical stability constraint for the parabolic equation. 

May 5, 2020, 11:44 

#5  
New Member
Colin
Join Date: Apr 2020
Posts: 7
Rep Power: 2 
Quote:
I don't normally do anything with CFD, but this is for a class project in a fluids class. I am satisfying the stability criteria, since my step size is .005 in both dx and dy, and then my underrelaxation technique takes care of the numerical stability by new = old + alpha*(old  calculated). So by going u(i+1,j) I need to change my for loop from i = 2:1:nx1 to 2:1:? if it is to be time stepped? Thanks [/B] 

May 5, 2020, 12:01 

#6  
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 5,665
Rep Power: 60 
Quote:
The problem is solved in an open domain, you fix the initial condition at x=0 and the BCs at y=0 and y>+Inf. You have to proceed from a position x to a position x+dx by integrating along the normal to wall direction y. There is no theoretical ending in this process, you integrate from x=0 to the position x required by your problem where you get the Re_x value. 

May 5, 2020, 12:10 

#7  
New Member
Colin
Join Date: Apr 2020
Posts: 7
Rep Power: 2 
Quote:
I am going the entire length of the plate to look at the growth of the boundary layer. from x = 0 to x = .25m, where i is each grid point I fill with a u and v velocity value I wish I could put a photo in of my grid to show 

Tags 
cfd, convergence, flat plate 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
LES outflow condition for a subsonic flat plate boundary layer  chris x  Main CFD Forum  3  August 25, 2018 05:09 
Compressible Flat Plate with coupled hydrodynamic/thermal boundary layer  Obad  OpenFOAM PreProcessing  0  September 14, 2016 14:45 
Wrong flow in ratating domain problem  Sanyo  CFX  17  August 15, 2015 06:20 
Bump in flat plate boundary layer  DCK  STARCCM+  2  July 23, 2015 01:55 
Flat plate boundary layer problem  student  Main CFD Forum  3  May 21, 2007 13:10 