Real Solving on Bivariate Systems with Sturm Sequences and SLV Maple - - PowerPoint PPT Presentation

Real Solving on Bivariate Systems with Sturm Sequences and SLV MapleTM library

Dimitris Diochnos

University of Illinois at Chicago

• Dept. of Mathematics, Statistics, and Computer Science

September 27, 2007

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 1 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 2 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 3 / 65

Notation and Conventions

Complexity:

◮

• OB implies that we are ignoring (poly-)logarithmic factors.

Length function: L()

◮ Given ν ∈ Z, L(ν) implies the bitsize of integer ν. ◮ Given A ∈ Z[x] L(A) implies the maximum bitsize of the

coefficients.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 4 / 65

Operations on Lists.

Sign Variations: Given a list of signs compute sign-swaps. Ignore zeros. Example VAR([+, +, −, 0, −, 0, 0, +]) = 2 Intermediate Points: Given a list of (sorted) rational numbers compute rationals in between. Compute two more bounding rationals for the entire sequence.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 5 / 65

Polynomial GCD Computation

Euclid’s Algorithm.

◮ Works fine when F, G ∈ Q[x].

What if we want to work in Z[x]?

◮ Pseudo-divisions are required.

k · F = Q · G + λ · R where F, G, Q, R ∈ Z[x] and k, λ ∈ Z.

Pseudo-Euclidean:    (k, λ) = (lead (G)δ+1 , 1) R = rem (F, G) Q = quo (F, G)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 6 / 65

Polynomial GCD Computation

Euclid’s Algorithm.

◮ Works fine when F, G ∈ Q[x].

What if we want to work in Z[x]?

◮ Pseudo-divisions are required.

k · F = Q · G + λ · R where F, G, Q, R ∈ Z[x] and k, λ ∈ Z.

Pseudo-Euclidean:    (k, λ) = (lead (G)δ+1 , 1) R = rem (F, G) Q = quo (F, G)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 6 / 65

Polynomial GCD Computation

Euclid’s Algorithm.

◮ Works fine when F, G ∈ Q[x].

What if we want to work in Z[x]?

◮ Pseudo-divisions are required.

k · F = Q · G + λ · R where F, G, Q, R ∈ Z[x] and k, λ ∈ Z.

Pseudo-Euclidean:    (k, λ) = (lead (G)δ+1 , 1) R = rem (F, G) Q = quo (F, G)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 6 / 65

Sturm Sequences and Signed PRSs.

Corollary

Every sequence (Ai) = (A, A′, . . .) where λiAi = kiAi−2 + Ai−1Qi−1 where ki, λi ∈ R, kiλi < 0 and A1 = A square-free is a Sturm sequence in [Il, IR] where A(IL)A(IR) = 0.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 7 / 65

Definitions on Signed PRSs.

λiAi = kiAi−2 + Ai−1Qi−1

Pseudo-Euclidean: L(Ai) = O((1 + √ 2)i) ⇒ Impractical. Primitive-Part: 8 < : (ki, λi) = (lead (Ai−1)δi +1 , content (prem (Ai, Ai−1))) Ai = rem (λiAi−2, Ai−1) /ki Qi = quo (λiAi−2, Ai−1) Time Bound: e OB(p2q2τ). Output Bound: L(Ai) = O(pτ). SubResultant: ˘ Ai = prem (Ai−2, Ai−1) / lead (Ai−2) Time Bound: e OB(p2qτ). Output Bound: L(Ai) = O(pτ). Sturm-Habicht: Similar to SubResultant sequences.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 8 / 65

Example on PRSs.

Given: f = x8 + x6 − 3x4 − 3x3 + 8x2 + 2x − 5 g = 3x6 + 5x4 − 4x2 − 9x + 21 Euclidean: 15x4 − 3x2 + 9 15795x2 + 30375x − 59535 1254542875143750x − 1654608338437500 −12593338795500743100931141992187500 Primitive Part: 5x4 − x2 + 3 13x2 + 25x − 49 4663x − 6150 −1 SubResultant: 15x4 − 3x2 + 9 65x2 + 125x − 245 9326x − 12300 −260708

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 9 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials Bounds Real Algebraic Numbers Solving in One Variable

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 10 / 65

Bounding Roots

