필터 지우기
필터 지우기

Solve integral equation efficiently

조회 수: 13 (최근 30일)
jcd
jcd 2023년 8월 3일
댓글: jcd 2023년 8월 3일
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
Torsten
Torsten 2023년 8월 3일
편집: Torsten 2023년 8월 3일
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 2023년 8월 3일
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.

댓글을 달려면 로그인하십시오.

채택된 답변

Torsten
Torsten 2023년 8월 3일
편집: Torsten 2023년 8월 3일
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
Torsten
Torsten 2023년 8월 3일
편집: Torsten 2023년 8월 3일
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 2023년 8월 3일
It does! Thank you again.

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Solver Outputs and Iterative Display에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by