Solving general polynomial systems

solve(args...; options...)::Result

The solve function takes many different arguments and options depending on your specific situation, but in the it always returns a Result containing the result of the computations. In the following we show the different inputs solve takes and at the end we list all possible options.

Total Degree Homotopy

solve(F; options...)

Solve the system F using a start system computed from the degrees of the entries of F. The number of paths to track is equal to the total degree d1⋯dn, where di is the degree of the ith entry of F. F can be

• Vector{<:MultivariatePolynomials.AbstractPolynomial} (e.g. constructed by @polyvar)
• A composition of polynomial systems constructed by compose.
• AbstractSystem (the system has to represent a homogeneous polynomial system.)

Example

We can solve the system $F(x,y) = (x^2+y^2+1, 2x+3y-1)$ in the following way:

julia> @polyvar x y;
julia> solve([x^2+y^2+1, 2x+3y-1])
Result with 2 solutions
==================================
• 2 non-singular solutions (0 real)
• 0 singular solutions (0 real)
• 2 paths tracked
• random seed: 661766

Polyhedral Homotopy

solve(F; start_system = :polyhedral, only_torus=false, options...)

Solve the system F using a start system computed from the Newton Polytopes of the entries F. The number of paths to track is equal to the mixed volume of the Newton Polytopes of the entries of F. The mixed volume is at most the total degree of F. F can be

• Vector{<:MultivariatePolynomials.AbstractPolynomial} (e.g. constructed by @polyvar)
• A composition of polynomial systems constructed by compose.
• AbstractSystem (the system has to represent a homogeneous polynomial system.)

If only_torus == true then only solutions in the algebraic torus $(ℂ\setminus \{0\})^n$ will be computed.

Example

We can solve the system $F(x,y) = (x^2+y^2+1, 2x+3y-1)$ in the following way:

julia> @polyvar x y;
julia> solve([x^2+y^2+1, 2x+3y-1]; start_system = :polyhedral)
Result with 2 solutions
==================================
• 2 non-singular solutions (0 real)
• 0 singular solutions (0 real)
• 2 paths tracked
• random seed: 222880

Homogeneous Systems

If F has is homogeneous, we return results in projective space

Examples

julia> @polyvar x y z;
julia> solve([x^2+y^2+z^2, 2x+3y-z])
Result with 2 solutions
==================================
• 2 non-singular solutions (0 real)
• 0 singular solutions (0 real)
• 2 paths tracked
• random seed: 291729

If your polynomial system is already homogeneous, but you would like to consider it as an affine system you can do

@polyvar x y z
solve([x^2+y^2+z^2, 2x+3y-z], homvar=z)

This yields the same result as solve([x^2+y^2+1, 2x+3y-1]).

Multihomogeneous Systems

By exploiting the multi-homogenous structure of a polynomial system it is possible to decrease the number of paths necessary to track.

@polyvar x y
# Use variable groups to only track 2 paths instead of 4
solve([x*y - 6, x^2 - 5], variable_groups=[(x,), (y,)])

To check whether a certain variable grouping is beneficial you can use the bezout_number function.

Start Target Homotopy

solve(G, F, start_solutions; options...)

This constructs the homotopy $H(x,t) = tG(x)+(1-t)F(x)$ to compute solutions of the system F. start_solutions is a list of solutions of G which are tracked to solutions of F.

Example

@polyvar x y
G = [x^2-1,y-1]
F = [x^2+y^2+z^2, 2x+3y-z]
solve(G, F, [[1, 1], [-1, 1]])

Parameter Homotopy

solve(F, startsolutions; parameters, start_parameters, target_parameters, start_gamma=nothing, target_gamma=nothing)

Solve the parameter homotopy

