Question

can somebody please help me with matlab, i need someone to help me fix this code.
clear all;
clc;

den_b = 8550;
den_f = 1100;
visc_f = 0.1;
D_b = 0.0064;
g = 9.81;
u = 0.1;
step = 0.01;
n = 100;
tol = 10^-4;

Re = @(u)(D_b*abs(u)*den_f)/visc_f
bouyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45
CirArea =(pi/4).*D_b^2
f =@(u) bouyancy-(Drag(u).*CirArea*den_f.*u*abs(u)/2)

for ii = 1:n
u(ii+1) = u(ii)-f(u(ii))/finitediff(u(ii)) 11 years ago

Answer 1

@Lion.k
What does the program do? You don't have a lot of documentation.

Answer 2

basically i have to use newtons method to find the new velocity.

Answer 3

Ok, so it's a numerical methods problem!

Answer 4

What seems to be bugging you?

Answer 5

well, in my command window it doesn't evaluate the fuction but instead prints the function out.

Answer 6

which function?

Answer 7

den_b = 8550; %density ball
den_f = 1100; % density of fluid
visc_f = 0.1; %visocsity
D_b = 0.0064; % diameter of bal
g = 9.81; % gravity
u = 0.1; % velocity
step = 0.01;
n = 100;
tol = 10^-4;
Re = @(u)(D_b*abs(u)*den_f)/visc_f % function of Reynold number
bouyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g %function of bouyacy
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45 %function of drag(Cd)
CirArea =(pi/4).*D_b^2 % area of circ
f =@(u) bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2) % the main function. a combnation of all.

Answer 8

Is it the same code that you had above?

Answer 9

Did you write finitediff?

Answer 10

yup but with documentation, so you'd get clear idea.
main function that consists of all is this one

Answer 11

finitediff function file to find the derivative of
f =@(u) bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)

Answer 12

