InputMode:=Word

PrecisionDigits:=64

"This is an example of calculating the length of a Bezier curve."

"This is a Bezier curve."

v1:=[x=a3*t^3+a2*t^2+a1*t+x0,y=b3*t^3+b2*t^2+b1*t+y0]

"The curve starts at location (x0,y0) when t=0."

SUBST(v1,t,0)

"The curve finishes at location (x3,y3) when t=1."

v2:=SUBST(v1,[x,y,t],[x3,y3,1])

"The slope and the rate of change of slope at t=0 is determined"

" by the control point (x1,y1)."

v3:=VECTOR(RHS(k),k,DIF(v1,t))

v4:=SUBST([3*(x1-x0),3*(y1-y0)]=v3,t,0)

"The slope and the rate of change of slope at t=1 is determined"

" by the control point (x2,y2)."

v5:=SUBST([3*(x3-x2),3*(y3-y2)]=v3,t,1)

"Solve for a1, a2, a3, b1, b2, and b3 in terms of "

"x0, x1, x2, x3, y0, y1, y2, and y3."

v6:=v2 SUB 1 AND v2 SUB 2 AND v4 AND v5

v7:=(SOLUTIONS(v6,[a1,a2,a3,b1,b2,b3])) SUB 1

"INT(int1,t,0,1) is the length of a Bezier curve. "

"It is an elliptic integral."

'DIF(s,t,1)='SQRT(DIF(x,t,1)^2+DIF(y,t,1)^2)

length='INT(SQRT(DIF(x,t,1)^2+DIF(y,t,1)^2),t,0,1)

int1:=SQRT(v3 . v3)

"Pick a actual Bezier curve to check the calculations."

"In this example Bezier curve [[x0,y0],[x1,y1],[x2,y2],[x3,y3]]=v8."

v8:=[[6,8],[1,10],[7,3],[4,4]]

"Solve for [a1,a2,a3,b1,b2,b3] in this example."

v10:=SUBST(v7,[x0,y0,x1,y1,x2,y2,x3,y3],APPEND(v8))

int2:=INT(SUBST(int1,[a1,a2,a3,b1,b2,b3],v10),t,0,1)

"refint is a check value calulated to 64 digits."

refint:=1809305842082148221456654802739005603491766996504376040840721629/2500~
00000000000000000000000000000000000000000000000000000000000

"Do numeric integration to show this check value is correct."

APPROX(int2-refint,16)

"Factor the polynomial to put the integral in a standard form."

v11:=SUBST(v3 . v3,[a1,a2,a3,b1,b2,b3],v10)

v12:=FACTORS(v11,Complex)

"km is the scale factor."

km:=SQRT(v12 SUB 5 SUB 1)

"mm=m in the elliptic integral I(m)"

mm:=[1,1,1,1]

am:=VECTOR(-SUBST(v12 SUB i SUB 1,t,0),i,1,4)

bm:=[-1,-1,-1,-1]

"Set i=1, n=4, and h=4."

LOAD("c2001.mth")

v13:=RE179(1,mm,4,4)

"This is the symbolic result for a general case where [a1, a2, a3, a4]"

"and [b1, b2, b3, b4] are not known."

v14:=RE175(v13,am,bm)

"AHAT(e1+e2+e3+e4), AHAT(e2+e3+e4), and AHAT(2*e1+e2+e3+e4)"

v15:=AHAT_xy([1,1,1,1],4,4,am,bm,1,0)

v16:=AHAT_xy([0,1,1,1],4,4,am,bm,1,0)

v17:=AHAT_xy([2,1,1,1],4,4,am,bm,1,0)

"Substitute the values of AHAT(e1+e2+e3+e4), AHAT(e2+e3+e4), "

"and AHAT(2*e1+e2+e3+e4)"

"into the expression v14."

v19:=SUBST(v14,[AHAT(e1+e2+e3+e4),AHAT(e2+e3+e4),AHAT(2*e1+e2+e3+e4)],[v15,v1~
6,v17])

"This is the results in terms of basic integrals I(-e1), I(e0), and I(e1)."

"Find the basic integrals in the table."

LOAD("c2001a.mth")

IS([IA(-e1),IA(0),IA(e1)],4)

"The DERIVE functions I428(a,b,x,y), I414(a,b,x,y), "

"and I429(a,b,x,y) are in c1998.mth."

"They are defined by equation 4.28, 4.14, and 4.29 in the paper."

"There is a revision to these equations in a latter paper"

"Reduction Theorems for Elliptic Integrands with Square Root "

"of two Quadratic Factors."

"Evaluate the basic integrals using the new notation."

"Now subsitute the values of the basic integrals into the expression v19 and"

"compare that to the check value."

v20:=APPROX(ISHAT(v19,am,bm,1,0),128)

v21:=APPROX(v20,128)

"Find the final integral."

v18:=APPROX(km*v21,128)

"Now compare v18 to the check value."

LOAD("carlson.mth")

APPROX(v18-refint,128)

"This result matches the check value to 62 places."