Let ̺ ∈ C be a root of the polynomial A(x) = adxd + · · · + a1x + a0. Cauchy, Mignotte: |̺| ≤ 1 + max{|a0|, |a1|, . . . , |ad−1|} |ad| |̺| ≤ max

• d|ad−1|

|ad| ,

• d|ad−2|

|ad| ,

3

• d|ad−3|

|ad| , . . . ,

d

• d|a0|

|ad|

• Zassenhaus:

|̺| ≤ 2 · max

k

• d−k
• ak

ad

• Complexity:

OB(dτ).

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 11 / 65

Real Algebraic Numbers

Definition: Let A ∈ Z[x]. Each ̺ ∈ R such that A(̺) = 0 is a Real Algebraic Number and we will write ̺ ∈ Ralg. Representation: Isolating Intervals. ̺ ∼ = [ square-free(A), [ IL, IR ] ] such that IL, IR ∈ Q, IL ≤ ̺ ≤ IR and ̺ unique root in interval [IL, IR]. Basic Operations: Sign-At. Comparison.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 12 / 65

Real Algebraic Numbers

Definition: Let A ∈ Z[x]. Each ̺ ∈ R such that A(̺) = 0 is a Real Algebraic Number and we will write ̺ ∈ Ralg. Representation: Isolating Intervals. ̺ ∼ = [ square-free(A), [ IL, IR ] ] such that IL, IR ∈ Q, IL ≤ ̺ ≤ IR and ̺ unique root in interval [IL, IR]. Basic Operations: Sign-At. Comparison.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 12 / 65

Real Algebraic Numbers

Definition: Let A ∈ Z[x]. Each ̺ ∈ R such that A(̺) = 0 is a Real Algebraic Number and we will write ̺ ∈ Ralg. Representation: Isolating Intervals. ̺ ∼ = [ square-free(A), [ IL, IR ] ] such that IL, IR ∈ Q, IL ≤ ̺ ≤ IR and ̺ unique root in interval [IL, IR]. Basic Operations: Sign-At. Comparison.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 12 / 65

Polynomial Evaluation at a Ralg point.

Input: α = [A, [L, R]] and f ∈ R[x]; A, f square-free Output: sign(f(α))

1

Compute SubResultant sequence.

2

Evaluate on endpoints.

3

Yield result with Sign-Variations. Complexity: OB(pqτ + p min{p, q}2τ).

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 13 / 65

Comparison between Ralg numbers.

Input: α = [A, [L, R]] and β = [B, [L, R]] Output: Decide α ≤ β. Idea: Compute SIGN-AT(A(β)). Note: We know the sign of A′(α). Complexity: OB(pqτ + p min{p, q}2τ).

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 14 / 65

Solving in One Variable

Input: A ∈ Z[x], square-free. Output: A list S of intervals that contain the real roots of A. Subdivision Method. Complexity: [Emiris,Mourrain,Tsigaridas - 2006] Time: OB(p6 + p4τ 2) (+ multiplicities) Output: O(pτ)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 15 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials Extensions. Resultant

4

Real Solving on Bivariate Systems

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 16 / 65

Bivariate SignAt.

Input: α ∼ = [A, [AL, AR]], β ∼ = [B, [BL, BR]] and f ∈ Z[x, y]. Output: Compute the sign of f(α, β).

1

Compute a Sturm sequence of (A, f).

2

Evaluate the sequence on each of endpoints AL, AR of α.

3

Perform SIGN-AT at y = β on each polynomial in sequence.

4

Count sign-variations on sequences.

5

Yield result with (VL − VR) · A′(α) Complexity: Sign can be computed in OB(n2

1n3 2τ).

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 17 / 65

Working in Extension-Fields.

GCD Computation [Hoeij and Monagan - 2002] Other required operations are reduced to univariate SIGN-AT.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 18 / 65

Resultant

Theorem

Given f, g ∈ Z[x] such that f = n

i=1 aixx and g = m j=1 bjxj with

anbm = 0, there is a unique (up to sign) irreducible polynomial res(f, g) ∈ Z[an, . . . , a0, bm, . . . , b0] which is zero iff f, g have a common factor. It is homogeneous and deg(res(f, g)) = deg(f) + deg(g) = n + m. This polynomial is called resultant.

Theorem

