Items tagged with fsolve


The issue I am currently having is that, while analyticity (and physics) indicates a certain function must have roots, fsolve is having trouble finding them. In fact, I have even found roots manually in a certain region myself, simply inputting into the function various values until I found them. However, fsolve does not seem to want to find these roots, and I believe it is a numerics issue: when I changed the digits around, for extremely low values of Digits, it would find a root (even though it was incorrect). Further, this problem arose elsewhere in the domain of interest for other values of Digits (in particular, for Digits:=5, fsolve failed in a region it had not failed before).

The region of interest is the "peak" of the output of poleR(M0, 0.935, mK), which should be somewhere around M0 = 0.95 or so. However, because fsolve cannot find the roots, the plot cannot be made.

Anyone have any ideas as to why fsolve cannot find the roots? I was also experiencing issues with some of these functions having multiple roots, which itself is weird as well (note that I am working over the complex plane).

Attached is the document.

Any ideas?

How can this error be corrected '' error, (in fsolve) fsolve cannot solve on 0=0 ''. See the worksheet



I have a problem solving two equations.  They are as follows:

s := 1/(273.16+50); s1 := 1/(273.16+145); s3 := 1/(273.16+250); s2 := 1/(273.16+197.5); gamma0 := 0.1e-3; gamma1 := .5; gamma2 := 0.15e-2; beta := -3800:

c := 300; n := 200; tau1 := 99; tau2 := 120;

Delta := solve(1-exp(-(gam0*tau1+(1/2)*gam1*tau1^2)*exp(beta*s1)) = 1-exp(-(gam0*a+(1/2)*gam1*a^2)*exp(beta*s2)), a);
a := Delta[1];

Theta := solve(1-exp(-(gam0*(a+tau2-tau1)+(1/2)*gam1*(a+tau2-tau1)^2)*exp(beta*s2)) = 1-exp(-(gam0*b+(1/2)*gam1*b^2)*exp(beta*s3)), b);
b := Theta[1];

n1 := int((gam1*t+gam0)*exp(beta*s1)*exp(-(gam0*t+(1/2)*gam1*t^2)*exp(beta*s1)), t = 0 .. tau1);
n22 := (n-n1)*(int((gam1*t+gam0)*exp(beta*s2)*exp(-(gam0*t+(1/2)*gam1*t^2)*exp(beta*s2)), t = a1 .. a1+tau2-tau1));
n2 := eval(n22, a1 = a);
n33 := (n-n1-n2)*(Int((gam1*t+gam0)*exp(beta*s3)*exp(-(gam0*t+(1/2)*gam1*t^2)*exp(beta*s3)), t = b1 .. c));
n3 := eval(n33, a1 = a);
n4 := n-n1-n2-n3;

g1 := -n1*(Int((1/(gam1*t+gam0)-t*exp(beta*s1))*(gamma2*t^2+gamma1*t+gamma0)*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1)), t = 0 .. tau1))-n2*(Int((1/(gam0+gam1*(a+t-tau1))-(a+t-tau1)*exp(beta*s2))*(gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2)*exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2)), t = tau1 .. tau2))-n3*(Int((1/(gam0+gam1*(b+t-tau2))-(b+t-tau2)*exp(s3))*(gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3)), t = tau2 .. c))+(n-n1-n2-n3)*(1/(gam0+gam1*(b+c-tau2))-(b+c-tau2)*exp(s3))*(gamma0+gamma1*(b+c-tau2)+gamma2*(b+c-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+c-tau2)+(1/2)*gamma1*(b+c-tau2)^2+(1/3)*gamma2*(b+c-tau2)^3)*exp(beta*s3));

