On Computing Ranges of Polynomials Some Improvements to the Bernstein Approach

1
On Computing Ranges of Polynomials - Some
Improvements to the Bernstein Approach

• P. S. V. Nataraj
• Shashwati Ray
• IIT Bombay
• Presented at the Second Taylor Model Workshop,
  Miami, Dec. 2003
• Indian Institute of Technology, Bombay 2003

2
SCOPE OF THE PRESENTATION

• Introduction
• Bernstein forms
• Definition
• Properties
• Basis conversion
• Range calculations
• Important theorems
• Subdivision
• Subdivision direction selection strategies
• New propositions
• Comparison with Globsol
• Examples and conclusion
• Future work

3
INTRODUCTION

• Polynomials are encountered in Engineering and
  Scientific Applications.
• Required to find bounds for the range of
  polynomials over an interval.
• Our method based on expansion of multivariate
  polynomial into Bernstein form.

4
ADVANTAGES OF BERNSTEIN FORM

• Avoids function evaluations may be otherwise
  costly for polynomials of high degrees.
• Nice features about range enclosures with
  Bernstein form
• Information when enclosure is exact.
• Sharpness of enclosure can be improved by
  elevating the Bernstein degree or by subdivision
  of the domain

5
BERNSTEIN FORM

• Any polynomial in power form is represented as
• Bernstein form of representation is
• Bernstein coefficients are given by

6
BERNSTEIN FORM

• Bernstein functions on 0,1 are defined by
• On a general interval the
  Bernstein functions are defined as
• Each set of coefficients ai or bin can be
  computed from the other (i.e., power form to
  Bernstein form).

7
BASIS CONVERSION

• On unit interval a polynomials equivalent power
  and Bernstein forms
• Using matrix multiplication p(x) XA BXB
• X is the variable (row) matrix.
• A is the coefficient (column) matrix
• BX is the Bernstein basis (row) matrix
• B is the Bernstein coefficient (column) matrix
• After certain computations, BXB XUXB
• UX is a lower triangular matrix

8
BASIS CONVERSION

• BX XWXVXUX
• WX is an upper triangular matrix
• VX is a diagonal matrix
• XA BXB
• XA XWXVXUX B
• B(UX)-1 (VX)-1 (WX)-1A

9
BASIS CONVERSION

• For a bivariate case
• X1 and X2 are variable matrices
• By analogy with univariate case

10
RANGE CALCULATION IMPORTANT THEOREMS

• Theorem 1 The minimum and maximum Bernstein
  coefficients give an enclosure of the range of
  polynomial on the given interval.
• Theorem 2 Vertex condition If the minimum
  and maximum Bernstein coefficients of polynomial
  in Bernstein form occur at the vertices of the
  Bernstein coefficient array, then the enclosure
  is the exact range.
• Theorem 3 Bernstein approximations converge to
  the range and the convergence is at least linear
  in the order of approximations.

11
SUBDIVISION

• Minimum maximum Bernstein coefficients enclose
  the range of the polynomial
• Vertex condition not satisfied subdivide in the
  chosen component direction
• Calculate Bernstein coefficient of p(x) on all
  these subboxes
• Check vertex condition on each subbox
• Find range of polynomial

12
SWEEP PROCEDURE
13
SUBDIVISION DIRECTION SELECTION
STRATEGIESEXISTING RULES

• Rule A
• Direction chosen is cyclic
• Rule B (derivative based)
• Find an upper bound for the absolute value of
  the partial derivative (from its Bernstein form)
  in the rth direction
• Compute the same for each direction
• The direction for which this is maximum is
  selected as the optimal direction for bisection

14
SUBDIVISION DIRECTION SELECTION
STRATEGIESEXISTING RULES

• Rule C (derivative based)
• Find an upper bound for the absolute value of
  the partial derivative (from its Bernstein form)
  in the rth direction
• Find the product of this upper bound and the
  degree of rth variable
• Compute the same for each direction
• The direction for which this is maximum is
  selected as the optimal direction for bisection

15
SUBDIVISION DIRECTION SELECTION
STRATEGIESEXISTING RULES

• Rule D (derivative based)
• Find an upper bound for the absolute value of
  the partial derivative (from its Bernstein form)
  in the rth direction
• Find the product of this upper bound and the edge
  length
• Compute the same for each direction
• The direction for which this is maximum is
  selected as the optimal direction for bisection

16
SUBDIVISION DIRECTION SELECTION
STRATEGIESEXISTING RULES

• Rule E (max width)
• Bisection is done along the component direction
  of maximal width
• Rule F (randomized)
• Pick a box randomly from the current list
• Bisect it in all component directions
• Find the maximal width of the polynomial over
  both these boxes
• Choose the component direction where this width
  is least

17
DISADVANTAGES OF THE EXISTING METHODS OF
BERNSTEIN COEFFICIENTS COMPUTATIONS

