there are a few general mathematical considerations that i think will help, but one of them implies that an analytic solution is unlikely

if we write f(x) for the LHS of the equation, so that f(x) = a*x^b - c*x^d, and we only consider x>0 (as we have been told), then we can see that f(x) = 0 for exactly one value of x (as long as a != b and c != d; if either of these statements fails, the analysis is much easier)

it is easy to show that as x--> infinity, f(x) --> 0 also (as long as b and d are negative, which i see is not always the case - however, for b 0 for x infinity as x --> 0, there clearly must be exactly one solution for a give e>0 - In that case also, f(x) < 0 for x>g, and this is the only solution

if f'(g) > 0, then f(x) < 0 for x < g, so there is no solution for any e > 0 in the interval [0, g) - In that case, f(x) > 0 for x>g, so we can see that there may be a solution for e (so f(x) starts at 0 when x=g, rises to something, then falls gradually to 0 as x-->infinity)

here is the hard part - it is also easy to show that f'(x) = 0 for exactly one point, so either there are no solutions (the difference in the two exponential expressions never gets big enough), or there are two solutions (the difference gets larger than e), or there is exactly one (when the maximum difference is exactly e)

f'(x) is easy to compute, and the unique x for which f'(x) = 0 is easy to compute, and it should be easy to show that the unique x for which f'(x) = 0 has x > g

a little long-winded, but i hope this helps

Have you tried to fit a model to your numerical solution. May be in the least square sense. Try fitting the model and then examine the residue. See whether you can reduce the residue.

If you want to use the equation to approximate a group of data that represents the solution for 'x', you can use he Matlab optimization tool box. The tool box can easily solve this problem, for details you can visit: www.mathworks.com/products/optimization

Good luck!

You can get good approximate to the radius of the disk of the roots domain by the numerical radius of the companion matrix for this see for example Bounds and majorization for the real parts of the zeros of polynomials by F.Kittaneh

if you find numerical convergence is slow, may be you can try MCMC

Since b,d are negative my guess x should be positive otherwise you could get complex solutions. Consideration of approximation x is big positive gives a x_0^b=e if b>d. Further you may consider iterations ax_(n+1)^b=e+cx_n^d, n=0,1,2,3,...

What are the ranges of your parameters a,b,c,d,e?

Is there any constraint on their values?

Here is a simple comment: there are simple transformations when bd one of the coefficients has to become negative

demetris' question about parameter and variable ranges is an important part of making progress in analyzing this equation

I think it's better to rewrite the equation as a *x^(-b)-c*(x^(-d))-e=0 with b,d positive exponents and a,c,e positive constants. With this formulation, you can transform your equation in e*x^b+cx^(b-d)-a=0 if b>= d or in e*x^d-a*x^(d-b)+c=0.With this form we have polynomial equations for which there exist analytic solution if b and d are smaller or equal to 4 (like Cardan formula).

Thanks a lot for the feedback, guys! In the majority of cases I have simulated, the values of b and d are negative (using my original notation), but the exponent d can assume positive values. The range of b is -0.67 to -0.01, and the range of d is -1.4 to 9. The other constants a, c, and e all have positive values

Is this equation somewhat important. Is there any applications behind the equation. I can develop an approximate solution using the same reasoning of estimation of polynomial roots by perturbation theory.

Article Perturbation theorems for estimating magnitudes of roots of ...

there are a few general mathematical considerations that i think will help, but one of them implies that an analytic solution is unlikely

if we write f(x) for the LHS of the equation, so that f(x) = a*x^b - c*x^d, and we only consider x>0 (as we have been told), then we can see that f(x) = 0 for exactly one value of x (as long as a != b and c != d; if either of these statements fails, the analysis is much easier)

it is easy to show that as x--> infinity, f(x) --> 0 also (as long as b and d are negative, which i see is not always the case - however, for b 0 for x infinity as x --> 0, there clearly must be exactly one solution for a give e>0 - In that case also, f(x) < 0 for x>g, and this is the only solution

if f'(g) > 0, then f(x) < 0 for x < g, so there is no solution for any e > 0 in the interval [0, g) - In that case, f(x) > 0 for x>g, so we can see that there may be a solution for e (so f(x) starts at 0 when x=g, rises to something, then falls gradually to 0 as x-->infinity)

here is the hard part - it is also easy to show that f'(x) = 0 for exactly one point, so either there are no solutions (the difference in the two exponential expressions never gets big enough), or there are two solutions (the difference gets larger than e), or there is exactly one (when the maximum difference is exactly e)

f'(x) is easy to compute, and the unique x for which f'(x) = 0 is easy to compute, and it should be easy to show that the unique x for which f'(x) = 0 has x > g

a little long-winded, but i hope this helps

Find my suggested answer to the attachment file.

All you need is an estimation for the 'mean value' of the size x.

Then you can take a Taylor series expansion around that point and build approximate analytic solutions for y=1/x and then for x.

Good luck!

Thanks again, I appreciate your support on this. Mehmet Pakdemirli, the equation form is important in this case because it resulted from some mechanistic considerations on how fish growth depends on consumption rates and activity (both power functions of body size x), and all constants a-e have biological meaning.

Christopher Landauer, what if we constrain the analysis to f '(g)

recall f(x) = a*x^b - c*x^d for x > 0, constants a, c > 0 > b, d,

and g is the unique solution of f(x) = 0, that is,

g = (a/c) ^ (1/(d-b))

(if d=b, the entire analysis is simpler and was not provided),

and we are looking for a solution to f(x) = e

now if f'(g) < 0, then there is exactly one solution for f(x) = e for e > 0,

and that solution has x0, it might be better numerically to divide the original equation by x^(min(b,d)) to get an equation in x with positive exponents that is similar in form to the original - then when x gets small, the relevant functional values do not get large (which is always useful for computer-based analyses)

since all the derivatives of f(x) are easy to compute, you could try to use a second or third order extension to Newton's method - in fact, since f(x) is smooth, these approaches should work pretty well

if we write f^(k) for the kth derivative, k>=0 an integer, then

f^(k)(x) = a * (b^k) * (x^(b-k)) - c * (d^k) * (x^(d-k))

(as is easily seen from a few examples and proven by a simple induction)

expanding a taylor series around x=g a few terms gives a polynomial equation for a value of x that should more closely approximate the solution to f(x) = e (yes, this will be a polynomial even if b and d are not integers, since the taylor expansion takes the form f(g+h) = polynomial in h + error term, and it is h that we are looking for)

this approach is limited only by your patience and the algebraic capabilities of your math clerk (which i hope is a computer program 8-))

The most stable method is bisection, but it is very slow in convergence.

Fast method is constant point method (Newton' s is a special case of it), but you have to build it for your specific function:

f(x)=0 x = g(x)

can be implemented with many ways and not all of them converge to a solution.

Taylor expansion is the best choice for you because x, the size of fish, is known, at least for its mean value x0, so taking the corresponding series around x=x0 will give you a very good approximation of f(x).

(Or f(y) if you do the transformation y=1/x)

Another solution is Christopher' s suggestion, also based on Taylor series.

Perhaps you could give us a realistic equation to study it.

You can transform this equation to : X^p + X^q + C = 0 or : X^p + -X^q + C = 0

This you can Solve. ( C is a constante ) ( p and q are now positive.)

If only 2 or 3 significant digits are sufficient, then the really simple method is to use well known gnuplot program (designed first of all to make 2D or 3D plots). Just define the function f(x) being equal to the left side of your equation, then set the values of all constants and plot the result. Zoom the neighborhood of zero crossing(s), the need will be. This doesn't work (easily) for complex solutions.

Besides: since your constant are most likely of experimental origin, are you sure you know their values exactly? That's probably your next problem ...

Hi,

Perhaps computing the best approximating function could be used instead of Taylor expansion, as suggested before. The advantage of this is that you have a function that is approximate in a whole region and not around a single point. This "best" approximating function could be defined in the 2-norm (or another one if you want to).

I computed for an approximating quadratic polynomial but the coefficients are really extensive. Perhaps a mix with a numerical treatment could be done. In the attached file you can see the first part of this procedure, prior to the integration over the domain of interest.

Good luck,

Nicolas

Of course we can do the familiar in statisticians linear regression, based on the minimization of the l2 norm. The problem in the regression analysis is that you can never be sure that your model is close to reality, even if it localy seems to be identical with your data! That' s the reason why we are trying to ecape form regression...

Octav Olteanu, University Politehnica of Bucharest

Let us look for positive solutions. The idea is to use the convexity of a suitable function. If d/b>1, you make the change of variable y=x^b. The equation in y is: f(y)=c*y^(d/b)-a*y+e=0. Because of strictly convexity of f, also observing that f(0)>0, f(infinity)>0, the conclusion is that the equation has at most two positive roots. This happens if and only if f(y min)0, there are no positive solutions. Standard computations lead to

y min = (a*b/c*d)^(b/(d-b)). The condition f(y min)

I will give the general form in wich you can transform it with standaard algebra:

X^p + X^q + C = 0

An general example: X^4 + X^3 + -t^4 = 0 of degree 4

X^4 + X^3 + -5 = 0 X = 1.296153312363 X = - 1.823974456169 Brian

t = 5^(0.25) = 1.495348781 in this particular example in wich C = -5

Sept 1 develop an simple serie for X. ( study a simple way to do this )

X = t + -1/4 + 3/32t + -1/32*t^2 + 15/2048*t^3 + -77/65536*t^5

X = 1.495348781 + -0.25 + 0.062694403 + -0.013975424 + 0.00219045231

+ -0.00015714417 = 1.29610106846

Two interative corrections will give X to 11 digits. Wich I consider to be enhough.

X = 1.29610106846 + -(-0.00071833546/13.7488005) = 1.29615331514

X = 1.29615331514 + -(0.00000003627/13.75026008) = 1.29615331236

This methode is apliable to all cases and of any degree.

I hope this will help you.

With friendly Greetings Ruben Hogenhout

Standard algebra with y=1/x^d, p=b/d, C=c/a and E=e/a transforms your equation to y^p=Cy+E where again p, C and E are positive constants. Assuming p>1 (the case p

a*X^b-c*X^d=0 implies

a*X^b=c*X^d implies

loga+blogx=logc+dlogx

this implies

loga-logc=(b+d)logx

log(a/c)=(b+d)logx

log[(a/c)^1/(b+d)]=logx

So X=(a/c)^(1/(b+d))

Dear Marek, Sorry there was a small mistake in sign in the previous solution. The correct solution is as follows.

a*X^b-c*X^d=0 implies

a*X^b=c*X^d implies

loga+blogx=logc+dlogx

this implies

loga-logc=(d-b)logx

log(a/c)=(d-b)logx

log[(a/c)^(1/(d-b))]=logx

So X=(a/c)^(1/(d-b))

I've spotted this mistake, but I didn't comment it. Telepathy ??

I was ashamed not to present this EXACT solution myself.

Op's original equation has an e>0. :(

Correcting a little inexactness above, y=x^d, p=b/d, C=c/a and E=e/a transforms your equation to y^p=Cy+E where again p, C and E are positive constants. We are looking for positive solutions y. The case p=1 is trivial (only a positive solution for C> C E^(1/p), i.e. E^(1-1/p) >> C or E >> C^(p/(p-1)), resp. ]

c) If E/C \approx y then y \approx (2 E)^(1/p) [a posteriori: if E \approx 2^(1/p) C E^(1/p), i.e., E \approx 2^(1/(p-1)) C^(p/(p-1)) ]

Here, a posteriori means using the approximate solution to obtain conditions only involving the parameters C, E, and p of the problem.

But the above is not the whole story:

Now, for positive y, the RHS Cy+E is greater than E > 0 and increasing linearly. For p>1, the LHS increases monotonously with initially zero slope from y=0 and faster than linearly for increasing y. Thus, in this case, we have a single solution y>0. This solution depends on the ratio of E/C, i.e. on cases a),b), and c).