Given f, g ∈ K[x] we can compute the resultant res(f, g) via the Sylvester matrix Syl(f, g). More specifically, we have: res(f, g) = det(Syl(f, g)). Complexity: OB(nmτ).

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 19 / 65

Resultant

Definition (Sylvester Matrix)

The Sylvester matrix of f and g is the (n + m) × (n + m) matrix defined as follows: Syl(f, g) =                  an an−1 an−2 . . . . . . a0 an an−1 an−2 . . . . . . a0 ... ... ... ... ... ... an an−1 an−2 . . . . . . a0 bm bm−1 bm−2 . . . b0 bm bm−1 bm−2 . . . b0 ... ... ... ... ... ... ... ... ... ... bm bm−1 bm−2 . . . b0                 

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 20 / 65

Resultant

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 21 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems Introduction GRID M RUR G RUR

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 22 / 65

The Problem

Given f, g ∈ Z[x, y] we want to compute all real solutions of the system f = g = 0. Previous Work

◮ Assuming df, dg ≤ N and L(f), L(g) ≤ N previous complexity

bounds were: Isolating Intervals: OB(N30) [Arnon, McCallum - 1988] Thom’s Encoding: OB(N16) [Gonz´ alez-Vega, El Kahoui - 1996]

Our Work [Diochnos, Emiris, Tsigaridas - 2007]

◮ 3 Algorithms using projection and Isolating Intervals ⋆ GRID, M RUR, G RUR ⋆ Algorithms differentiate on matching ◮ Running Times:

OB(N14) and OB(N12)

⋆ Input Size: e

OB(N3)

⋆ Output Size: e

OB(N4)

⋆ Projection Complexity: e

OB(N12)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 23 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems Introduction GRID M RUR G RUR

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 24 / 65

GRID Algorithm

1

Project solutions on axes and find roots.

2

Generate all candidate pairs.

3

Check with BIVARIATESIGNAT for solutions.

Theorem (GRID Complexity)

Isolating all real roots of the system using GRID has complexity

• OB(n14 + n13τ), provided that τ = O(n3)
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 25 / 65

GRID Example

Assume we want to solve the system: f(x, y) = y − x3 + 2x − 1 g(x, y) = y − x2 The system has solutions:    (+0.45, +0.20) (−1.25, +1.55) (+1.80, +3.25)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 26 / 65

GRID Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 27 / 65

GRID Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 28 / 65

GRID Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 29 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems Introduction GRID M RUR G RUR

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 30 / 65

Generic Position

[Gonz´ alez-Vega, El Kahoui and Basu, Pollack, Roy]

Definition (Generic Position)

Two polynomials f, g ∈ R[x, y] are in generic position if (f(̺, β1) = g(̺, β1) = 0) ∧ (f(̺, β2) = g(̺, β2) = 0) → β1 = β2. Leading coefficients w.r.t. y are not zero when specializing x.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 31 / 65

Solutions in Generic Position

Theorem (Solution Existence)

Let A, B square-free and co-prime polynomials in generic position. If SRj(x, y) = srj(x)yj + srj,j−1(x)yj−1 + . . . + srj,0(x), then if (̺, β) is a real solution of the system A = B = 0, then there exists k, such that sr0(̺) = . . . = srk−1(̺) = 0, srk(̺) = 0 and β = − 1

k srk,k−1(̺) srk(̺)

What we need is: The minimum k such that srk(̺) = psck(̺) = 0.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 32 / 65

How to Compute K?

Naive: Evaluate all pscj(x) at x = ̺ and pick k = minj(pscj(̺)) = 0. Better: Assume that it works for j = k > 0. It follows that pscj(̺) = 0 ∀j ∈ {0, 1, . . . , k − 1}. Hence ̺ appears on the gcd(psc0, . . . , psck−1). We define the polynomials: Φ0(x) = psc0(x) gcd(psc0(x), psc′

0(x))

Φj(x) = gcd(Φj−1(x), pscj(x)) as well as Γj = Φj−1(x) Φj(x) , j ∈ {1, . . . , n − 1} k = j + 1 Complexity: Index k can be computed in time

• OB(n8τ).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 33 / 65

How to Compute K?

Naive: Evaluate all pscj(x) at x = ̺ and pick k = minj(pscj(̺)) = 0. Better: Assume that it works for j = k > 0. It follows that pscj(̺) = 0 ∀j ∈ {0, 1, . . . , k − 1}. Hence ̺ appears on the gcd(psc0, . . . , psck−1). We define the polynomials: Φ0(x) = psc0(x) gcd(psc0(x), psc′

0(x))

