Rosetta Code

Bernstein basis polynomials

Revision as of 20:01, 26 May 2023 by PureFox (talk | contribs) (Added Wren)

The Bernstein basis polynomials of degree are defined as

Bernstein basis polynomials is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Any real polynomial can written as a linear combination of such Bernstein basis polynomials. Let us call the coefficients in such linear combinations Bernstein coefficients.

The goal of this task is to write subprograms for working with degree-2 and degree-3 Bernstein coefficients. A programmer is likely to have to deal with such representations. For example, OpenType fonts store glyph outline data as as either degree-2 or degree-3 Bernstein coefficients.

The task is as follows:

  1. Write a subprogram (which we will call Subprogram (1)) to transform the ordinary monomial coefficients of a real polynomial, of degree 2 or less, to Bernstein coefficients for the degree-2 basis.
  2. Write Subprogram (2) to evaluate, at a given point, a real polynomial written in degree-2 Bernstein coefficients. For this use de Casteljau's algorithm.
  3. Write Subprogram (3) to transform the monomial coefficients of a real polynomial, of degree 3 or less, to Bernstein coefficients for degree 3.
  4. Write Subprogram (4) to evaluate, at a given point, a real polynomial written in degree-3 Bernstein coefficients. Use de Casteljau's algorithm.
  5. Write Subprogram (5) to transform Bernstein coefficients for degree 2 to Bernstein coefficients for degree 3.
  6. For the following polynomials, use Subprogram (1) to find Bernstein coefficients for the degree-2 basis: , . Display the results.
  7. Use Subprogram (2) to evaluate and at and at . Display the results. Optionally also evaluate the polynomials in the monomial basis and display the results. (For the monomial basis it is best to use a Horner scheme.)
  8. Use Subprogram (3) to find degree-3 Bernstein coefficients for and , and additionally also for . Display the results.
  9. Use Subprogram (4) to evaluate , , and at and . Display the results. Optionally also evaluate them in the monomial basis and display the results. (You may have done that already for and .)
  10. Use Subprogram (5) to get degree-3 Bernstein coefficients for and from their degree-2 Bernstein coefficients. Display the results.

ALGOL 60 and Python implementations are provided as initial examples. The latter does the optional monomial-basis evaluations.

Here is ALGOL 60 more nicely printed: 

procedure tobern2 (a0, a1, a2, b0, b1, b2);
 value a0, a1, a2; comment pass by value;
 real a0, a1, a2;
 real b0, b1, b2;
comment Subprogram (1): transform monomial coefficients a0, a1, a2 to Bernstein coefficients b0, b1, b2;
begin
 b0 := a0;
 b1 := a0 + ((1/2) × a1);
 b2 := a0 + a1 + a2
end tobern2;

real procedure evalbern2 (b0, b1, b2, t);
 value b0, b1, b2, t;
 real b0, b1, b2, t;
comment Subprogram (2): evaluate, at t, the polynomial with Bernstein coefficients b0, b1, b2. Use de Casteljau’s algorithm;
begin
 real s, b01, b12, b012;
 s := 1 - t;
 b01 := (s × b0) + (t × b1);
 b12 := (s × b1) + (t × b2);
 b012 := (s × b01) + (t × b12);
 evalbern2 := b012
end evalbern2;

procedure tobern3 (a0, a1, a2, a3, b0, b1, b2, b3);
 value a0, a1, a2, a3;
 real a0, a1, a2, a3;
 real b0, b1, b2, b3;
comment Subprogram (3): transform monomial coefficients a0, a1, a2, a3 to Bernstein coefficients b0, b1, b2, b3;
begin
 b0 := a0;
 b1 := a0 + ((1/3) × a1);
 b2 := a0 + ((2/3) × a1) + ((1/3) × a2);
 b3 := a0 + a1 + a2 + a3
end tobern3;  

real procedure evalbern3 (b0, b1, b2, b3, t);
 value b0, b1, b2, b3, t;
 real b0, b1, b2, b3, t;
