| 624 | } |
| 625 | |
| 626 | int quadraticSolve(const poly p, number &s1, number &s2, |
| 627 | const number tolerance) |
| 628 | { |
| 629 | poly q = pCopy(p); |
| 630 | int result; |
| 631 | |
| 632 | if (q == NULL) result = -1; |
| 633 | else |
| 634 | { |
| 635 | int degree = pGetExp(q, 1); |
| 636 | if (degree == 0) result = 0; /* constant polynomial <> 0 */ |
| 637 | else |
| 638 | { |
| 639 | number c2 = nInit(0); /* coefficient of var(1)^2 */ |
| 640 | number c1 = nInit(0); /* coefficient of var(1)^1 */ |
| 641 | number c0 = nInit(0); /* coefficient of var(1)^0 */ |
| 642 | if (pGetExp(q, 1) == 2) |
| 643 | { nDelete(&c2); c2 = nCopy(pGetCoeff(q)); q = q->next; } |
| 644 | if ((q != NULL) && (pGetExp(q, 1) == 1)) |
| 645 | { nDelete(&c1); c1 = nCopy(pGetCoeff(q)); q = q->next; } |
| 646 | if ((q != NULL) && (pGetExp(q, 1) == 0)) |
| 647 | { nDelete(&c0); c0 = nCopy(pGetCoeff(q)); q = q->next; } |
| 648 | |
| 649 | if (degree == 1) |
| 650 | { |
| 651 | c0 = nInpNeg(c0); |
| 652 | s1 = nDiv(c0, c1); |
| 653 | result = 1; |
| 654 | } |
| 655 | else |
| 656 | { |
| 657 | number tmp = nMult(c0, c2); |
| 658 | number tmp2 = nAdd(tmp, tmp); nDelete(&tmp); |
| 659 | number tmp4 = nAdd(tmp2, tmp2); nDelete(&tmp2); |
| 660 | number discr = nSub(nMult(c1, c1), tmp4); nDelete(&tmp4); |
| 661 | if (nIsZero(discr)) |
| 662 | { |
| 663 | tmp = nAdd(c2, c2); |
| 664 | s1 = nDiv(c1, tmp); nDelete(&tmp); |
| 665 | s1 = nInpNeg(s1); |
| 666 | result = 2; |
| 667 | } |
| 668 | else if (nGreaterZero(discr)) |
| 669 | { |
| 670 | realSqrt(discr, tolerance, tmp); /* sqrt of the discriminant */ |
| 671 | tmp2 = nSub(tmp, c1); |
| 672 | tmp4 = nAdd(c2, c2); |
| 673 | s1 = nDiv(tmp2, tmp4); nDelete(&tmp2); |
| 674 | tmp = nInpNeg(tmp); |
| 675 | tmp2 = nSub(tmp, c1); nDelete(&tmp); |
| 676 | s2 = nDiv(tmp2, tmp4); nDelete(&tmp2); nDelete(&tmp4); |
| 677 | result = 3; |
| 678 | } |
| 679 | else |
| 680 | { |
| 681 | discr = nInpNeg(discr); |
| 682 | realSqrt(discr, tolerance, tmp); /* sqrt of |discriminant| */ |
| 683 | tmp2 = nAdd(c2, c2); |
no test coverage detected