g2 := -n1*(Int((t/(gam1*t+gam0)-(1/2)*t^2*exp(beta*s1))*(gamma2*t^2+gamma1*t+gamma0)*exp(beta*s1)*exp(-(gamma0*t+(1/2)*gamma1*t^2+(1/3)*gamma2*t^3)*exp(beta*s1)), t = 0 .. tau1))-n2*(Int(((a+t-tau1)/(gam0+gam1*(a+t-tau1))-(1/2)*(a+t-tau1)^2*exp(beta*s2))*(gamma0+gamma1*(a+t-tau1)+gamma2*(a+t-tau1)^2)*exp(beta*s2)*exp(-(gamma0*(a+t-tau1)+(1/2)*gamma1*(a+t-tau1)^2+(1/3)*gamma2*(a+t-tau1)^3)*exp(beta*s2)), t = tau1 .. tau2))-n3*(Int(((b+t-tau2)/(gam0+gam1*(b+t-tau2))-(1/2)*(b+t-tau2)^2*exp(s3))*(gamma0+gamma1*(b+t-tau2)+gamma2*(b+t-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+t-tau2)+(1/2)*gamma1*(b+t-tau2)^2+(1/3)*gamma2*(b+t-tau2)^3)*exp(beta*s3)), t = tau2 .. c))+(n-n1-n2-n3)*((b+c-tau2)/(gam0+gam1*(b+c-tau2))-(1/2)*(b+c-tau2)^2*exp(s3))*(gamma0+gamma1*(b+c-tau2)+gamma2*(b+c-tau2)^2)*exp(beta*s3)*exp(-(gamma0*(b+c-tau2)+(1/2)*gamma1*(b+c-tau2)^2+(1/3)*gamma2*(b+c-tau2)^3)*exp(beta*s3));

solve({g1 = 0, g2 = 0}, {gam0, gam1});

Warning, solutions may have been lost.

What do I do wrong?

Thanks for advice in advance.


I have a rather complex expression that I want to find the zero for as a function of two other parameters, i.e. I have a function

Denom := (s,M,g) -> stuff

that I want to find the zero of for a variety of values of M and g. In some cases the solution will be complex, which is entirely acceptable. However, the real part of the solution should never be negative, and yet that is the kind of result I am getting.