Can you tell me which is your target variable (that you need to iterate using Newton's method), and how do you plan to find its derivative.

Answer 13

u

Answer 14

So you have a separate file for finitediff, or is it defined in the program above?

Answer 15

i made a function file to caluculate its derivative, since the derivative cannot be found directly. sinc the function is embedded within a function.

Answer 16

function [d]= finitediff(f,x)
h = sqrt(eps);
d = (f(x+h) - f(x))/h;
end

Answer 17

So you are trying to refine
u1=u0-f(u0)/f'(u0) ?

Answer 18

yeah.

Answer 19

So your starting value for u = 0.1?

Answer 20

f is basically :
Re = @(u)(D_b*abs(u)*den_f)/visc_f % function of Reynold number
bouyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g %function of bouyacy
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45 %function of drag(Cd)
CirArea =(pi/4).*D_b^2 % area of circ.

Answer 21

yeah

Answer 22

the answer that i'm suppose to get is -0.584

Answer 23

Perhaps I missed it, but I don't see the loop to iterate u. Can you help me find it?

Answer 24

for ii = 1:n
u(ii+1) = u(ii)-f(u(ii))/finitediff(u(ii))<tol
end

Answer 25

the thing is.. f is not evaluating.

Answer 26

do you need the values of u other than for printing?

Answer 27

yea.. to plot the graph

Answer 28

You used finitediff(...) which Matlab complains is undefined.
Did you compile it elsewhere and forgot to link it?

Answer 29

*its

Answer 30

its there.. its working for me. but f is not working.

Answer 31

you need to save the finitediff function and the script in same folder . though i think the reason its not working is because of f.

Answer 32

finitediff(@cos,0) if you try this, the function is working. its the script which is not.

Answer 33

Are you allowed to put the whole thing in one file?

Answer 34

like... self-contained?

Answer 35

Is this a course on fluid mechanics or numerical methods?

Answer 36

tbh i haven't tried. i don't know how to

Answer 37

its mathematical modelling

Answer 38

ok.
What I suggest you do, is not to start big.
Can you start with a simple Newton's method test, like finding square root, to make sure the logic is good. Once you get it going, it will be relative simple to replace the function and the derivative by something else.
This way, it will be easier for you or anybody to debug.

Answer 39

Also, with Newtons, I always print out the results of every iteration.

Answer 40

The reason for that is also Newton's problem, ... it does not always converge, and is very dependent on the initial value. But I don't think you're there yet.

Answer 41

the newtons method is not the problem. its the function.
Re = @(u)(D_b*abs(u)*den_f)/visc_f % function of Reynold number
bouyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g %function of bouyacy
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45 %function of drag(Cd)
CirArea =(pi/4).*D_b^2 % area of circ.

there is something wrong with these, and i can't figure it out.

it can't derive the function directly. i i will have to solve it using finitediff function. but the function itself is not function. how will i apply the newtons method?

Answer 42

Is u the velocity? Does it ever go negative?

Answer 43

yea it will go negative when its going down.

Answer 44

Is it a scalar or a vector?

Answer 45

vector

Answer 46

So u has two components?

Answer 47

i did newtons method example: here clear;
clc;

x0 = 3;
x = x0;
C = -1;

f = @(x) 10*x.^2 -10*x + 0.1*sinh(x) + C;%function
d = @(x) 20*x + 0.1*cosh(x)-10;%derivative
n = 5;
% newtons method
for ii = 1:n
x = x - f(x)/d(x)%calculating new root
end
i've applied the same method in my code. but in this question i can derive the function derectly but for my code, i cannot since the function is embedded within a function.

Answer 48

actually no u is scalar.

Answer 49

`plot(urray,yarray)`
should be
`plot(uarray,yarray) `

?

Answer 50

So does u ever go negative?
The reason I ask is because you have abs(u), which may give difficulties in finding derivative.

Answer 51

u will go negative.

Answer 52

So a backflow?

Answer 53

Oh yes, the answer is negative, right?

Answer 54

thats why in the equation its has abs

Answer 55

yup u tuns out to be -0.584

Answer 56

@UnkleRhaukus i know that's a typo. but thats not the problem. the funtion "f" is not evaluating.

Answer 57

when you plug in u as -0.584 , then Re is 41.1.

Answer 58

Is finitediff the code you posted later?

Answer 59

yea. function file

Answer 60

function [d]= finitediff(f,x)
h = sqrt(eps);
d = (f(x+h) - f(x))/h;
end

Answer 61

The finitediff takes two arguments
so the iterations should be

```
for ii = 1:n
u(ii+1) = u(ii)-f(u(ii))/finitediff(f,u(ii)) < tol;
end
```

Answer 62

so function f was passed as an argument, did you verify that it evaluates the f you want?

Answer 63

Actually, I don't think the loop is even in the second version of the code.

Answer 64

Can you post ALL the code in one post and show the result that you get?

Answer 65

With the finitediff function in a separate file
the rest of the code i have is :

```
%%% L0 %%%
% %
%%% %%%

%% Clear
clear all;
clc;

%% Variables
den_b = 8550; % density ball
den_f = 1100; % density of fluid
visc_f = 0.1; % viscosity
D_b = 0.0064; % diameter of bal
g = 9.81; % gravity
u = 0.1; % velocity

step = 0.01; % step size
n = 100; % number of iterations
tol = 10^-4; % total something?

%% function of Reynold number
Re = @(u)(D_b*abs(u)*den_f)/visc_f;

%% function of bouyacy
buoyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g;

%% function of drag(Cd)
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45;

%% area of circle
CirArea =(pi/4).*D_b^2;

% the main function. a combnation of all.
f =@(u) buoyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2);

u = zeros(n);
for ii = 1:n
u(ii+1) = u(ii)-f(u(ii))/finitediff(f,u(ii)) < tol;
end

uarray = 0.1:0.1:100;
yarray = f(uarray);
plot(uarray,yarray)
```

and it plots this:

Answer 66

@UnkleRhaukus Thanks!

Answer 67

Re =

@(u)(D_b*abs(u)*den_f)/visc_f


bouyancy =

-0.0100


Drag =

@(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45


CirArea =

3.2170e-05


f =

@(u)bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)

Undefined function or variable 'buoyancy'.

Error in @(u)buoyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)
i get this

