How do I complete my code to plot the Moody Chart?
26 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
function [f] = frictionFactor(Re, ed) %Re = Reynolds Number, ed = eps/d, relative roughness
colebrook = @(f) 1/sqrt(f)+2*log10((ed/3.7)+(2.51)/(Re*sqrt(f)));
if Re > 4000 %turbulent
f = fzero(colebrook, [0.008, 0.1]);
elseif Re < 2000 %laminar
f = 64/Re;
else %transitional
f = (((Re-2000)/(4000-2000))*(0.1-0.008))+0.008;
end
array = linspace(0.000001, 0.05, 21);
for i = array
colebrook(i)
end
end
2 Kommentare
Walter Roberson
am 1 Mär. 2018
It is not clear to me why you calculate f and then ignore it when you for i = array ?
I somehow suspect that the values in array are intended to represent different Re values that you want to evaluate colebrook with after figuring out what f value you want to use ?
Antworten (2)
Walter Roberson
am 1 Mär. 2018
function [f, rough] = frictionFactor(Re, ed) %Re = Reynolds Number, ed = eps/d, relative roughness
colebrook_fed = @(f, ed) 1/sqrt(f)+2*log10((ed/3.7)+(2.51)/(Re*sqrt(f)));
if Re > 4000 %turbulent
f = fzero(@(f) colebrook_fed(f, ed), [0.008, 0.1]);
elseif Re < 2000 %laminar
f = 64/Re;
else %transitional
f = (((Re-2000)/(4000-2000))*(0.1-0.008))+0.008;
end
array = linspace(0.000001, 0.05, 21);
narray = length(array);
rough = zeros(1, narray);
for K = 1 : narray
i = array(K);
rough(K) = colebrook_fed(f, i);
end
0 Kommentare
Nikolaj Maack Bielefeld
am 28 Mär. 2020
Bearbeitet: Nikolaj Maack Bielefeld
am 28 Mär. 2020
You could also have used the symbolic math implicit plot command:
close all; clear; clc
% symbolic math: y = f, x = Re
syms x y
% relative roughness
relrough = [0 2e-7 1e-6 5e-6 10e-6 50e-6 100e-6 200e-6 400e-6 600e-6 ...
800e-6 1e-3 2e-3 4e-3 6e-3 8e-3 10e-3 15e-3 20e-3 30e-3 40e-3 50e-3];
% Colebrook equation
eqn = 1/sqrt(y) == -2*log10(relrough/3.7+2.51/(x*sqrt(y)));
% implicit plot with x- and y-limits
fimplicit(eqn,[2000 10e8 0.006 0.1])
set(gca, 'XScale', 'log') % logarithmic x-axis
set(gca, 'YScale', 'log') % logarithmic y-axis
title('Moody Chart')
xlabel('Reynolds Number')
ylabel('Friction Factor')
2 Kommentare
Walter Roberson
am 28 Mär. 2020
Because of the three different ranges of values, your lower bound on x for the fimplicit should be 4000. You would need to also hold on and plot the other two parts (you could plot both together if you used piecewise)
Nikolaj Maack Bielefeld
am 29 Mär. 2020
Yes, I agree Walter.
I just posted such a plot in this thread:
Siehe auch
Kategorien
Mehr zu Line Plots finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!