Solving an implicit equation with fsolve and validating with fplot

10 Ansichten (letzte 30 Tage)
Below is my code to try and solve this function:
D/(log10(RX*cfe)-C))^2 = cfe
I'm not sure if I'm using fsolve correctly. I get a value for cfe as 0.0016 but when I try to some validation by plotting two sides of the equation the curves cross at .0017377. That's about a 7% difference in values so I am questioning if I'm using fsove correctly.
I'm actually unsure of how to implement the initial guess which I call cfi.
I thought maybe adjusting the options would make the code run longer and would get a closer answer to the plotting method, but I'm not even sure if I'm doing that right.
Any direction would be helpful.
% AL, BE, E, and TWTD are all constants, which makes D and C also constants.
D = 0.242*(asin(AL)+asin(BE))/sqrt(E);
C = 1.26*log10(TWTD);
RX = 3.5e6 % Constant
options = optimoptions('fsolve','FunctionTolerance',1.0000e-1);
cfi = .455/(log10(RX))^2.58; % Initial guess
f = @(cfe,cfi) (D/(log10(RX*cfi)-C))^2 - cfe;
fsol = fsolve(@(cfe) f(cfe,cfi),0,options);
%------------------------------------------------------
syms cfe
eqnLeft = 0.242*(asin(AL)+asin(BE))/(log10(RX*cfe)-C);
eqnRight = sqrt(cfe)*sqrt(E);
fplot([eqnLeft eqnRight])
title([texlabel(eqnLeft) ' = ' texlabel(eqnRight)])

Akzeptierte Antwort

Ryan McBurney
Ryan McBurney am 1 Feb. 2021
Fixed my problem with the following adjustments to my code.
cfi = .455/(log10(RX))^2.58;
f = @(cfe,RX,Me) ((D)/(log10(RX*cfe)-C))^2 - cfe;
fsol = fsolve(@(cfe) f(cfe,RX,Me),cfi);
I had a suspicion I was using the initial guess the wrong way and I believe that is what is giving me a correct answer. The answer now matches the fplot result.

Weitere Antworten (1)

Alan Weiss
Alan Weiss am 31 Jan. 2021
Bearbeitet: Alan Weiss am 31 Jan. 2021
You have set an enormous value of the FunctionTolerance option. Don't set any options and I expect that you will get a better answer.
If this advice does not help, then please give numeric values for your constants so that we might have a chance at running your code and seeing what is happening.
Alan Weiss
MATLAB mathematical toolbox documentation
  1 Kommentar
Ryan McBurney
Ryan McBurney am 1 Feb. 2021
Here is the full code. Eventually I want to make Me and REN a vector of inputs.
clear all; clc;
Me = 3.5;
REN = 1.3e6;
LCFT = 12;
g = 1.4;
RX = REN*LCFT;
E = (g-1)/2 * Me^2;
TWTD = 1 + 0.89 * E;
A = E/TWTD;
B = (1+E)/(TWTD) - 1;
AL = (2*A-B)/sqrt(B^2+4*A);
BE = B/sqrt(B^2+4*A);
D = 0.242*(asin(AL)+asin(BE))/sqrt(E);
C = 1.26*log10(TWTD);
cfi = .455/(log10(RX))^2.58; % Initial guess
f = @(cfe,cfi) (D/(log10(RX*cfi)-C))^2 - cfe;
fsol = fsolve(@(cfe) f(cfe,cfi),0);
%-------------------------------------------
syms cfe
eqnLeft = 1/(log10(RX*cfe)-C);
eqnRight = sqrt(cfe)*sqrt(E)/(0.242*(asin(AL)+asin(BE)));
fplot([eqnLeft eqnRight])

Melden Sie sich an, um zu kommentieren.

Kategorien

Mehr zu Mathematics finden Sie in Help Center und File Exchange

Produkte

Community Treasure Hunt

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

Start Hunting!

Translated by