Colebrook Equation Help with matrix errors

33 visualizaciones (últimos 30 días)
Mac
Mac el 14 de Abr. de 2016
Respondida: Walter Roberson el 15 de Abr. de 2016
I'm working on solving the colebrook equation for "f", given specific parameters, and plotting Re vs f. However I'm having trouble with matrix dimensions and running fsolve. I'm given Re is a range from 1000 to 10 million and E is a set of values 0,0.001,0.002,0.005,0.010,0.020,0.050. Can anyone help?
Here is my colebrook function:
function F = colebrook(f,Re,E)
F = 1./sqrt(f)+4*log((E./3.7)+(1.256./Re*sqrt(f)));
end
And my main script:
options = optimoptions('fsolve','Display','none','PlotFcn',@colebrook);
fun = @colebrook;
Re = 1000:10000000;
E = [0 0.001 0.002 0.005 0.010 0.020 0.050];
f = fsolve(@(f)colebrook(f,Re,E),0.1,options);
And my most recent error code:
Error using +
Matrix dimensions must agree.
Error in colebrook (line 3)
F = 1./sqrt(f)+4*log((E./3.7)+(1.256./Re*sqrt(f)));
Error in @(f)colebrook(f,Re,E)
Error in fsolve (line 219)
fuser = feval(funfcn{3},x,varargin{:});
Error in solution (line 5)
f = fsolve(@(f)colebrook(f,Re,E),0.1,options);
Caused by:
Failure in initial user-supplied objective function evaluation. FSOLVE cannot continue.

Respuestas (1)

Walter Roberson
Walter Roberson el 15 de Abr. de 2016
F = 1./sqrt(f)+4*log((E./3.7)+(1.256./Re*sqrt(f)));
Your initial conditions are scalar, so your f will be scalar, so 1./sqrt(f) will be scalar
Your E is a vector of length 7, so 4*log((E./3.7) will be a vector of length 7.
Your Re is a vector of length 9999000, so (1.256./Re*sqrt(f))) will be a vector of length 9999000 .
You cannot add a vector of length 7 to a vector of length 9999000.
If you are trying to solve over every combination of Re and E then you will need to do something like:
[gE, gRe] = ndgrid(E, Re);
F = arrayfun(@(G, RE) fsolve(@(f) colebrook(f, G, RE), 0.1), gE, gRe);
Unfortunately, the function is nonlinear, and needs to be solved individually.
You should consider, though, that the solution increases with Re and decreases with E, so you might be more efficient looping, moving to increasing value of one parameter and using the previous result as a guess.

Categorías

Más información sobre Surface and Mesh Plots en Help Center y File Exchange.

Community Treasure Hunt

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

Start Hunting!

Translated by