comment Subprogram (4): evaluate, at t, the polynomial with Bernstein coefficients b0, b1, b2, b3. Use de Casteljau's algorithm;
begin
 real s, b01, b12, b23, b012, b123, b0123;
 s := 1 - t;
 b01 := (s × b0) + (t × b1);
 b12 := (s × b1) + (t × b2);
 b23 := (s × b2) + (t × b3);
 b012 := (s × b01) + (t × b12);
 b123 := (s × b12) + (t × b23);
 b0123 := (s × b012) + (t × b123);
 evalbern3 := b0123
end evalbern3;

procedure bern2to3 (q0, q1, q2, c0, c1, c2, c3);
 value q0, q1, q2;
 real q0, q1, q2;
 real c0, c1, c2, c3;
comment Subprogram (5): transform the quadratic Bernstein coefficients q0, q1, q2 to the cubic Bernstein coefficients c0, c1, c2, c3;
begin
 c0 := q0;
 c1 := ((1/3) × q0) + ((2/3) × q1);
 c2 := ((2/3) × q1) + ((1/3) × q2);
 c3 := q2
end bern2to3;

ALGOL 60

Works with: GNU MARST version 2.7
procedure tobern2 (a0, a1, a2, b0, b1, b2);
  value a0, a1, a2; comment pass by value;
  real a0, a1, a2;
  real b0, b1, b2;
comment Subprogram (1): transform monomial coefficients
        a0, a1, a2 to Bernstein coefficients b0, b1, b2;
begin
  b0 := a0;
  b1 := a0 + ((1/2) * a1);
  b2 := a0 + a1 + a2
end tobern2;  

real procedure evalbern2 (b0, b1, b2, t);
  value b0, b1, b2, t;
  real b0, b1, b2, t;
comment Subprogram (2): evaluate, at t, the polynomial with
        Bernstein coefficients b0, b1, b2. Use de Casteljau's
        algorithm;
begin
  real s, b01, b12, b012;
  s := 1 - t;
  b01 := (s * b0) + (t * b1);
  b12 := (s * b1) + (t * b2);
  b012 := (s * b01) + (t * b12);
  evalbern2 := b012
end evalbern2;

procedure tobern3 (a0, a1, a2, a3, b0, b1, b2, b3);
  value a0, a1, a2, a3;
  real a0, a1, a2, a3;
  real b0, b1, b2, b3;
comment Subprogram (3): transform monomial coefficients
        a0, a1, a2, a3 to Bernstein coefficients b0, b1,
        b2, b3;
begin
  b0 := a0;
  b1 := a0 + ((1/3) * a1);
  b2 := a0 + ((2/3) * a1) + ((1/3) * a2);
  b3 := a0 + a1 + a2 + a3
end tobern3;  

real procedure evalbern3 (b0, b1, b2, b3, t);
  value b0, b1, b2, b3, t;
  real b0, b1, b2, b3, t;
comment Subprogram (4): evaluate, at t, the polynomial
        with Bernstein coefficients b0, b1, b2, b3. Use
        de Casteljau's algorithm;
begin
  real s, b01, b12, b23, b012, b123, b0123;
  s := 1 - t;
  b01 := (s * b0) + (t * b1);
  b12 := (s * b1) + (t * b2);
  b23 := (s * b2) + (t * b3);
  b012 := (s * b01) + (t * b12);
  b123 := (s * b12) + (t * b23);
  b0123 := (s * b012) + (t * b123);
  evalbern3 := b0123
end evalbern3;

procedure bern2to3 (q0, q1, q2, c0, c1, c2, c3);
  value q0, q1, q2;
  real q0, q1, q2;
  real c0, c1, c2, c3;
comment Subprogram (5): transform the quadratic Bernstein
        coefficients q0, q1, q2 to the cubic Bernstein
        coefficients c0, c1, c2, c3;
begin
  c0 := q0;
  c1 := ((1/3) * q0) + ((2/3) * q1);
  c2 := ((2/3) * q1) + ((1/3) * q2);
  c3 := q2
end bern2to3;