As an example (illustrated in the worksheet, when attempting to find the zero for M = 3 and g = 0.2, fsolve gives me s = -6.1 -1.4i. However, when I plot the function with the parameters input already, I can clearly see a zero at s = 9 with no imaginary component. Why won't fsolve find this zero? How can I make it do so?

See the bottom of the attached worksheet for the main problem.

Could anyone assist in rectifying this error ''Error, (in fsolve) {f[1], f[2], f[3], f[4], f[5], f[6], f[7], f[8], f[9], f[10], f[11], theta[11]} are in the equation, and are not solved for''. Here is the worksheet

hello, i went solve these equation ,with a & b take any value


thank you

I have a nonlinear system with 4 equations and 4 unknowns. I am using fsolve. I know that there are multiple solutions for each variable but am only getting one. I need the others. what do I do??

This is my code:

R__1 := Matrix([[1, 0] , [0, 1] ]);

R__2 := Matrix([[1/2, sqrt(3)/2] , [-sqrt(3)/2, 1/2] ]);

R__3 := Matrix([[-1/2, sqrt(3)/2] , [-sqrt(3)/2, -1/2] ]);

R__4 := Matrix([[-1, 0] , [0, -1] ]);

R__5 := Matrix([[-1/2, -sqrt(3)/2] , [sqrt(3)/2, -1/2] ]);


d__1 := Vector( [ 0, 5.4] );

d__2 := Vector( [ 6.4, 4.539] );

d__3 := Vector( [ 11, 4.078] );

d__4 := Vector( [ 15.5, 2.079] );

d__5 := Vector( [ 19, 1.039] );


a := Vector( [ a__x, a__y] );
















Thanks for any tips you may be able to offer


Is it possible in Maple 15 to solve an equation with a parameter for a given set of parameters? How can this be passed to the solve function, should I use some kind of list?  After obtaining the solution how can I assign the solutions to variables such as x1 for the first value of the parameter, x2 for the second value of the parameters and so on. Furthermore, is this possible with the fsolve command?



Hello I am a Maple 15 user and I am using the command fsolve to solve for the intersection of two curves over a specified interval in x, namely from 0 to the lim defined in the Maple document. The specified interval contains asymptotes and when I specify the full interval only one of the three solutions is returned even if I can see that there are three distinct solutions by looking at the plot of RHS and LHS. Should I use another technique to find the solution or is my implementation of fsolve command wrong?

Thanks in advance



n1 := 1:

n2 := 1.50:

n3 := 1.40:

lambda := 1.3:

k0 := 2*Pi/lambda:

d := 3:

x0 := k0*d:

arg1 := sqrt(x0^2*(n2^2-n1^2)):

arg2 := sqrt(x0^2*(n2^2-n3^2)):

lim := FindMinimalElement([arg1, arg2]):

sqr1 := sqrt(x0^2*(n2^2-n1^2)-x^2):

sqr2 := sqrt(x0^2*(n2^2-n3^2)-x^2):

LHS := tan(x):

RHS := (sqr1+sqr2)/(x*(1-sqr1*sqr2/x^2)):

plot([LHS, RHS], x = 0 .. lim, y = -6 .. 6)


fsolve(RHS = LHS, x = (1/2)*Pi .. 3*Pi*(1/2))



fsolve(RHS = LHS, x = 3*Pi*(1/2) .. 9*Pi*(1/4))



fsolve(RHS = LHS, x = 9*Pi*(1/4) .. lim)






I have to numerically solve this equation many (tens of thousands) of times for theta_n (as I vary different parameters, especially n):

tan(theta_n) - C_a*Z_0*(-v^2*(Pi*n-theta_n)^2/x_l^2+omega_a^2)*x_l/(v*(Pi*n-theta_n)) = 0

theta_n should be in [-pi/2, pi/2). It seems like solving this is the slowest single component of the chain of calculations (that follow this).

Currently I do it with this function:

get_theta_n_array:=proc(max_n::integer, omega_a::float, v::float, x_l::float, C_a::float, Z_0::float)
    local theta_n_array:=Array(
        select(x->Re(x)<evalf(Pi/2) and Re(x)>=-evalf(Pi/2), [seq(
        #NOTE: careful with fsolve - in some cases returns unevaluated equation
        fsolve(tan(theta_n) - C_a*Z_0*(-v^2*(Pi*n-theta_n)^2/x_l^2+omega_a^2)*x_l/(v*(Pi*n-theta_n)) = 0, theta_n) , n=1..max_n)]
        , real)
    if ArrayNumElems(theta_n_array) <> max_n then
        printf("Bad Array Dimensions! Got too many or not enough solutions.");
        theta_n_array:="CHECK: get_theta_n_array()": #dirrrrty hack that will ring an alarm bell if array is not the right size

And call it like so (for say n=1000)

result:=get_theta_n_array(1000, 100e9, 1e8, 0.3, 20e-14, 50.0);

This will take say 3.5s on my PC.

Does anyone have any ideas how to speed this up? I would hope this to take at least an order of magnitude less time. I played with DirectSearch lib but that was not faster.

Also, I should note that this is the only portion of my code that is not thread safe (because of the fsolve call), which leads to "extra" slowdowns because I have to use Grid:-Map, instead of Thread:-Map when parallelizing, and more importantly because I can't compile the rest of the code (Grid:-Map is not compatible with compiled functions).

Let me know if you have any ideas...


i want to solve an equation by fsolve but i cant assign a value as an input for next step!

please help me

s := fsolve(G), x = -1 .. 1     

s := .1449607418, x = -1 .. 1  


Error, invalid input: subs received .1449607418, which is not valid for its 1st argument                  


I want to solve system of equation but it has unknow parameter.

Then I test system of equation. It hasn't unknowparameter.

eq1 := x^2+y^2 = 4

eq2 := y-x^2 = 0

fsolve({eq1, eq2}, {x, y})

{x = -1.249621068, y = 1.561552813}

So I get answer by using fsolve.


Then I try to put unknow parameter in system of equation.

eq3 := x^2+ky^2 = 4

eq4 := ay-hx^2 = 0

fsolve({eq3, eq4}, {x, y})

Error, (in fsolve) {ay, hx, ky} are in the equation, and are not solved for

I don't get answer and open link. The link hasn't similar this problem.

To motivate some ideas in my research, I've been looking at the expected number of real roots of random polynomials (and their derivatives).  In doing so I have noticed an issue/bug with fsolve and RootFinding[Isolate].  One of the polynomials I came upon was

f(x) = -32829/50000-(9277/50000)*x-(37251/20000)*x^2-(6101/6250)*x^3-(47777/20000)*x^4+(291213/50000)*x^5.

We know that f(x) has at least 1 real root and, in fact, graphing shows that f(x) has exactly 1 real root (~1.018).  However, fsolve(f) and Isolate(f) both return no real roots.  On the other hand, Isolate(f,method=RC) correctly returns the root near 1.018.  I know that fsolve's details page says "It may not return all roots for exceptionally ill-conditioned polynomials", though this system does not seem especially ill-conditioned.  Moreover, Isolate's help page says confidently "All significant digits returned by the program are correct, and unlike purely numerical methods no roots are ever lost, although repeated roots are discarded" which is clearly not the case here.  It also seems interesting that the RealSolving package used by Isolate(f,method=RS) (default method) misses the root while the RegularChains package used by Isolate(f,method=RC) correctly finds the root.

 All-in-all, I am not sure what to make of this.  Is this an issue which has been fixed in more recent incarnations of fsolve or Isolate?  Is this a persistent problem?  Is there a theoretical reason why the root is being missed, particularly for Isolate?

Any help or insight would be greatly appreciated.

Hi! I'm trying to solve a system of four non-linear equations in Maple 17 but it doesn't work.

Equations are: F, Fw, Ft, Fk and varibles are T,w,k,ki.

Parametars are Mn=10, Ms=2, alfa=0.2 and Tf=(k*T^(alfa)/Mn)^(1/alfa).

Solutions must be positive. 

This is maple script:

> assume (w>0):
> assume (T>0):
> assume (k>0):
> assume (ki>0):
> assume (Mn>0):
> assume (Ms>0):
> assume (Tf>0):
> assume (alfa>0):

> Gp:=1/exp(sqrt(w*I));

> Cf:=((T*w*I+1)/(Tf*w*I+1));

> Cpi:=(k*w*I+ki)/(w*I);

> L:=Cf*Cpi*Gp;
> L:=evalc(L):
> F:=subs(Ms=2,Tf=(k*T^(alfa)/Mn)^(1/alfa),evalc(abs(1+L)^2-1/Ms^2));
> F:=subs(Mn=10,alfa=0.2,evalc(F));

> Fw:=diff(F,w):
> Fk:=diff(F,k):
> Ft:=diff(F,T):
> fsolve({F,Fw,Fk,Ft},{w,k,ki,T});

Thanks in advance for any help. Dragoslav

i := I;






Hz := k*(z^2-2*r*a+r^2)/((z-1)*(z^2-2*b*z+1));

Hzw := eval(Hz, z = exp(i*w)); assume(a > 0);

Habs := simplify(abs(Hzw)^2);

p1 := eval(Habs, w = Pi) = (10^((-3.3018)*(1/20)))^2;

p2 := eval(Habs, w = (1/2)*Pi) = (10^((-.1758)*(1/20)))^2;

p3 := eval(Habs, w = (1/4)*Pi) = (10^(6.425*(1/20)))^2;

p4 := eval(Habs, w = (1/8)*Pi) = (10^(54.5699*(1/20)))^2;

solve({p1, p2, p3, p4})


4 5 6 7 8 9 10 Page 6 of 13