where $p₁$ (=start_parameters) and $p₀$ (=target_parameters) are vectors of parameter values for $F$ and $γ₁$ (=start_gamma) and $γ₀$ (=target_gamma) are complex numbers. If start_parameters or target_parameters is nothing, it is assumed that γ₁ and γ₀ are $1$. The input parameters specifies the variables of F which should be considered as parameters. Necessarily we have length(parameters) == length(p₁) == length(p₀).

solve(F, startsolutions; parameters, p₁, p₀, γ₁=nothing, γ₀=nothing)

This is a unicode variant where γ₁=start_parameters, γ₀=target_parameters, γ₁=start_gamma, γ₀=target_gamma.

Example

We want to solve a parameter homotopy $H(x,t) := F(x; t[1, 0]+(1-t)[2, 4])$ where

and let's say we are only intersted in tracking of $[1,1]$. This can be accomplished as follows

@polyvar x[1:2] a[1:2]
F = [x[1]^2-a[1], x[1]*x[2]-a[1]+a[2]]
startsolutions = [[1, 1]]
p₁ = [1, 0]
p₀ = [3im, 0.5+2im]
solve(F, startsolutions; parameters=a, start_parameters=p₁, target_parameters=p₀)
# If you like unicode this is also possible
solve(F, startsolutions; parameters=a, p₁=p₁, p₀=p₀)

Abstract Homotopy

solve(H::AbstractHomotopy, start_solutions; options...)

Solve the homotopy H by tracking the each solution of $H(⋅, t)$ (as provided by start_solutions) from $t=1$ to $t=0$. Note that H has to be a homotopy between homogeneous polynomial systems. If it should be considered as an affine system indicate which is the index of the homogenization variable, e.g. solve(H, startsolutions, homvar=3) if the third variable is the homogenization variable.

Options

General options:

• seed::Int: The random seed used during the computations.
• show_progress=true: Whether a progress bar should be printed to standard out.
• threading=true: Enable or disable multi-threading.
• path_result_details=:default: The amount of information computed in each path result. Possible values are :minimal (minimal details), :default (default) and :extensive (all information possible).
• homvar::Union{Int,MultivariatePolynomials.AbstractVariable}: This considers the homogeneous system F as an affine system which was homogenized by homvar. If F is an AbstractSystem homvar is the index (i.e. Int) of the homogenization variable. If F is an AbstractVariables (e.g. created by @polyvar x) homvar is the actual variable used in the system F.
• system::AbstractSystem: A constructor to assemble a AbstractSystem. The default is SPSystem. This constructor is only applied to the input of solve. The constructor is called with system(polynomials, variables) where polynomials is a vector of MultivariatePolynomials.AbstractPolynomials and variables determines the variable ordering. If you experience significant compilation times, consider to change system to FPSystem.
• homotopy::AbstractHomotopy: A constructor to construct a AbstractHomotopy for the totaldegree and start target homotopy. The default is StraightLineHomotopy. The constructor is called with homotopy(start, target) where start and target are homogeneous AbstractSystems.
• affine_tracking::Bool=true: Indicate whether path tracking should happen in affine space.
• projective_tracking::Bool=false: Indicate whether path tracking should happen in projective space. The flag affine_tracking is dominant.
• path_jumping_check::Bool=true: Enable a check whether one of the paths jumped to another one.

Path tracking specific options:

