Drexel dragonThe Math ForumDonate to the Math Forum



Search All of the Math Forum:

Views expressed in these public forums are not endorsed by Drexel University or The Math Forum.


Math Forum » Discussions » sci.math.* » sci.math.num-analysis.independent

Topic: Solving a nonlinear system of equations OR finding all roots of a polynomial?
Replies: 4   Last Post: Mar 12, 2012 8:55 AM

Advanced Search

Back to Topic List Back to Topic List Jump to Tree View Jump to Tree View   Messages: [ Previous | Next ]
deltaquattro@gmail.com

Posts: 77
Registered: 7/21/06
Solving a nonlinear system of equations OR finding all roots of a polynomial?
Posted: Mar 8, 2012 4:59 AM
  Click to see the message monospaced in plain text Plain Text   Click to reply to this topic Reply

Hi,

I'm trying to solve

(a-L)^2/a*2/(a+b)-p=0
(b-U)^2/b*2/(a+b)-p=0

in Visual Basic (!), under the following constraints:

1. 0<p<1 (it's a probability)
2. a>L>0, b>U>0

A possible example is L=U=1, p = 0.022750132: the only solution which satisfies my constraints is [1.1776225370286649805,1.1776225370286649805]. Initially I thought that the system had to have a simple, closed solution, since it comes out of a simple problem with similar triangles. Little did I know that apparently simple math problems, often pose the hardest of challenges :) Luckily, the great folks at sci.math.symbolic (thanks, Martin, Waldek and all the others, for the invaluable help!!) showed me that the system can be reduced to an equation in b:

2*b^4*(p - 1)*(p - 2) + 2*b^3*(L*p*(p - 2) - U*(p^2 - 8*p + 8))+ b^2*(L^2*p^2 + 8*L*U*p + U^2*(p^2 - 16*p + 24))- 4*b*U^2*(L*p + 2*U*(2 - p)) + 2*U^4*(2 - p) = 0

whose solution I can substitute in the following equation for a:

a = (2*b^3*(p - 1)*(p - 2) + 2*b^2*(L*p*(p - 2) - U*(p^2 - 8*p + 8)) + b*(L^2*p^2 + 8*L*U*p + 4*U^2*(5 - 3*p)) - 4*U^2*(L*p + U*(2 - p)))/(U^2*p*(p - 2))

I tried solving f(b)=0 using Newton's method. Unfortunately, it doesn't always converge to the one root satisfying my constraints, out of the four possible. I guess that depends on the starting point: since b must be positive, I choose as a starting point one of the following three, depending on which one is positive:

1) b0=U (which is the solution in the case p=0)
2) b0=(-Sqrt(p) * Sqrt((4 - 2 * p) * U ^ 2 + (4 - 2 * p) * L * U + p * L ^ 2) - (2 - p) * U - p * L) / (3 * p - 2)
3) b0=(Sqrt(p) * Sqrt((4 - 2 * p) * U ^ 2 + (4 - 2 * p) * L * U + p * L ^ 2) + (p - 2) * U - p * L) / (3 * p - 2)

but that's not enough, and sometimes I don't end up with the good root. I guess I should find all real solutions with some methods for finding polynomial roots. Can you point me to a method which a) is simple to code (I'm not exactly an expert in the field, you know), b) is guaranteed to work (since this must go in an automatic process) and c) doesn't rely on complex arithmetics? VBA doesn't support that. However, I could try to implement a complex number class, if that's strictly necessary.

Of course, if you think that I'm completely off and that there's a far better way to solve the initial problem, feel free to show me :)

Best Regards

Sergio Rossi




Point your RSS reader here for a feed of the latest messages in this topic.

[Privacy Policy] [Terms of Use]

© Drexel University 1994-2014. All Rights Reserved.
The Math Forum is a research and educational enterprise of the Drexel University School of Education.