• Conventional method is very inefficient for
  higher dimensional polynomials
• Even the Matrix multiplication method cannot be
  generalized to higher dimensional polynomials

18
DISADVANTAGES OF THE EXISTING METHODS OF
BISECTION DIRECTION SELECTION

• The existing rules generally need large
  computation times.
• Certain existing rules are grossly inefficient,
  especially for higher dimensional polynomials.

19
DISADVANTAGES OF THE EXISTING RANGE COMPUTATION
ALGORITHMS

• Fortran 95 can not readily cater to more than six
  dimensional arrays.
• Bernstein coefficients generated are stored in
  multidimensional arrays
• For getting sharp enclosures, subdivision usually
  creates large amount of data
• All these are found to heavily slow down
  computations

20
NEW PROPOSITIONS

• We propose new tools for
• Bernstein Coefficient Computations
• Subdivision Direction Selection
• Acceleration devices

21
PROPOSITION 1Computation of Bernstein
Coefficients

• Matrix Method
• This gives

22
PROPOSITION 1Computation of Bernstein
coefficients

• Extending the above to trivariate case
• where transpose means converting 3rd coordinate
  direction to 2nd , 2nd coordinate direction to
  1st , 1st coordinate direction to 3rd
• The same logic can be extended to the l-variate
  case

23
PROPOSED MATRIX METHOD

• For a 3-dimensional case instead of considering
  the polynomial coefficient matrix A as a 3
  dimensional array, it can be considered as a
  matrix with 0 to n1 rows and 0 to (n21)(n31)
  1 columns
• After the third transpose and reshape the A
  matrix would come in its original form
• The same is shown as figures in the report.

24
Proposed Algorithm Bernstein Matrix

• Inputs
• Degree of each variable
• Widths of each interval component
• Lower end point of each interval component
• Polynomial coefficients a in matrix form
• Output Bernstein coefficient a in matrix form
• BEGIN Algorithm
• Compute UX, VX and WX for each variable, j
  1,2,, n
• Compute their inverses (UX)-1 , (VX )-1 and
  (WX)-1 for j 1,2,, n
• Compute their product Mj for j 1,2,, n
• FOR j 1,2,, n Set a Mj a
• b transpose(a)
• a reshape(b)
• END FOR
• END Algorithm   

25
Comparison and results (Bernstein coefficients
evaluation)

• Comparison is done between conventional method
  and the proposed Matrix method
• All problems are taken from Verscheldes PHC pack
  (database of polynomial systems)
• Both methods generate the same Bernstein
  coefficients
• Proposed method gives results for all the
  examples
• Conventional method takes much longer time it
  does not give results for some examples even
  after a long time.

26
Time in seconds taken by the Matrix and the
Conventional methods to compute Bernstein
coefficients
27
Some Definitions

• Solution box is the box where vertex property
  of the Bernstein coefficients is satisfied.
• Solution patch is given by enclosure of the
  Bernstein coefficients over the solution box.
• Current range is the hull of all solution
  patches.

28
PROPOSITION 2Subdivision Direction Selection

• Proposed Rule 1 Select patch with minimum hull
• Find the patches with smallest and largest
  Bernstein coefficients from the entire list of
  patches
• Find the respective distance of smallest and
  largest Bernstein coefficient to the infimum and
  supremum of the current range
• Select that patch as the one which gives the
  largest distance, and bisect this selected patch
  in all component directions

29
EXAMPLEDistance

• Let a 3-d polynomial be given as
• The polynomial coefficient and the Bernstein
  coefficient matrices are given as
• Let at some stage the Current range be given by
  1,20

30
EXAMPLE

• Let the two existing patches in the entire list
  be
• Minimum and maximum value of Bernstein
  coefficients are 0.5 and 7.70 respectively
• The respective distances are (1-0.5) 0.5 and
  (20-7.70) 12.30
• Select the patch which has 7.70 as the Bernstein
  coefficient

31
Subdivision Direction Selection

• Select patch with minimum hull
• Hull of the Bernstein coefficients over both the
  generated subboxes is found out in each direction
  
• The selected bisection direction is the one which
  gives the minimal hull.

32
PROPOSITION 3Subdivision Direction Selection

• Proposed Rule 2 Randomized box with minimum
  hull
• Choose a box randomly from the current list.
• The chosen box is bisected along each and every
  component direction.
• Hull of the Bernstein coefficients over both the
  generated subboxes is found out in each
  direction.
• The selected bisection direction is the one which
  gives the minimal hull.

33
RESULTS AND CONCLUSIONS Subdivision Direction
Selection

• All problems are taken from Verscheldes PHC pack
  (database of polynomial systems)
• Comparison of the bisection direction selection
  rules, on the basis of the following performance
  metrics
• Number of solution boxes where vertex
  property/simplified vertex condition is satisfied