Φj(x) = gcd(Φj−1(x), pscj(x)) as well as Γj = Φj−1(x) Φj(x) , j ∈ {1, . . . , n − 1} k = j + 1 Complexity: Index k can be computed in time

• OB(n8τ).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 33 / 65

How to Compute K?

Naive: Evaluate all pscj(x) at x = ̺ and pick k = minj(pscj(̺)) = 0. Better: Assume that it works for j = k > 0. It follows that pscj(̺) = 0 ∀j ∈ {0, 1, . . . , k − 1}. Hence ̺ appears on the gcd(psc0, . . . , psck−1). We define the polynomials: Φ0(x) = psc0(x) gcd(psc0(x), psc′

0(x))

Φj(x) = gcd(Φj−1(x), pscj(x)) as well as Γj = Φj−1(x) Φj(x) , j ∈ {1, . . . , n − 1} k = j + 1 Complexity: Index k can be computed in time

• OB(n8τ).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 33 / 65

How to Compute K?

Naive: Evaluate all pscj(x) at x = ̺ and pick k = minj(pscj(̺)) = 0. Better: Assume that it works for j = k > 0. It follows that pscj(̺) = 0 ∀j ∈ {0, 1, . . . , k − 1}. Hence ̺ appears on the gcd(psc0, . . . , psck−1). We define the polynomials: Φ0(x) = psc0(x) gcd(psc0(x), psc′

0(x))

Φj(x) = gcd(Φj−1(x), pscj(x)) as well as Γj = Φj−1(x) Φj(x) , j ∈ {1, . . . , n − 1} k = j + 1 Complexity: Index k can be computed in time

• OB(n8τ).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 33 / 65

How to Compute K?

Naive: Evaluate all pscj(x) at x = ̺ and pick k = minj(pscj(̺)) = 0. Better: Assume that it works for j = k > 0. It follows that pscj(̺) = 0 ∀j ∈ {0, 1, . . . , k − 1}. Hence ̺ appears on the gcd(psc0, . . . , psck−1). We define the polynomials: Φ0(x) = psc0(x) gcd(psc0(x), psc′

0(x))

Φj(x) = gcd(Φj−1(x), pscj(x)) as well as Γj = Φj−1(x) Φj(x) , j ∈ {1, . . . , n − 1} k = j + 1 Complexity: Index k can be computed in time

• OB(n8τ).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 33 / 65

Matching Solutions in M RUR - FIND.

k computation yields first k such that psck(x) = 0 at x = ̺. SRk(x = ̺) ∈ (Z[̺])[y] has k roots. Generic Position ⇒ At most 1 root. Hence 1 root with multiplicity k of the form: R(̺) = −1 k srk,k−1(̺) psck(̺) . Binary search for the unique interval such that: Qj < R(̺) < Qj+1 Complexity: The matching procedure takes time OB(n7τ) for each candidate on x-axis.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 34 / 65

Matching Solutions in M RUR - FIND.

k computation yields first k such that psck(x) = 0 at x = ̺. SRk(x = ̺) ∈ (Z[̺])[y] has k roots. Generic Position ⇒ At most 1 root. Hence 1 root with multiplicity k of the form: R(̺) = −1 k srk,k−1(̺) psck(̺) . Binary search for the unique interval such that: Qj < R(̺) < Qj+1 Complexity: The matching procedure takes time OB(n7τ) for each candidate on x-axis.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 34 / 65

Matching Solutions in M RUR - FIND.

k computation yields first k such that psck(x) = 0 at x = ̺. SRk(x = ̺) ∈ (Z[̺])[y] has k roots. Generic Position ⇒ At most 1 root. Hence 1 root with multiplicity k of the form: R(̺) = −1 k srk,k−1(̺) psck(̺) . Binary search for the unique interval such that: Qj < R(̺) < Qj+1 Complexity: The matching procedure takes time OB(n7τ) for each candidate on x-axis.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 34 / 65

Matching Solutions in M RUR - FIND.