• accuracy=1e-7: The accuracy required during the path tracking.
• corrector::AbstractCorrector: The corrector used during in the predictor-corrector scheme. The default is NewtonCorrector.
• initial_step_size=0.1: The size of the first step.
• max_corrector_iters=2: The maximal number of correction steps in a single step.
• max_lost_digits::Real: The tracking is terminated if we estimate that we loose more than max_lost_digits in the linear algebra steps. This threshold depends on the precision argument.
• max_refinement_iters=5: The maximal number of correction steps used to refine the final value.
• max_steps=1_000: The maximal number of iterations the path tracker has available. Note that this changes to 10_000 for parameter homotopies.
• max_step_size=Inf: The maximal step size.
• min_step_size=1e-14: The minimal step size.
• precision::PrecisionOption=PRECISION_FIXED_64: The precision used for evaluating the residual in Newton's method.
• predictor::AbstractPredictor: The predictor used during in the predictor-corrector scheme. The default is Heun()`.
• refinement_accuracy=1e-8: The precision required for the final value.
• simple_step_size_alg=false: Use a more simple step size algorithm.
• steps_jacobian_info_update::Int=1: Every n-th step a linear system will be solved using a QR factorization to obtain an estimate for the condition number of the Jacobian.
• terminate_ill_conditioned::Bool=true: Indicates whether the path tracking should be terminated for ill-conditioned paths. A path is considerd ill-conditioned if the condition number of the Jacobian is larger than ≈1e14 or if it is larger than 1emax_lost_digits.

Endgame specific options:

• accuracy_eg::Float64=min(accuracy, 1e-5)): It is possible to change the accuracy during the path tracking. Usually you want lower the accuracy.
• cond_eg_start::Float64=1e4: The endgame is only started if the condition of the Jacobian is larger than this threshold.
• max_winding_number::Int=12: This limits the maximal number of loops taken in applying Cauchy's formula.
• min_cond_at_infinity::Float64=1e7: A path is declared as going to infinity only if it's Jacobian is also larger than this threshold.
• samples_per_loop::Int=12: To compute singular solutions Cauchy's integral formula is used. The accuracy of the solutions increases with the number of samples per loop.
• t_eg_start::Float64=0.1: The endgame starts only if t is smaller than this threshold.
• tol_val_inf_accurate::Float64=1e-4: A valuation which would result in a path declared as going to infinity is only accepted if the estimated accuracy of the valuation is less than this threshold.
• tol_val_finite_accurate::Float64=1e-3: A valuation which would result in a proper solution is only accepted if the estimated accuracy of the valuation is less than this threshold. This is only affects solutions where the path has at some point near 0 a condition number larger than cond_eg_start.

It is recommended to also take a look at the PathTracker documentation for some context.

Overdetermined system specific options:

• overdetermined_min_accuracy=1e-5: The minimal accuracy a non-singular solution needs to have to be considered a solution of the original system.
• overdetermined_min_residual=1e-3: The minimal residual a singular solution needs to have to be considered a solution of the original system.

Result{V<:AbstractVector}

The result of solve. This is a wrapper around the results of each single path (PathResult) and it contains some additional informations like a random seed to replicate the result.

seed(result)

The random seed used in the computation.

nonsingular(result::Results; conditions...)

Return all PathResults for which the solution is non-singular. This is just a shorthand for results(R; only_nonsingular=true, conditions...). For the possible conditions see results.

singular(R::Results; tol=1e10, multiple_results=false, kwargs...)

Return all PathResults for which the solution is singular. A solution is labeled singular if the condition number is greater than singular_tol, or if the winding number is > 1. If multiple_results=false only one point from each cluster of multiple solutions is returned. If If multiple_results=true all singular solutions in R are returned. For the possible kwargs see results.

results(result; only_real=false, real_tol=1e-6, only_nonsingular=false,
            onlysigular=false, singular_tol=1e10, onlyfinite=true, multiple_results=false)

Return all PathResults for which the given conditions apply.

Example

R = solve(F)

# This gives us all PathResults considered non-singular and real (but still as a complex vector).
real_solutions = results(R, only_real=true, only_nonsingular=true)

mapresults(f::Function, result; conditions...)

Apply the function f to all PathResults for which the given conditions apply. For the possible conditions see results.

Example

# This gives us all solutions considered real (but still as a complex vector).
real_solutions = mapresults(solution, R, only_real=true)

mapresults(f, result::MonodromyResult; only_real=false, real_tol=1e-6)

Apply the function f to all entries of MonodromyResult for which the given conditions apply.

Example

# This gives us all solutions considered real (but still as a complex vector).
real_solutions = mapresults(solution, R, only_real=true)

solutions(result; conditions...)

Return all solution (as Vectors) for which the given conditions apply. For the possible conditions see results.

Example

julia> @polyvar x y
julia> result = solve([(x-2)y, y+x+3]);
julia> solutions(result)
[[2.0+0.0im, -5.0+0.0im], [-3.0+0.0im, 0.0+0.0im]]

solutions(result::MonodromyResult; only_real=false, real_tol=1e-6)

Return all solutions (as SVectors) for which the given conditions apply.

Example

real_solutions = solutions(R, only_real=true)

solutions(MS::MonodromySolver)

Get the solutions of the loop.

real_solutions(result; tol=1e-6, conditions...)

Return all real solution (as Vectors of reals) for which the given conditions apply. For the possible conditions see results. Note that only_real is always true and real_tol is now tol.

Example

julia> @polyvar x y
julia> result = solve([(x-2)y, y+x+3]);
julia> real_solutions(result)
[[2.0, -5.0], [-3.0, 0.0]]

real_solutions(res::MonodromyResult; tol=1e-6)

Returns the solutions of res whose imaginary part has norm less than 1e-6.

finite(result::AffineResults; conditions...)

Return all PathResults for which the solution is finite. This is just a shorthand for results(R; onlyfinite=true, conditions...). For the possible conditions see results.

real(result, tol=1e-6)

Get all results where the solutions are real with the given tolerance tol. See is_real for details regarding the determination of 'realness'.

at_infinity(result::AffineResult)

Get all results where the solutions is at infinity.

failed(result)

Get all results where the path tracking failed.

multiplicities!(result::Result; tol=1e-6)

Compute the multiplicities of the solutions in result with respect to the given tolerance.

multiplicities(V::Results; tol=1e-6)

Returns a Vector of Vector{PathResult}s grouping the PathResults whose solutions appear with multiplicities greater 1 in 'V'. Two solutions are regarded as equal, when their pairwise distance is less than 'tol'.

nresults(result; only_real=false, real_tol=1e-6, only_nonsingular=false, singular_tol=1e10, onlyfinite=true)

The number of solutions which satisfy the corresponding predicates.

Example

result = solve(F)
# Get all non-singular results where all imaginary parts are smaller than 1e-8
nresults(result, only_real=true, real_tol=1e-8, only_nonsingular=true)

nfinite(result)

The number of finite solutions.

nreal(result; tol=1e-6)

The number of real solutions where all imaginary parts of each solution are smaller than tol.

nreal(res::MonodromyResult; tol=1e-6)

Counts how many solutions of res have imaginary part norm less than 1e-6.

nsingular(result; singular_tol=1e10, multiplicitytol=1e-5, counting_multiplicities=false, kwargs...)

The number of singular solutions. A solution is considered singular if its windingnumber is larger than 1 or the condition number is larger than tol. If counting_multiplicities=true the number of singular solutions times their multiplicities is returned.

nnonsingular(result; tol=1e-10)

The number of non-singular solutions.

nat_infinity(result)

The number of solutions at infinity.

nafailed(result)

The number of failed paths.

bezout_number(F::MPPolys; variable_groups=[variables(F)], homvars=nothing, parameters=nothing)
bezout_number(multidegrees, groups::VariableGroups)

Compute the multi-homogeneous bezout number associated to the given system and variable groups.

mixed_volume(F::Vector{<:MP.AbstractPolynomialLike}; show_progress=true, algorithm=:regeneration)
mixed_volume(𝑨::Vector{<:Matrix}; show_progress=true, algorithm=:regeneration)

Compute the mixed volume of the given polynomial system F resp. represented by the support 𝑨. There are two possible values for algorithm:

• :total_degree: Use the total degree homotopy algorithm described in Section 7.1
• :regeneration: Use the tropical regeneration algorithm described in Section 7.2