/* [wxMaxima batch file version 1] [ DO NOT EDIT BY HAND! ]*/
/* [ Created with wxMaxima version 0.8.4 ] */

/* [wxMaxima: title   start ]
Approximating a Small Circular Arc with a Cubic Bezier Curve
   [wxMaxima: title   end   ] */

/* [wxMaxima: comment start ]
Some graphical systems (e.g. PostScript) lack primitives for circles and circular arcs. Many systems do have
support for Bezier curves, however. We will show how to determine the control points of a cubic Bezier curve
that approximates a small circular arc. The Bezier curve may be positioned, scaled, and rotated arbitrarily.
This enables drawing a full circle with a cubic Bezier spline.
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
Definitions
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
We start by defining the equation for a general Bezier curve from Berstein polynomials.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
bernstein(i, n, t) := binomial(n, i)*t^i*(1-t)^(n-i);
bezier(n) := sum(P[i]*bernstein(i,n,t), i, 0, n);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
Specifically, we are interested in a cubic Bezier curve.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
C(t):=bezier(3);
C(t);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Geometrical Configuration
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
Consider a unit circular arc A of sweep < 90°, bisected by the x-axis. Let %phi=%theta/2.
The cubic Bezier curve is defined by four points, P[0..3]. P[0] and P[3] are endpoints, while
P[1] and P[2] are their respective control points.
   [wxMaxima: comment end   ] */

/* [wxMaxima: subsect start ]
Endpoints
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
Force the endpoints of C to the endpoints of A. Note the symmetry due to bisection.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
eq0 : x[0]=cos(%phi);
eq1 : y[0]=sin(%phi);
eq2 : x[3]=x[0];
eq3 : y[3]=-y[0];
/* [wxMaxima: input   end   ] */

/* [wxMaxima: subsect start ]
Midpoint
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
Force the midtime of the Bezier curve (at t=1/2) to pass through the midpoint of the arc at point (1, 0).
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
eq4 : at(subst(x,P,C(t)),t=1/2)=1;
eq5 : at(subst(y,P,C(t)),t=1/2)=0;
/* [wxMaxima: input   end   ] */

/* [wxMaxima: subsect start ]
Slope
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
The tangent slope at each endpoint is equal to the negative reciprocol of the the slope of the segment connecting
it to the arc's center. We ensure that this slope is equal to the the slope of the segment connecting each endpoint
to its control point.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
m0 : y[0]/x[0]$
eq6 : -1/m0 = (y[0]-y[1])/(x[0]-x[1]);
m3 : y[3]/x[3]$
eq7 : -1/m3 = (y[3]-y[2])/(x[3]-x[2]);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: subsect start ]
Solution
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
We now solve a system of eight variables (X[0..3] and Y[0..3]) and eight equations. The solution is simplified
through the application of a trignometric identity.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
ans : solve(
    [eq0, eq1, eq2, eq3, eq4, eq5, eq6, eq7],
    [x[0], x[1], x[2], x[3], y[0], y[1], y[2], y[3]]
)$
map(factor, map(trigsimp, ans[1]));
/* [wxMaxima: input   end   ] */

/* Maxima can't load/batch files which end with a comment! */
"Created with wxMaxima"$