k computation yields first k such that psck(x) = 0 at x = ̺. SRk(x = ̺) ∈ (Z[̺])[y] has k roots. Generic Position ⇒ At most 1 root. Hence 1 root with multiplicity k of the form: R(̺) = −1 k srk,k−1(̺) psck(̺) . Binary search for the unique interval such that: Qj < R(̺) < Qj+1 Complexity: The matching procedure takes time OB(n7τ) for each candidate on x-axis.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 34 / 65

Regarding Genericity

But is the algorithm general? YES!

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 35 / 65

Regarding Genericity

But is the algorithm general? YES!

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 35 / 65

Performing Shear.

We are performing the transformation (X, Y) − → (X + α · Y, Y) for some α ∈ Z. Complexity: The deterministic computation of α has complexity

• OB(n9τ).

Bad practical performance. Random shear is the choice.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 36 / 65

M RUR Algorithm

1

Project solutions on axes and find roots.

2

Find Intermediate Points on y-axis.

3

For each x solution, compute kx.

4

Find suitable interval on y-axis to match solutions.

Theorem (M RUR Complexity)

Isolating all real roots of the system using M RUR has complexity

• OB(n12 + n10τ 2).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 37 / 65

M RUR Example

Assume we want to solve the system: f(x, y) = y − x3 + 2x − 1 g(x, y) = y − x2 The system is not in generic position. The shear (x, y) − → (x + 3y, y) makes it generic: f(x, y) = −27y3 − 27xy2 + (7 − 9x2)y − x3 + 2x − 1 g(x, y) = −9y2 + (1 − 6x)y − x2 The sheared system has solutions:    (−0.15, +0.20) (−5.91, +1.55) (−7.94, +3.25)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 38 / 65

M RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 39 / 65

M RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 40 / 65

M RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 41 / 65

M RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 42 / 65

M RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 43 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems Introduction GRID M RUR G RUR

5

Implementation

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 44 / 65

G RUR Algorithm

1

Project solutions on axes and find roots.

2

Find Intermediate Points on y-axis.

3

For each candidate solution ̺ on x-axis:

◮ Compute H(y) = gcd(˜

f(̺, y), ˜ g(̺, y)) ∈ (Z[̺])[y].

◮ Check for solutions of H(y) on candidate intervals along the y-axis.

Theorem (G RUR Complexity)

Isolating all real roots of the system using G RUR has complexity

• OB(n12 + n10τ 2).
• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 45 / 65

G RUR Example

Assume we want to solve the system: f(x, y) = y − x3 + 2x − 1 g(x, y) = y − x2 The system has solutions:    (+0.45, +0.20) (−1.25, +1.55) (+1.80, +3.25)

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 46 / 65

G RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 47 / 65

G RUR Example

Assume ̺ ∈ Ralg such that resy(f, g)(̺) = 0. We want to compute H = gcd(f(̺, y), g(̺, y)) ∈ (Z[̺])[y]. Since resy(f, g)(α) = 0 we have ̺3 − ̺2 − 2̺ + 1 = 0. Or equivalently: ̺3 = ̺2 + 2̺ − 1 Hence f(̺, y) = y − ̺3 + 2̺ − 1 = y − ̺2 = g(̺, y). As a result gcd(f(̺, y), g(̺, y)) = y − ̺2.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 48 / 65

G RUR Example

Assume ̺ ∈ Ralg such that resy(f, g)(̺) = 0. We want to compute H = gcd(f(̺, y), g(̺, y)) ∈ (Z[̺])[y]. Since resy(f, g)(α) = 0 we have ̺3 − ̺2 − 2̺ + 1 = 0. Or equivalently: ̺3 = ̺2 + 2̺ − 1 Hence f(̺, y) = y − ̺3 + 2̺ − 1 = y − ̺2 = g(̺, y). As a result gcd(f(̺, y), g(̺, y)) = y − ̺2.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 48 / 65

G RUR Example

Assume ̺ ∈ Ralg such that resy(f, g)(̺) = 0. We want to compute H = gcd(f(̺, y), g(̺, y)) ∈ (Z[̺])[y]. Since resy(f, g)(α) = 0 we have ̺3 − ̺2 − 2̺ + 1 = 0. Or equivalently: ̺3 = ̺2 + 2̺ − 1 Hence f(̺, y) = y − ̺3 + 2̺ − 1 = y − ̺2 = g(̺, y). As a result gcd(f(̺, y), g(̺, y)) = y − ̺2.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 48 / 65

