Question

I have following equation:

f(N):  N = ((1+lam)^3 )/ ((1-lam)*(1+lam^2));

I need to create a function that finds lam for specified N.

Right now I'm doing it using simple loop:

lam = 0.9999;
n = f(lam);
pow = 0;
delta = 0.1;
while(abs(N - n)) > 0.1 & pow < 10000)
    lam = lam - 0.001;
    n = f(lam)
    pow = pow+1;
end

How can I solve it more accurate and without using loops?

Was it helpful?

Solution

If you have

N = ((1+lam)^3 )/ ((1-lam)*(1+lam^2))

then you know that

(1+lam)^3 = N*(1-lam)*(1+lam^2)

Suppose you were to expand these terms? Coalesce into one simple cubic equation, with real coefficients, equal to zero? Is there a function that will solve it for you?

The answer is yes. One solution might be to use fzero, but since the equation is just a cubic polynomial, roots is the answer unless you needed a symbolic solution. Use the symbolic toolbox for symbolic problems.

OTHER TIPS

Here is a solution for N=10 by Wolfram Alpha:

http://www.wolframalpha.com/input/?i=(1%2Bx^3)/((1-x)*(1%2Bx^2))%3D10

An algebraic solution will work for your particular case, because it's not terribly difficult. The problem is that, in general, non-linear equations require an iterative solution: start with a guess, step in a particular direction, and hopefully converge to a solution. You can't solve non-linear equations in general without iteration and looping.

Rearrange the equation to be 0 = f(x)/g(x) (where f and g are polynomials). Then solve for 0 = f(x). This should be easy enough as f will be cubic (http://en.wikipedia.org/wiki/Cubic_function#Roots_of_a_cubic_function). In fact, Matlab has the roots() function to do this.

Plotting suggest that for N positive, there is exactly one solution in the interval [-1,1). You should consider Newton's method, it will converge for a zero initial guess fairly quickly.

You can solve this equation in closed form, as discussed in other answers, but to be honest, closed-form solutions to polynomials of degree > 2 are not very useful in practice, because the results tend to be poorly conditioned.

For your particular polynomial, I agree with Alexandre that Newton's method is probably the way to go.

In the long run, though, I highly recommend writing (or reusing from the Internet) an implementation of the Jenkins-Traub root-finding algorithm. Wikipedia describes it as "practically a standard in black-box polynomial root-finders," and they're not exaggerating. It has served all of my polynomial-solving needs for years; in my experience it's more robust than Newton's method (no reliance on a good initial guess) and eigenvalue-based methods, and is quite fast to boot.

There is an algebraic solution to your problem for most values of N. Here is the solution, as solved by Wolfram Alpha:

if N+1!=0
   x = (20 N^3+18 N^2+3 sqrt(3) sqrt(16 N^6+32 N^5-20 N^4-72 N^3-9 N^2+54 N+27)-27 N-27)^(1/3)/(3 2^(1/3) (N+1))-(2^(1/3) (2 N^2+3 N))/(3 (N+1) (20 N^3+18 N^2+3 sqrt(3) sqrt(16 N^6+32 N^5-20 N^4-72 N^3-9 N^2+54 N+27)-27 N-27)^(1/3))+N/(3 (N+1))

Yes, it's ugly.

If you have one, an exact algebraic solution, even a big ugly one like this one, is always superior to a numerical solution. As duffymo indicated, solving a problem with numerical methods require iterations (so it's slow), and the solver can get stuck in local minima.

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top