For p < 1, the LHS has a vertical tangent at y=0 and increases slower than linearly with increasing y. Hence, for positive E there may be zero, one or two intersections of the LHS with the RHS.

A) If E is very small, there are two solutions, one close to y=0 (corresponding to case b)), the other close to the above mentioned value C^(1/(p-1)) of case a).

B) If E is large, then there are no solutions.

C) The borderline is obtained when the LHS is tangent to the RHS at the solution, i.e. for the simultaneous equations p y^(p-1) = C, y^p=Cy+E,

whence y^p=C/p y = Cy+E and

thus y= E/(C(1/p-1)) = (C/p) ^ (1/(p-1)). Thus, for E=C(1/p-1) * (C/p)^(1/(p-1)) = C^(-p/(1-p)) * (1-p)/p * (p)^(1/(1-p)), we obtain exactly one solution. Note, that this is a relation of the form E= K(p) C^(-p/(1-p)). Up to the p-dependent factor K(p), note the relation to the a posteriori conditions on E in a), b) and c).

The strategy is clear: In order to obtain approximate solutions of a three-term equation, you may discuss the cases, where one term is negligible with respect to the others.

Some maple file is attached to illustrate the numerical performance of the simple approximations as discussed above.

Hello Henrique,

The following is an attempt to give you a visual insight that may help you solve the problem.