Answer 68

@UnkleRhaukus the command window shows this, the graph is basically just plotting uarray against uarray.

Answer 69

1. Buoyancy was spelled wrong in the definition.
2. uarray that you plotted does not seem to be what you saved.

Answer 70

and in the definition of f.

Answer 71

i checked again: clear all;
clc;
den_b = 8550;
den_f = 1100;
visc_f = 0.1;
D_b = 0.0064;
g = 9.81;
u = 0.1;
step = 0.01;
n = 100;
tol = 10^-4;
Re = @(u)(D_b*abs(u)*den_f)/visc_f
bouyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45
CirArea =(pi/4).*D_b^2
f =@(u) bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)
%for ii = 1:n
%u(ii+1) = u(ii)-f(u(ii))/finitediff(u(ii))
%end
uarray = -100:0.1:100;
yarray = f(uarray);
plot(uarray,yarray)

Answer 72

Re = @(u)(D_b*abs(u)*den_f)/visc_f
bouyancy = (den_f-den_b)*(pi/6).*D_b.^3.*g
Drag = @(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45
CirArea =(pi/4).*D_b^2
f =@(u) bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)
these are the problem.

Answer 73

I don't understand the problem,
What isn't it doing, that you want it to do?

Answer 74

evaluate. like give a number. but instead it prints

Answer 75

prints what?

Answer 76

Can you post the output like you did earlier?

Answer 77

Re =

@(u)(D_b*abs(u)*den_f)/visc_f


bouyancy =

-0.0100


Drag =

@(u)((2.25./Re(u).^0.31)+(0.358.*Re(u).^0.06)).^3.45


CirArea =

3.2170e-05


f =

@(u)bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)

>>

Answer 78

Did you put comments on the loop?

Answer 79

loop is newtons method.
but for newtons method to work, f needs to be right. but the function is guving any output.

Answer 80

i meant function is not giving an output

Answer 81

For f(u) to execute, you need to give a value to u as well.

Answer 82

You just defined it, not executing it except in the loop.

Answer 83

yea u is 0.1 initial guess.

Answer 84

f= @(u)
i did

Answer 85

I am not familiar with Matlab, but that looks like a definition to me, not execution.
If you put f(0.1) on the following line, you should see numbers.

Answer 86

... or f(u), since u is defined.

Answer 87

f is a function made of other functions, which are inside of others.
so for that to work and accept a value i needed create function handles(@(u))

and then type up the function

Answer 88

f(u) will online work in the loop for newtons method.

Answer 89

Exactly, that's what I call definition. The u is just a dummy variable, takes no value when you defined it.
f(u) will execute and output a value, just like in the loop.
... unless I don't get what your problem is.

Answer 90

I mean @(u) is just a dummy variable, it would have worked with @(q), as long as you're consistent.

Answer 91

So you don't expect numbers to spill out on definition.

Answer 92

i get what your saying.
but to understand me i'll give show u somethung,
if i subsitute u = 0.1 in this eqution should get
Re = D_b*abs(u)*den_f)/visc_f.
and i get 7.04

Answer 93

but i get Re= D_b*abs(u)*den_f)/visc_f. IN COMMAND WINDOW. as soon as i add a @

Answer 94

What I thought you expected numbers to come out after:
```
f =

@(u)bouyancy-(Drag(u).*CirArea*den_f.*u.*abs(u)./2)
```
I guess I misunderstood.
So what is the problem bugging you?

Answer 95

display(Re(0.1));

ans =

7.0400