begin
  real p0m, p1m, p2m;
  real q0m, q1m, q2m;
  real r0m, r1m, r2m, r3m;

  real p0b2, p1b2, p2b2;
  real q0b2, q1b2, q2b2;

  real p0b3, p1b3, p2b3, p3b3;
  real q0b3, q1b3, q2b3, q3b3;
  real r0b3, r1b3, r2b3, r3b3;

  real pc0, pc1, pc2, pc3;
  real qc0, qc1, qc2, qc3;

  real x, y;

  p0m := 1; p1m := 0; p2m := 0;
  q0m := 1; q1m := 2; q2m := 3;
  r0m := 1; r1m := 2; r2m := 3; r3m := 4;

  tobern2 (p0m, p1m, p2m, p0b2, p1b2, p2b2);
  tobern2 (q0m, q1m, q2m, q0b2, q1b2, q2b2);
  outstring (1, "Subprogram (1) examples:\n");
  outstring (1, "  mono ( ");
  outreal (1, p0m); outstring (1, ", ");
  outreal (1, p1m); outstring (1, ", ");
  outreal (1, p2m); outstring (1, ") --> bern ( ");
  outreal (1, p0b2); outstring (1, ", ");
  outreal (1, p1b2); outstring (1, ", ");
  outreal (1, p2b2); outstring (1, ")\n");
  outstring (1, "  mono ( ");
  outreal (1, q0m); outstring (1, ", ");
  outreal (1, q1m); outstring (1, ", ");
  outreal (1, q2m); outstring (1, ") --> bern ( ");
  outreal (1, q0b2); outstring (1, ", ");
  outreal (1, q1b2); outstring (1, ", ");
  outreal (1, q2b2); outstring (1, ")\n");

  outstring (1, "Subprogram (2) examples:\n");
  x := 0.25;
  y := evalbern2 (p0b2, p1b2, p2b2, x);
  outstring (1, "  p ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 7.50;
  y := evalbern2 (p0b2, p1b2, p2b2, x);
  outstring (1, "  p ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 0.25;
  y := evalbern2 (q0b2, q1b2, q2b2, x);
  outstring (1, "  q ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 7.50;
  y := evalbern2 (q0b2, q1b2, q2b2, x);
  outstring (1, "  q ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");

  tobern3 (p0m, p1m, p2m, 0, p0b3, p1b3, p2b3, p3b3);
  tobern3 (q0m, q1m, q2m, 0, q0b3, q1b3, q2b3, q3b3);
  tobern3 (r0m, r1m, r2m, r3m, r0b3, r1b3, r2b3, r3b3);
  outstring (1, "Subprogram (3) examples:\n");
  outstring (1, "  mono ( ");
  outreal (1, p0m); outstring (1, ", ");
  outreal (1, p1m); outstring (1, ", ");
  outreal (1, p2m); outstring (1, ", ");
  outreal (1, 0); outstring (1, ") --> bern ( ");
  outreal (1, p0b3); outstring (1, ", ");
  outreal (1, p1b3); outstring (1, ", ");
  outreal (1, p2b3); outstring (1, ", ");
  outreal (1, p3b3); outstring (1, ")\n");
  outstring (1, "  mono ( ");
  outreal (1, q0m); outstring (1, ", ");
  outreal (1, q1m); outstring (1, ", ");
  outreal (1, q2m); outstring (1, ", ");
  outreal (1, 0); outstring (1, ") --> bern ( ");
  outreal (1, q0b3); outstring (1, ", ");
  outreal (1, q1b3); outstring (1, ", ");
  outreal (1, q2b3); outstring (1, ", ");
  outreal (1, q3b3); outstring (1, ")\n");
  outstring (1, "  mono ( ");
  outreal (1, r0m); outstring (1, ", ");
  outreal (1, r1m); outstring (1, ", ");
  outreal (1, r2m); outstring (1, ", ");
  outreal (1, r3m); outstring (1, ") --> bern ( ");
  outreal (1, r0b3); outstring (1, ", ");
  outreal (1, r1b3); outstring (1, ", ");
  outreal (1, r2b3); outstring (1, ", ");
  outreal (1, r3b3); outstring (1, ")\n");

  outstring (1, "Subprogram (4) examples:\n");
  x := 0.25;
  y := evalbern3 (p0b3, p1b3, p2b3, p3b3, x);
  outstring (1, "  p ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 7.50;
  y := evalbern3 (p0b3, p1b3, p2b3, p3b3, x);
  outstring (1, "  p ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 0.25;
  y := evalbern3 (q0b3, q1b3, q2b3, q3b3, x);
  outstring (1, "  q ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 7.50;
  y := evalbern3 (q0b3, q1b3, q2b3, q3b3, x);
  outstring (1, "  q ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 0.25;
  y := evalbern3 (r0b3, r1b3, r2b3, r3b3, x);
  outstring (1, "  r ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");
  x := 7.50;
  y := evalbern3 (r0b3, r1b3, r2b3, r3b3, x);
  outstring (1, "  r ( "); outreal (1, x);
  outstring (1, ") = "); outreal (1, y);
  outstring (1, "\n");

  bern2to3 (p0b2, p1b2, p2b2, pc0, pc1, pc2, pc3);
  bern2to3 (q0b2, q1b2, q2b2, qc0, qc1, qc2, qc3);
  outstring (1, "Subprogram (5) examples:\n");
  outstring (1, "  bern ( ");
  outreal (1, p0b2); outstring (1, ", ");
  outreal (1, p1b2); outstring (1, ", ");
  outreal (1, p2b2); outstring (1, ") --> bern ( ");
  outreal (1, pc0); outstring (1, ", ");
  outreal (1, pc1); outstring (1, ", ");
  outreal (1, pc2); outstring (1, ", ");
  outreal (1, pc3); outstring (1, ")\n");
  outstring (1, "  bern ( ");
  outreal (1, q0b2); outstring (1, ", ");
  outreal (1, q1b2); outstring (1, ", ");
  outreal (1, q2b2); outstring (1, ") --> bern ( ");
  outreal (1, qc0); outstring (1, ", ");
  outreal (1, qc1); outstring (1, ", ");
  outreal (1, qc2); outstring (1, ", ");
  outreal (1, qc3); outstring (1, ")\n")
end
Output:
Subprogram (1) examples:
  mono ( 1 , 0 , 0 ) --> bern ( 1 , 1 , 1 )
  mono ( 1 , 2 , 3 ) --> bern ( 1 , 2 , 6 )
Subprogram (2) examples:
  p ( 0.25 ) = 1 
  p ( 7.5 ) = 1 
  q ( 0.25 ) = 1.6875 
  q ( 7.5 ) = 184.75 
Subprogram (3) examples:
  mono ( 1 , 0 , 0 , 0 ) --> bern ( 1 , 1 , 1 , 1 )
  mono ( 1 , 2 , 3 , 0 ) --> bern ( 1 , 1.66666666667 , 3.33333333333 , 6 )
  mono ( 1 , 2 , 3 , 4 ) --> bern ( 1 , 1.66666666667 , 3.33333333333 , 10 )
Subprogram (4) examples:
  p ( 0.25 ) = 1 
  p ( 7.5 ) = 1 
  q ( 0.25 ) = 1.6875 
  q ( 7.5 ) = 184.75 
  r ( 0.25 ) = 1.75 
  r ( 7.5 ) = 1872.25 
Subprogram (5) examples:
  bern ( 1 , 1 , 1 ) --> bern ( 1 , 1 , 1 , 1 )
  bern ( 1 , 2 , 6 ) --> bern ( 1 , 1.66666666667 , 3.33333333333 , 6 )

Julia

Translation of: Python
""" rosettacode.org/wiki/Bernstein_basis_polynomials """

""" Subprogram (1). """
monomial_to_bernstein_degree2(mc) = ((a0, a1, a2) = mc; (a0, a0 + ((1/2) * a1), a0 + a1 + a2))

""" Subprogram (2). """
function evaluate_bernstein_degree2(bernstein_coefficients, t)
    (b0, b1, b2) = bernstein_coefficients
    # de Casteljau’s algorithm.
    s = 1 - t
    b01 = (s * b0) + (t * b1)
    b12 = (s * b1) + (t * b2)
    b012 = (s * b01) + (t * b12)
    return b012
end

""" Subprogram (3). """
function monomial_to_bernstein_degree3(monomial_coefficients)
    (a0, a1, a2, a3) = monomial_coefficients
    return (a0, a0 + ((1/3) * a1),
            a0 + ((2/3) * a1) + ((1/3) * a2),
            a0 + a1 + a2 + a3)
end

""" Subprogram (4). """
function evaluate_bernstein_degree3(bernstein_coefficients, t)
    (b0, b1, b2, b3) = bernstein_coefficients
    # de Casteljau’s algorithm.
    s = 1 - t
    b01 = (s * b0) + (t * b1)
    b12 = (s * b1) + (t * b2)
    b23 = (s * b2) + (t * b3)
    b012 = (s * b01) + (t * b12)
    b123 = (s * b12) + (t * b23)
    b0123 = (s * b012) + (t * b123)
    return b0123
end

""" Subprogram (5). """
function bernstein_degree2_to_degree3(bernstein_coefficients)
    (b0, b1, b2) = bernstein_coefficients
    return (b0, ((1/3) * b0) + ((2/3) * b1),
            ((2/3) * b1) + ((1/3) * b2), b2)
end

evaluate_monomial_degree2(mc, t) = ((a0, a1, a2) = mc; a0 + (t * (a1 + (t * a2))))

evaluate_monomial_degree3(mc, t) = ((a0, a1, a2, a3) = mc; a0 + (t * (a1 + (t * (a2 + (t * a3))))))

#
# For the following polynomials, use Subprogram (1) to find
# coefficients in the degree-2 Bernstein basis:
#
#   p(x) = 1
#   q(x) = 1 + 2x + 3x²
#
# Display the results.
#
pmono2 = (1, 0, 0)
qmono2 = (1, 2, 3)
pbern2 = monomial_to_bernstein_degree2(pmono2)
qbern2 = monomial_to_bernstein_degree2(qmono2)
println("Subprogram (1) examples:")
println("  mono ", pmono2, " --> bern ", pbern2)
println("  mono ", qmono2, " --> bern ", qbern2)

#
# Use Subprogram (2) to evaluate p(x) and q(x) at x = 0.25, 7.50.
# Display the results. Optionally also display results from evaluating
# in the original monomial basis.
println("Subprogram (2) examples:")
for x in (0.25, 7.50)
    println("  p(", x, ") = ", evaluate_bernstein_degree2(pbern2, x),
          " ( mono: ", evaluate_monomial_degree2(pmono2, x), ")")
end
for x in (0.25, 7.50)
    println("  q(", x, ") = ", evaluate_bernstein_degree2(qbern2, x),
          " ( mono: ", evaluate_monomial_degree2(qmono2, x), ")")
end

#
# For the following polynomials, use Subprogram (3) to find
# coefficients in the degree-3 Bernstein basis:
#
#   p(x) = 1
#   q(x) = 1 + 2x + 3x²
#   r(x) = 1 + 2x + 3x² + 4x³
#
# Display the results.
#
pmono3 = (1, 0, 0, 0)
qmono3 = (1, 2, 3, 0)
rmono3 = (1, 2, 3, 4)
pbern3 = monomial_to_bernstein_degree3(pmono3)
qbern3 = monomial_to_bernstein_degree3(qmono3)
rbern3 = monomial_to_bernstein_degree3(rmono3)
println("Subprogram (3) examples:")
println("  mono ", pmono3, " --> bern ", pbern3)
println("  mono ", qmono3, " --> bern ", qbern3)
println("  mono ", rmono3, " --> bern ", rbern3)

#
# Use Subprogram (4) to evaluate p(x), q(x), and r(x) at x = 0.25,
# 7.50.  Display the results. Optionally also display results from
# evaluating in the original monomial basis.
println("Subprogram (4) examples:")
for x in (0.25, 7.50)
    println("  p(", x, ") = " , evaluate_bernstein_degree3(pbern3, x),
          " ( mono: " , evaluate_monomial_degree3(pmono3, x), ")")
end
for x in (0.25, 7.50)
    println("  q(", x, ") = ", evaluate_bernstein_degree3(qbern3, x),
          " ( mono: ", evaluate_monomial_degree3(qmono3, x), ")")
end
for x in (0.25, 7.50)
    println("  r(", x, ") = ", evaluate_bernstein_degree3(rbern3, x),
          " ( mono: ", evaluate_monomial_degree3(rmono3, x), ")")
end
#
# For the following polynomials, using the result of Subprogram (1)
# applied to the polynomial, use Subprogram (5) to get coefficients
# for the degree-3 Bernstein basis:
#
#   p(x) = 1
#   q(x) = 1 + 2x + 3x²
#
# Display the results.
#
println("Subprogram (5) examples:")
pbern3a = bernstein_degree2_to_degree3(pbern2)
qbern3a = bernstein_degree2_to_degree3(qbern2)
println("  bern ", pbern2, " --> bern ", pbern3a)
println("  bern ", qbern2, " --> bern ", qbern3a)
Output:

Same as Python example.

Python

#!/bin/env python3

# Subprogram (1).
def monomial_to_bernstein_degree2(monomial_coefficients):
    (a0, a1, a2) = monomial_coefficients
    return (a0, a0 + ((1/2) * a1), a0 + a1 + a2)

# Subprogram (2).
def evaluate_bernstein_degree2(bernstein_coefficients, t):
    (b0, b1, b2) = bernstein_coefficients
    # de Casteljau’s algorithm.
    s = 1 - t
    b01 = (s * b0) + (t * b1)
    b12 = (s * b1) + (t * b2)
    b012 = (s * b01) + (t * b12)
    return b012

# Subprogram (3).
def monomial_to_bernstein_degree3(monomial_coefficients):
    (a0, a1, a2, a3) = monomial_coefficients
    return (a0, a0 + ((1/3) * a1),
            a0 + ((2/3) * a1) + ((1/3) * a2),
            a0 + a1 + a2 + a3)

# Subprogram (4).
def evaluate_bernstein_degree3(bernstein_coefficients, t):
    (b0, b1, b2, b3) = bernstein_coefficients
    # de Casteljau’s algorithm.
    s = 1 - t
    b01 = (s * b0) + (t * b1)
    b12 = (s * b1) + (t * b2)
    b23 = (s * b2) + (t * b3)
    b012 = (s * b01) + (t * b12)
    b123 = (s * b12) + (t * b23)
    b0123 = (s * b012) + (t * b123)
    return b0123

# Subprogram (5).
def bernstein_degree2_to_degree3(bernstein_coefficients):
    (b0, b1, b2) = bernstein_coefficients
    return (b0, ((1/3) * b0) + ((2/3) * b1),
            ((2/3) * b1) + ((1/3) * b2), b2)

def evaluate_monomial_degree2(monomial_coefficients, t):
    (a0, a1, a2) = monomial_coefficients
    # Horner’s rule.
    return a0 + (t * (a1 + (t * a2)))

def evaluate_monomial_degree3(monomial_coefficients, t):
    (a0, a1, a2, a3) = monomial_coefficients
    # Horner’s rule.
    return a0 + (t * (a1 + (t * (a2 + (t * a3)))))

#
# For the following polynomials, use Subprogram (1) to find
# coefficients in the degree-2 Bernstein basis:
#
#   p(x) = 1
#   q(x) = 1 + 2x + 3x²
#
# Display the results.
#
pmono2 = (1, 0, 0)
qmono2 = (1, 2, 3)
pbern2 = monomial_to_bernstein_degree2(pmono2)
qbern2 = monomial_to_bernstein_degree2(qmono2)
print("Subprogram (1) examples:")
print("  mono", pmono2, " --> bern", pbern2)
print("  mono", qmono2, " --> bern", qbern2)

#
# Use Subprogram (2) to evaluate p(x) and q(x) at x = 0.25, 7.50.
# Display the results. Optionally also display results from evaluating
# in the original monomial basis.
print("Subprogram (2) examples:")
for x in (0.25, 7.50):
    print("  p(", x, ") = ", evaluate_bernstein_degree2(pbern2, x),
          " ( mono: ", evaluate_monomial_degree2(pmono2, x), ")")
for x in (0.25, 7.50):
    print("  q(", x, ") = ", evaluate_bernstein_degree2(qbern2, x),
          " ( mono: ", evaluate_monomial_degree2(qmono2, x), ")")

#
# For the following polynomials, use Subprogram (3) to find
# coefficients in the degree-3 Bernstein basis:
#
#   p(x) = 1
#   q(x) = 1 + 2x + 3x²
#   r(x) = 1 + 2x + 3x² + 4x³
#
# Display the results.
#
pmono3 = (1, 0, 0, 0)
qmono3 = (1, 2, 3, 0)
rmono3 = (1, 2, 3, 4)
pbern3 = monomial_to_bernstein_degree3(pmono3)
qbern3 = monomial_to_bernstein_degree3(qmono3)
rbern3 = monomial_to_bernstein_degree3(rmono3)
print("Subprogram (3) examples:")
print("  mono", pmono3, " --> bern", pbern3)
print("  mono", qmono3, " --> bern", qbern3)
print("  mono", rmono3, " --> bern", rbern3)

#
# Use Subprogram (4) to evaluate p(x), q(x), and r(x) at x = 0.25,
# 7.50.  Display the results. Optionally also display results from
# evaluating in the original monomial basis.
print("Subprogram (4) examples:")
for x in (0.25, 7.50):
    print("  p(", x, ") = ", evaluate_bernstein_degree3(pbern3, x),
          " ( mono: ", evaluate_monomial_degree3(pmono3, x), ")")
for x in (0.25, 7.50):
    print("  q(", x, ") = ", evaluate_bernstein_degree3(qbern3, x),
          " ( mono: ", evaluate_monomial_degree3(qmono3, x), ")")
for x in (0.25, 7.50):
    print("  r(", x, ") = ", evaluate_bernstein_degree3(rbern3, x),
          " ( mono: ", evaluate_monomial_degree3(rmono3, x), ")")

#
# For the following polynomials, using the result of Subprogram (1)
# applied to the polynomial, use Subprogram (5) to get coefficients
# for the degree-3 Bernstein basis:
#
#   p(x) = 1
#   q(x) = 1 + 2x + 3x²
#
# Display the results.
#
print("Subprogram (5) examples:")
pbern3a = bernstein_degree2_to_degree3(pbern2)
qbern3a = bernstein_degree2_to_degree3(qbern2)
print("  bern", pbern2, " --> bern", pbern3a)
print("  bern", qbern2, " --> bern", qbern3a)
Output:
Subprogram (1) examples:
  mono (1, 0, 0)  --> bern (1, 1.0, 1)
  mono (1, 2, 3)  --> bern (1, 2.0, 6)
Subprogram (2) examples:
  p( 0.25 ) =  1.0  ( mono:  1.0 )
  p( 7.5 ) =  1.0  ( mono:  1.0 )
  q( 0.25 ) =  1.6875  ( mono:  1.6875 )
  q( 7.5 ) =  184.75  ( mono:  184.75 )
Subprogram (3) examples:
  mono (1, 0, 0, 0)  --> bern (1, 1.0, 1.0, 1)
  mono (1, 2, 3, 0)  --> bern (1, 1.6666666666666665, 3.333333333333333, 6)
  mono (1, 2, 3, 4)  --> bern (1, 1.6666666666666665, 3.333333333333333, 10)
Subprogram (4) examples:
  p( 0.25 ) =  1.0  ( mono:  1.0 )
  p( 7.5 ) =  1.0  ( mono:  1.0 )
  q( 0.25 ) =  1.6874999999999998  ( mono:  1.6875 )
  q( 7.5 ) =  184.75000000000034  ( mono:  184.75 )
  r( 0.25 ) =  1.7499999999999998  ( mono:  1.75 )
  r( 7.5 ) =  1872.25  ( mono:  1872.25 )
Subprogram (5) examples:
  bern (1, 1.0, 1)  --> bern (1, 1.0, 1.0, 1)
  bern (1, 2.0, 6)  --> bern (1, 1.6666666666666665, 3.333333333333333, 6)

Wren

Translation of: Algol 60
Library: Wren-fmt
import "./fmt" for Fmt

var toBern2 = Fn.new { |a| [a[0], a[0] + a[1] / 2, a[0] + a[1] + a[2]] }

// uses de Casteljau's algorithm
var evalBern2 = Fn.new { |b, t|
    var s = 1 - t
    var b01 = s * b[0] + t * b[1]
    var b12 = s * b[1] + t * b[2]
    return s * b01 + t * b12
}

var toBern3 = Fn.new { |a|
    var b = List.filled(4, 0)
    b[0] = a[0]
    b[1] = a[0] + a[1] / 3
    b[2] = a[0] + a[1] * 2/3 + a[2] / 3
    b[3] = a[0] + a[1] + a[2] + a[3]
    return b
}

// uses de Casteljau's algorithm
var evalBern3 = Fn.new { |b, t|
    var s = 1 - t
    var b01  = s * b[0] + t * b[1]
    var b12  = s * b[1] + t * b[2]
    var b23  = s * b[2] + t * b[3]
    var b012 = s * b01  + t * b12
    var b123 = s * b12  + t * b23
    return s * b012 + t * b123
}

var bern2to3 = Fn.new { |q|
    var c = List.filled(4, 0)
    c[0] = q[0]
    c[1] = q[0] / 3   + q[1] * 2/3
    c[2] = q[1] * 2/3 + q[2] / 3
    c[3] = q[2]
    return c
}

var pm = [1, 0, 0]
var qm = [1, 2, 3]
var rm = [1, 2, 3, 4]
var x
var y

System.print("Subprogram(1) examples:")
var pb2 = toBern2.call(pm)
var qb2 = toBern2.call(qm)
Fmt.print("mono$n --> bern$n", pm, pb2)
Fmt.print("mono$n --> bern$n", qm, qb2)

System.print("\nSubprogram(2) examples:")
x = 0.25
y = evalBern2.call(pb2, x)
Fmt.print("p($4.2f) = $j", x, y)
x = 7.5
y = evalBern2.call(pb2, x)
Fmt.print("p($4.2f) = $j", x, y)
x = 0.25
y = evalBern2.call(qb2, x)
Fmt.print("q($4.2f) = $j", x, y)
x = 7.5
y = evalBern2.call(qb2, x)
Fmt.print("q($4.2f) = $j", x, y)

System.print("\nSubprogram(3) examples:")
pm.add(0)
qm.add(0)
var pb3 = toBern3.call(pm)
var qb3 = toBern3.call(qm)
var rb3 = toBern3.call(rm)
Fmt.print("mono$n --> bern$n", pm, pb3)
Fmt.print("mono$n --> bern$n", qm, qb3)
Fmt.print("mono$n --> bern$n", rm, rb3)

System.print("\nSubprogram(4) examples:")
x = 0.25
y = evalBern3.call(pb3, x)
Fmt.print("p($4.2f) = $j", x, y)
x = 7.5
y = evalBern3.call(pb3, x)
Fmt.print("p($4.2f) = $j", x, y)
x = 0.25
y = evalBern3.call(qb3, x)
Fmt.print("q($4.2f) = $j", x, y)
x = 7.5
y = evalBern3.call(qb3, x)
Fmt.print("q($4.2f) = $j", x, y)
x = 0.25
y = evalBern3.call(rb3, x)
Fmt.print("r($4.2f) = $j", x, y)
x = 7.5
y = evalBern3.call(rb3, x)
Fmt.print("r($4.2f) = $j", x, y)

System.print("\nSubprogram(5) examples:")
var pc = bern2to3.call(pb2)
var qc = bern2to3.call(qb2)
Fmt.print("mono$n --> bern$n", pb2, pc)
Fmt.print("mono$n --> bern$n", qb2, qc)
Output:
Subprogram(1) examples:
mono[1, 0, 0] --> bern[1, 1, 1]
mono[1, 2, 3] --> bern[1, 2, 6]

Subprogram(2) examples:
p(0.25) = 1
p(7.50) = 1
q(0.25) = 1.6875
q(7.50) = 184.75

Subprogram(3) examples:
mono[1, 0, 0, 0] --> bern[1, 1, 1, 1]
mono[1, 2, 3, 0] --> bern[1, 1.6666666666667, 3.3333333333333, 6]
mono[1, 2, 3, 4] --> bern[1, 1.6666666666667, 3.3333333333333, 10]

Subprogram(4) examples:
p(0.25) = 1
p(7.50) = 1
q(0.25) = 1.6875
q(7.50) = 184.75
r(0.25) = 1.75
r(7.50) = 1872.25

Subprogram(5) examples:
mono[1, 1, 1] --> bern[1, 1, 1, 1]
mono[1, 2, 6] --> bern[1, 1.6666666666667, 3.3333333333333, 6]
Retrieved from "https://rosettacode.org/wiki/Bernstein_basis_polynomials?oldid=343162"