NOTE: You must submit the listing of any function that you write and a log of yo

Question

NOTE: You must submit the listing of any function that you write and a log of your session in MATLAB with publish feature. 1. Colebrook equation and plot The frictional loses in pipe flows are computed using the Colebrook equation given below. It is also often presented in the form of a graph called Moody Chart. Your task is to write a MATLAB function and scripts to generate a plot of f us Re for various values of Take the range of Re E [ 1000, 10000000 and E.-|0 0.001 0.002 0.005 0.010 0.020 0.050] The lines for every value of should be on the same graph. The Colebrook equation is given below. f=|-410g, le/D) , 1.256]1-2 -4log10 3.7 Re 3.7 Revf fsolve Note that f appears implicitly in this equation and so you must use fsolve to solve for various values of f for a give set of Re,e/D. The graph should appear as shown below with all of the labels and titles as shown in the graphs. In doing this you will have to learn about the following commands on your own: xlim, xlabel, ylabel, title, grid, text.

Explanation / Answer

We will need two matlab functions .

First one(moody) to find the values of f and second one(MyMoody) to plot the graph.

Codes are given below:

function f = moody(ed,Re,verbose)
% moody Find friction factor by solving the Colebrook equation (Moody Chart)
%
% Synopsis: f = moody(ed,Re)
%
% Input: ed = relative roughness = epsilon/diameter
% Re = Reynolds number
%
% Output: f = friction factor
%
% Note: Laminar and turbulent flow are correctly accounted for
if Re<0
error(sprintf(’Reynolds number = %f cannot be negative’,Re));
elseif Re<2000
f = 64/Re; return % laminar flow
end
if ed>0.05
warning(sprintf(’epsilon/diameter ratio = %f is not on Moody chart’,ed));
end
if Re<4000, warning(’Re = %f in transition range’,Re); end
% --- Use fzero to find f from Colebrook equation.
% coleFun is an inline function object to evaluate F(f,e/d,Re)
% fzero returns the value of f such that F(f,e/d/Re) = 0 (approximately)
% fi = initial guess from Haaland equation, see White, equation 6.64a
% Iterations of fzero are terminated when f is known to whithin +/- dfTol
coleFun = inline(’(-4.0*log10( ed/3.7 + 1.256/( Re*sqrt(f)))^(-2) )’,...
’f’,’ed’,’Re’);
fi = 1/(1.8*log10(6.9/Re + (ed/3.7)^1.11))^2; % initial guess at f
dfTol = 5e-6;
f = fzero(coleFun,fi,optimset(’TolX’,dfTol,’Display’,’off’),ed,Re);
% --- sanity check:
if f<0, error(sprintf(’Friction factor = %f, but cannot be negative’,f)); end

------------------------------------SECOND FUNCTION STARTS HERE--------------------------------------------------------------------------

function myMoody
% myMoody Make a simple Moody chart
% --- Generate a log-spaced vector of Re values in the range 1000 <= Re < 10^8
Re = logspace(log10(1000),8,50);
ed = [0 0.001 0.002 0.005 0.010 0.020 0.050];
f = zeros(size(Re));
% --- Plot f(Re) curves for one value of epsilon/D at a time
% Temporarily turn warnings off to avoid lots of messages when
% 2000 < Re < 4000
warning(’off’)
figm = figure; hold(’on’);
for i=1:length(ed)
for j=1:length(Re)
f(j) = moody(ed(i),Re(j));
end
loglog(Re,f,’k-’);
end
ReLam = [100 2000];
fLam = 64./ReLam;
loglog(ReLam,fLam,’r-’);
xlabel(’Re’); ylabel(’f’,’Rotation’,0)
axis([100 1e9 5e-3 2e-1]);
hold(’off’);
grid(’on’)
warning(’on’);
% --- MATLAB resets loglog scale? Set it back
set(gca,’Xscale’,’log’,’Yscale’,’log’);
% --- Add text labels like "epsilon/D = ..." at right end of the plot
Remax = max(Re);
ReLabel = 10^( floor(log10(Remax)) - 1);
for i=2:length(ed)
fLabel = moody(ed(i),Remax);
text(ReLabel,1.1*fLabel,sprintf(’\epsilon/D = %5.4f’,ed(i)),’FontSize’,14);
end

Related Questions

Navigate

• Browse (All)
• Browse N
• Subjects

---

---

Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.