Your problem is to find approximate POSITIVE solutions (x > 0) of the following equation

ax^(-b) - cx^d = e

where a, c and e are positive constant coefficients in the unit interval

0 < a < 1

0 < c < 1

0 < e < 1

whereas b and d are constant exponents

0.01 < b < 0.67 (always positive)

-1.4 < d < 9.0

--------------------------------------------

First, reduce the number of independent parameters from 5 to 4.

Let E = e/a and let C = c/a

then, multiply the equation by x^b and divide by a to obtain

1 = Ex^b + Cx^(d+b)

-------------------------------------------------------------------------

Next, introduce the following change of variable: X = Ex^b and Y = Cx^(d+b)

Your equation becomes 1 = X + Y (straight line) where Y = f(X) is a power function of X.

Namely, f(X) = C(X/E)^P with P = 1 + (d/b).

Let us describe the general idea.

The problem now reduces to finding the intersection of the straight line 1 = X + Y with the curve Y = f(X). More importantly, this intersection (if any) will be necessarily inside the unit square of the (X, Y) space. In fact, since both X and Y must be positive and less than 1, then we can even calculate beforehand (when they exist) the interval [Xmin, Xmax] inside which this intersection will take place.

Consequently, you can either:

1) Draw a simple plot of the curve Y = f(X) with X running from Xmin to Xmax. This will help you locate visually the solutions X before transforming back to the original variable x = (X/E)^(1/b).