G RUR Example

Assume ̺ ∈ Ralg such that resy(f, g)(̺) = 0. We want to compute H = gcd(f(̺, y), g(̺, y)) ∈ (Z[̺])[y]. Since resy(f, g)(α) = 0 we have ̺3 − ̺2 − 2̺ + 1 = 0. Or equivalently: ̺3 = ̺2 + 2̺ − 1 Hence f(̺, y) = y − ̺3 + 2̺ − 1 = y − ̺2 = g(̺, y). As a result gcd(f(̺, y), g(̺, y)) = y − ̺2.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 48 / 65

G RUR Example

The solutions of the system come from the H = gcd ∈ (Z[̺])[y]. Consequently, by substituting y with Intermediate Points H must change sign at x = ̺. The Intermediate Points are: Q1 = −1 Q2 = 163

128

Q3 = 73

32

Q4 = 5 Checking interval [Q1, Q2] : H1(x) = −1 − x2 H2(x) = +163 − 128x2 Note that H1(x) < 0 ∀x ∈ R. Hence, we need to find a root ̺ ∈ Ralg such that H2(̺) > 0.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 49 / 65

G RUR Example

The solutions of the system come from the H = gcd ∈ (Z[̺])[y]. Consequently, by substituting y with Intermediate Points H must change sign at x = ̺. The Intermediate Points are: Q1 = −1 Q2 = 163

128

Q3 = 73

32

Q4 = 5 Checking interval [Q1, Q2] : H1(x) = −1 − x2 H2(x) = +163 − 128x2 Note that H1(x) < 0 ∀x ∈ R. Hence, we need to find a root ̺ ∈ Ralg such that H2(̺) > 0.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 49 / 65

G RUR Example

The solutions of the system come from the H = gcd ∈ (Z[̺])[y]. Consequently, by substituting y with Intermediate Points H must change sign at x = ̺. The Intermediate Points are: Q1 = −1 Q2 = 163

128

Q3 = 73

32

Q4 = 5 Checking interval [Q1, Q2] : H1(x) = −1 − x2 H2(x) = +163 − 128x2 Note that H1(x) < 0 ∀x ∈ R. Hence, we need to find a root ̺ ∈ Ralg such that H2(̺) > 0.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 49 / 65

G RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 50 / 65

G RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 51 / 65

G RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 52 / 65

G RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 53 / 65

G RUR Example

Working similarly . . .

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 54 / 65

G RUR Example

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 55 / 65

Outline

1

Introduction

2

Results in Univariate Polynomials

3

Results in Bivariate Polynomials

4

Real Solving on Bivariate Systems

5

Implementation SLV Library and Sample Usage In-Depth View

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 56 / 65

Sturm soLVer and Sample Usage

Open Source MAPLE Code.

◮ Homepage:

http://www.di.uoa.gr/˜erga/soft/SLV index.html

◮ 2700 lines of source code. ◮ Real Algebraic Numbers → Isolating Intervals Representation. ⋆ Univariate Sign-At - Sturm Sequences ⋆ Bivariate Sign-At - Sturm Sequences ◮ Real Solving on Univariate Polynomials ◮ Real Solving on Bivariate Systems

Sample Usage:

f := 1 + 2*x + xˆ2*y - 5*x*y + xˆ2: g := 2*x + y - 3: bivsols := SLV:-solveGRID ( f, g ): SLV:-display_2 ( bivsols ); < 2*xˆ2-12*x+1, [ 3, 7], 5.915475965 > , < xˆ2+6*x-25, [ -2263/256, -35/4], -8.830718995 > < x-1, [ 1, 1], 1 > , < x-1, [ 1, 1], 1 > < 2*xˆ2-12*x+1, [3/64, 3/32], .8452400565e-1 > , < xˆ2+6*x-25, [23179/8192, 2899/1024], 2.830943108 >

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 57 / 65

GCD in Extension Field, Filtering and Speedup.

GCD in Extension Field. [Hoeij and Monagan - 2002] Filtering

◮ Interval and floating point arithmetic. ⋆ Quadratic Interval Refinement [Abbott - 2006] ⋆ Iterations based on total degree of input polynomials. ◮ GCD computation. ◮ Finally exact algorithms and computation. ◮ M RUR pre-computation filtering.

