Manipulating an input variable for newton raphson numerical method

Hi

My code below solves colebrook eqn for friction factor using newton raphson method. The issue is that I dont know how to manipulate an input variable (e.g Rynolds no.) to see the resulting f factor ( x in this code）.

My function;

function nr(y,x0,n,varargin)
x = x0
%convert symbolic function to function_handle to find derivitve 
f = matlabFunction(y) 
y1= diff(y)
f1= matlabFunction(y1)
for c = 1:n 
x = x - f(x)/f1(x)
k = [x f(x)]; 
disp(k)
end

my script

format long
Re = 10000;
y = (x.^(-1/2))+(0.86*log(((0.0001)/3.7)+((2.51/Re)*x.^(-1/2))));
nr(y,0.01,10);

If for example I input Re = linspace(10000,60000,10) I get an error "matrix dimension must agree"

Additionally, I would like to make a matrix with each Re and final friction factor value (x)

Current output

>> samplescribt

x =

     0.010000000000000

f =

    function_handle with value:

      @(x)log(1.0./sqrt(x).*2.51e-4+2.702702702702703e-5).*(4.3e1./5.0e1)+1.0./sqrt(x)

y1 =

- 1/(2*x^(3/2)) - 199095708787547171/(1844674407370955161600*x^(3/2)*(4630132762501097/(18446744073709551616*x^(1/2)) + 1/37000))

f1 =

    function_handle with value:

      @(x)1.0./x.^(3.0./2.0).*(-1.0./2.0)-(1.0./x.^(3.0./2.0).*1.0793e-4)./(1.0./sqrt(x).*2.51e-4+2.702702702702703e-5)

x =

     0.018957804377517

     0.018957804377517   1.851220793655044

x =

     0.027612285982097

     0.027612285982097   0.447238114020280

x =

     0.031211167390416

     0.031211167390416   0.037929950196298

x =

     0.031575183433487

     0.031575183433487   0.000314699305965

x =

     0.031578254339019

     0.031578254339019   0.000000021974257

x =

     0.031578254553479

     0.031578254553479   0.000000000000002

x =

     0.031578254553479

     0.031578254553479                   0

x =

     0.031578254553479

     0.031578254553479                   0

x =

     0.031578254553479

     0.031578254553479                   0

x =

     0.031578254553479

     0.031578254553479                   0

>>

Thanks

Answer 1

[Darshan Ramakant Bhat]

I see that in the equation

y = (x.^(-1/2))+(0.86*log(((0.0001)/3.7)+((2.51/Re)*x.^(-1/2))));

you are dividing a scalar value by "Re". So when you make "Re = Re = linspace(10000,60000,10) " it will become a vector. You cannot divide a scalar by a vector, that is the reason you are getting the error. To manipulate the "Re" value one by one, use a loop. Below is the sample code:

format long
Re = linspace(10000,60000,10)
syms x;
for i =1:length(Re)
y = (x.^(-1/2))+(0.86*log(((0.0001)/3.7)+((2.51/Re(i))*x.^(-1/2))));
nr(y,0.01,10);
end

[Abdullah Althobaity]

Thanks Darshan, that solved the problem.

Regards