2) Write a code to get automatic solutions of this problem using a simple bisection algorithm starting from the initial interval [Xmin, Xmax] (when they exist).

The following figures might help you visualize this technique.

Figure 1 shows the very simple case P > 0. In that case, the curve will always start from (0, 0) and will grow monotonically yielding thus one unique valid intersection. The initial interval (either for plotting or for automatic searching) is easily derived as follows:

Xmin = 0 and Xmax is calculated such that f(Xmax) = 1. Thus Xmax = E( (1/C)^(1/P) ).

The special case P = 0 can be worked out analytically. In that case, X = 1 – C. Then, depending on C, this solution will be either accepted or rejected as X must be in [0, 1].

The tricky case P < 0 is show in figure 2 and figure 3. First, the curve will no longer start from (0,0). Instead, it will be a decreasing curve. Thus, now Xmax = 1 while Xmin must be pre-calculated as shown above by requiring f(Xmin) = 1.

Of course, if it turns out that Xmin > 1 then there will be no solution. However, even if Xmin is found to be less than 1 you may not necessarily have a solution. In fact, you may either have 2 solutions (figure 2) or no solution at all (figure 3). There is also a third possibility (not shown here) where the curve (such as that of figure 3) becomes tangent to the straight line, which will yield then a single valid solution.

I hope this will help