• Number of subdivisions
• Maximum patch list length
• Computation time in seconds
• With the proposed rules we are able to solve all
  the examples
• With rules A, B, C, D, E and F we are able to
  solve 90, 70, 70, 50, 90 and 90 of the
  examples, respectively

34
RESULTS AND CONCLUSIONS Subdivision Direction
Selection

• Comparison is done with all the existing rules
  (except rule D, which is unable to solve half the
  problems).
• Proposed rules are found to be the most
  efficient, in terms of every performance metric
• Reduction in number of solution boxes compared to
  Rule A is 82, Rule B is 83, Rule C is 83, Rule
  E is 81 and Rule F is 84
• Reduction in number of subdivisions compared to
  Rule A is 84, Rule B is 83, Rule C is 83, Rule
  E is 81 and Rule F is 84

35
RESULTS AND CONCLUSIONS Subdivision Direction
Selection

• Reduction in the maximum list length compared to
  Rule A is 85, Rule B is 51, Rule C is 46, Rule
  E is 34 and Rule F is 63
• Reduction in the computational time compared to
  Rule A is 98, Rule B is 93, Rule C is 93, Rule
  E is 91 and Rule F is 93
• Both the rules give the same amount of reductions

36
PROPOSITION 4Accelerating Devices for Range
Computation

• Cut off test
• Let the enclosure of Bernstein coefficients over
  a box be included in the current range.
• This box can be deleted from further processing.
• Any further subdivision of this box is not going
  to improve the current range.
• Therefore, unnecessary subdivisions are
  eliminated

37
PROPOSITION 5Accelerating Devices for Range
Computation

• Simplified Vertex Condition
• Let smallest Bernstein coefficient in a patch
  appear at the vertex and the upper bound of
  this patch be included in the current range
• This patch can be taken to be a solution patch
• No further subdivision is required, as it would
  not give any further improvement in the current
  range
• Same is true if the largest Bernstein coefficient
  in a patch is appearing at the vertex and the
  lower bound of the patch is included in the
  current range

38
PROPOSITION 6Accelerating Devices for Range
Computation

• Monotonicity Test
• If the polynomial is monotonic w.r.t any
  direction on a box
• If the box is in the interior
• This box can be rejected
• If the box is not in the interior
• This box is retained
• No subdivision of this box in that direction
• Avoids unnecessary subdivisions

39
PROPOSED ALGORITHMRange

• Inputs
• Maximum degree of each variable of the polynomial
• Coefficients of the polynomial in the form of a
  matrix
• Initial domain box
• Accuracy tolerances
• Outputs
• An enclosure of the range of the specified
  accuracy
• Number of solution boxes
• Number of subdivisions
• Maximum patch list length used to store the
  patches

40
PROPOSED ALGORITHMRange

• BEGIN Algorithm
• Evaluate Bernstein coefficients in matrix form
  put them as a patch in a list.
• Set flags for simplified vertex condition, cut
  off test and monotonicity test
• Take the first patch from the list and check for
  vertex/simplified vertex condition.
• If true, then update current range.
• Else, perform cut off test.
• If list is not empty, choose a component
  direction based on a bisection rule

41
PROPOSED ALGORITHMRange

• Take each patch from list, do the monotonicity
  test. If test is passed, proceed to next step
  else discard tested box and choose next patch
  from the list.
• Perform the necessary subdivision in the given
  direction and evaluate the new Bernstein
  coefficients
• Update the list of patches and remove the tested
  box
• Take the next patch from the list and perform the
  same operations, till the list is empty
• Get the range of the polynomial from the current
  range
• END Algorithm

42
COMPARISON WITH GLOBSOL

• We compare our results with the results obtained
  with Globsol
• All problems are taken from Verscheldes PHC
  pack, a database of polynomial systems.
• Comparison of the computation time required by
  both the algorithms in terms of values of
• Ratio
• Percent reduction

43
COMPARISON WITH GLOBSOL

44
Comparison of computation time by both the
algorithms
45
(No Transcript)
46
(No Transcript)
47
(No Transcript)
48
RESULTSComparison with Globsol

• In one of the examples Globsol was unable to give
  any result (arithmetic exception core dumped)
  indicated by in the table.
• Compared to Globsol, the improved Bernstein
  method gives an average percent reduction in
  computational time of 47.38

49
CONCLUSIONS

• The new algorithm Range is faster than Globsol in
  fourteen out of eighteen problems.
• Further, on the average, we obtain considerable
  reductions in computation time with the proposed
  Algorithm Range.

50
FUTURE WORK

• Modify the new algorithm after incorporating
  certain modifications in the Monotonicity test ,
  thus further speeding up the algorithm Range
• Integrate the code with COSY package of Berz et
  al.
• Apply the method to Engineering problems,
  especially control engineering !

51
THANK YOU !