Speedup Real Solving:

◮ GRID −

→ 4×

◮ M RUR −

→ 10×

◮ G RUR −

→ 1.1×

+ multiplicities:

◮ GRID and M RUR −

→ 10×

◮ G RUR −

→ 2×

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 58 / 65

Analyzing SLV algorithms

deg ≤ 5 → similar performance in most cases As the degree increases: Real Solving:

GRID GRUR = 7-10 and MRUR GRUR ≃ 38

+ multiplicities:

GRID GRUR = 60 and MRUR GRUR ≃ 12

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 59 / 65

SLV: Performance Fragments

phase of the median mean algorithm

GRID

projections 00.04 00.08

• univ. solving

07.08 26.77

• biv. solving

96.18 73.03 sorting 00.06 00.12

MRUR

projection 00.06 00.14

• univ. solving

15.55 17.47 StHa seq. 01.17 05.80

• inter. points

00.09 00.32 filter x-cand 26.68 23.81 compute K 02.04 07.06

• biv. solving

51.17 45.41

GRUR

projections 00.23 00.48

• univ. solving

39.83 41.68

• inter. points

00.54 01.11 rational biv. 14.83 15.89 Ralg biv. 65.30 40.53 sorting 00.22 00.32

Table: Procedures at a glance.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 60 / 65

SLV: Running Times

system deg R Average Time (msecs) f g sols

GRID MRUR GRUR

R1 3 4 2 5 9 5 R2 3 1 1 66 21 36 R3 3 1 1 1 2 1 M1 3 3 4 87 72 10 M2 4 2 3 4 5 4 M3 6 3 5 803 782 110 M4 9 10 2 218 389 210 D1 4 5 1 6 12 6 D2 2 2 4 667 147 128 C1 7 6 6 1, 896 954 222 C2 4 3 6 177 234 18 C3 8 7 13 580 1, 815 75 C4 8 7 17 5, 903 80, 650 370 C5 16 15 17 > 20′ 60, 832 3, 877 W1 7 6 9 2, 293 2, 115 247 W2 4 3 5 367 283 114 W3 8 7 13 518 2, 333 24 W4 8 7 17 5, 410 77, 207 280

Table: Running Times are averages over 10 runs.

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 61 / 65

Comparative Performance

Real Solving Bivariate Systems

• vs. GBRS

Performance: G RUR Similar performance with GBRS on average

SLV increased fluctuation vs. uniform treatment in GBRS

Robustness: SLV → reliable kernel on loops rs isolate → may cause problems with MAPLE kernel

• vs. SYNAPS (STURM, SUBDIV and NEWMAC)

Performance:

◮ deg ≤ 5 → SYNAPS slightly faster ◮ deg ↑ → SYNAPS slower or incomplete solution set ◮

STURM → determinants compute square-free part

◮

SUBDIV and NEWMAC → double precision

⋆ Precision is limited ⋆ May lose solutions ⋆ May generate false solutions

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 62 / 65

Comparative Performance

Topology

• vs. INSULATE [R. Seidel, N. Wolpert]

Performance:

◮ deg ≤ 4 → INSULATE faster. No more than 2× ◮ deg ↑→ SLV faster. ⋆

G RUR → 40× faster when deg ∼ 16

⋆

M RUR → 2-3× faster when deg ∼ 16

• vs. TOP [L. Gonz´

alez-Vega, I. Necula] Performance:

◮ No theory for good settings on the extra required parameter of TOP ◮ Critical Points ↑ → G RUR faster even for good settings for TOP ◮

G RUR 7-20× faster than TOP (param. = 60, 500)

Quality: TOP may produce incorrect results under bad initialization

◮ May lose solutions ◮ May generate false solutions

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 63 / 65

Summary

Three algorithms for real solving of bivariate polynomial systems:

◮ GRID →

OB(N14).

◮ M RUR and G RUR →

OB(N12).

The topology of a real plane algebraic curve can be computed in

• OB(N12).

Novel and robust MAPLE implementation with promising results. Outlook

◮ Express solutions of the sheared system in the original coordinate

system.

◮ Better implementation on PRSs. ◮ Fan-In / Fan-Out

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 64 / 65

THANK YOU!

• D. I. Diochnos (UIC - MSCS)

Real Solving on Bivariate Algebraic Systems Sep ’07 65 / 65