vpasolve does not solve simple equation

I have a simple function and equation to solve: A - 2/(9*(1 - x)^(1/2)) - B*x =0, where A>B>1. The function is clearly continuous, at x = 0 it is positive, and it goes to minus infinity as x goes to 1. This means that there must be a solution on the interval (0, 1). Indeed, Mathematica says the solution is x = 0.987436 for A = 6.742 and B = 4.82. Yet Matlab finds no solution:

syms x; A = 6.742; B = 4.82;

vpasolve((A - 2/(9*(1 - x)^(1/2)) - B*x) == 0, x, [-0.001 1])

ans =

Empty sym: 0-by-1

My good friend ChatGPT suggests that the problem is that the solution is too close to the singularity (i.e. 1). So it suggested restricting the search interval (which I really don’t want, because the solution can indeed be very close to 1), giving good initial values, but none of this works. Curiously, if I get rid of the last term B*x, everything is fine:

vpasolve((A - 2/(9*(1 - x)^(1/2))) == 0, x, [-0.001 1])

ans =

0.99891358069017265865646884840762

Interestingly, this value is indeed very close to 1, closer than the Mathematica solution to the original equation, but poses no problem for Matlab.

Questions:

1. What am I doing wrong, why can’t Matlab solve the original simple equation?
2. Can it be fixed and if so how?
3. Bonus: which AI chatbot is best to use for coding questions?

Thank you!! aron

[Stephen23]

"which AI chatbot is best to use for coding questions?"

Currently: Claude. ChatGPT is good for writing emails. But as John D'Errico stated, this could be different next week.

And AI is just a tool: in the hands of a good programmer it can make writing code faster, but it does not (yet) replace that good programmer. If you want to use AI to help write good code then you still need to be able to write and iteratively improve prompts, review the results critically, and know where it missed something.

AI does not currently replace experience.

Answer 1

[John D'Errico]

Well, look more closely at your question. Plot it.

syms x; A = 6.742; B = 4.82;
F = A - 2/(9*(1 - x)^(1/2)) - B*x

F = 

fplot(F,[0,1]),grid on

The solution is not really extremely close to the singularity, but vpasolve uses a random starting value by default. And most of the time for this problem, that random starting value gets it in trouble, because for most values of x in the interval [0,1], a locally linear approximation to F predicts the solution must lie considerably beyond 1. And in that domain, things get a bit sticky, complexly so.

vpasolve(F)

 
ans =
 
Empty sym: 0-by-1
 

It can often help to give a starting value to any nonlinear solver. And since you know the solution will be a little less than 1, try a start point in that vicinity. Even 0.95 will suffice.

vpasolve(F,0.95)

ans = 

0.98743613546390317836514353721343

As far as knowing which chatbot is best for anything, that is a rapidly moving target, especially if you want a bot that actually gives a valid answer instead of a fabrication.

Oh, if you dispose of the linear term? That very much changes the shape of the problem and you don't even need a nonlinear solver anymore. An analytical one will work now. You could even do it with paper and pencil.

fplot(A - 2/(9*(1 - x)^(1/2)),[0 1])

solve(A - 2/(9*(1 - x)^(1/2)))

ans = 

[Aron]

Hi John, many thanks! Of course, this is a simplified version, A and B are changing parameters in a loop, and the solutions can indeed be anywhere on (0,1) in principle, that is why I was wary of giving starting values, especially so close to the solution. I am very surprised that it seems so sensitive to starting values (e.g. 0.9 does not work, but 0.91does), considering that it is a very well-behaved function for a solution algorithm to tackle (i.e. continuous and monotone) – at least so I thought! Many thanks again John, I really appreciate your help!

[Torsten]

considering that it is a very well-behaved function for a solution algorithm to tackle

A function with a singularity (like yours) is never well-behaved for numerical methods - especially if you want a value that is near to this singularity,

[John D'Errico]

Nonlinear solvers are always tricky things. If you want robustness, I'd suggest something like this:

syms x; A = 6.742; B = 4.82;
F = A - 2/(9*(1 - x)^(1/2)) - B*x;

xtest = sort(rand(1,1000));
Ffun = matlabFunction(F);
Ftest = Ffun(xtest);

[~,ind] = min(abs(Ftest));
xstart = xtest(ind)

xstart = 0.9881

Now, any solver will suffice, because this will be a very good start point.

vpasolve(F,xstart)

ans = 

0.98743613546390317836514353721343

Of course, we could just use fzero. And that will work, even without resorting to a starting value.

fzero(Ffun,[eps,1-eps])

ans = 0.9874

So then, can I use the same bracket for vpasolve?

vpasolve(F,x,[eps,1-eps])

 
ans =
 
Empty sym: 0-by-1
 

And strangely, vpasolve is a little more sensitive. But if I use a different slightly narrower bracket, it works well enough, again, without needing to find a decent starting value.

vpasolve(F,x,[0.000000000001,0.999999999999])

ans = 

0.98743613546390317836514353721343

[Aron]

Wow, this is super useful, thank you very much John!

[Torsten]

syms x; A = 6.742; B = 4.82;
F = A*(9*(1 - x)^(1/2)) - 2 - B*x*(9*(1 - x)^(1/2));
vpasolve(F==0,x)

ans = 

0.98743613546390317836514353721343