Solve integral equation efficiently

6 views (last 30 days)
jcd
jcd on 3 Aug 2023
Commented: jcd on 3 Aug 2023
I'm trying to find the parameter u in the following equation:
where and are known constants. So far, I succeded finding u by creating a vector of, say, 20 values in a range I think is reasonable, plotting it and the narrowing the range. This is not optimal because it takes so long to compute even 20 values and I have to do this process for at least 50 other cases. Here is the code I made:
nu = 1.04 * 10^(-6); %m2/s
Aplus = 26; kappa = 0.41; z = 0.3;
u_tau = linspace(0.0000001,0.05,20); % u_tau correct answer is 0.0000448
meanU = 4.3556e-04;
uplus2 = nan(1, length(u_tau));
uplus = meanU./u_tau;
yplus = u_tau * z / nu;
Unrecognized function or variable 'z'.
for ii = 1:length(yplus)
syms y; fun = ( 2 ) / ( 1 + ( 1 + 4 * kappa^2 * y^2 * ( 1 - exp( -y / Aplus ) )^2 )^(1/2) );
uplus2(ii) = double(int(fun, y, 0, yplus(ii)));
end
err = abs(uplus-uplus2);
plot(u_tau, err)
I also tried Newton-Raphson and the Secant methods but they fail too (I think with the derivative computation for NR). I also tried this:
syms ut
fun1 = ( 2 ) / ( 1 + ( 1 + 4 * kappa^2 * (ut*z/nu)^2 * ( 1 - exp( -(ut*z/nu) / Aplus ) )^2 )^(1/2) );
up2 = int(fun1, ut, 0, (ut*z/nu));
assume(ut > 0 & ut < 0.05)
ut_sol = vpasolve(up2 - (meanU/ut) == 0, ut); % I also used solve()
but it gives me an incorrect answer of 1.24 which I know it's not possible becuase the correct answer is very small.
How can I solve this the most optimal way without having to look at a plot to determine the correct value (and also as fast as possible)? How can I make sure there is only one correct answer? Ideally I'd like to plot the zero crossings too but that is not urgent.
I use Matlab R2021b if it helps.
Thanks in advance.
  2 Comments
Torsten
Torsten on 3 Aug 2023
Edited: Torsten on 3 Aug 2023
You didn't specify z and mnU. And why do you call the upper limit of the integral and the integration variable both y ? Should the expression for y = u*z/nu also be substituted for the integration variable ?
jcd
jcd on 3 Aug 2023
Just put the value of z and fixed mnU. I did substitute y = u*z/nu in the second attempt I put in the question but I'm not sure if that is giving me the wrong answer.

Sign in to comment.

Accepted Answer

Torsten
Torsten on 3 Aug 2023
Edited: Torsten on 3 Aug 2023
nu = 1.04 * 10^(-6); %m2/s
Aplus = 26;
kappa = 0.41;
z = 0.3;
meanU = 4.3556e-04;
fun = @(x) 2 ./ ( 1 + sqrt( 1 + 4 * kappa^2 * x.^2 .* ( 1 - exp( -x / Aplus ) ).^2 ) );
f = @(y) meanU./(nu*y/z) - integral(fun,0,y);
y = fzero(f,20)
y = 12.9167
f(y)
ans = 0
u = nu*y/z
u = 4.4778e-05
  4 Comments
Torsten
Torsten on 3 Aug 2023
Edited: Torsten on 3 Aug 2023
Say you have the equation x^2 - 4 = 0 and you want to use fzero to get a solution. If you set x0 = -1.9, you get -2 and if you set x0 = 1.9, you get 2. Does this answer your question ?
Thus a good initial guess is very important to get a solution and - at best - the "correct" solution.
Using "fsolve" instead of "fzero" could be another option for your zero-crossing problem. First tests with your problem from above show that
y = fsolve(f,x0)
seems to have a larger domain of convergence for x0 than
y = fzero(f,x0)
jcd
jcd on 3 Aug 2023
It does! Thank you again.

Sign in to comment.

More Answers (0)

Categories

Find more on Mathematics in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!