Chapter 11. Length and Center of Mass of a Curve

Rovenski Vladimir, Haifa

> restart:

11.2 Calculation of Length and Center of Mass

Curves: graphs (a)-(c).

> restart: a:=-2: b:=2: f:=x -> x^2: # define your function

> p1:=plot(f, a..b, scaling=constrained): %;

> int(sqrt(1+diff(f(x), x)^2), x=a..b); L1:=evalf(%);

> My:=evalf(int(x*sqrt(1+diff(f(x), x)^2), x=a..b)); Mx:=evalf(int(f(x)*sqrt(1+diff(f(x), x)^2), x=a..b)); rc1:=[My/L1, Mx/L1];

> n:=6: A:=i->[a+i*(b-a)/n, f(a+i*(b-a)/n)]:

> p2:=plot([seq(A(i), i=0..n)], color=blue): %;

> plots[display]([p1, p2], scaling=constrained);

> L2:=evalf(sum(linalg[norm](A(i+1)-A(i), 2), i=0..n-1)); d:=abs(L1-L2)/L1;

> rc2:=evalm(sum(0.5*(A(i+1)+A(i))*linalg[norm](A(i+1)-A(i), 2),i=0..n-1)/L2); evalm(rc2-rc1);

Parametrized plane curves (d)-(f).

> a:=0: b:=2*Pi: X:=t->t-sin(t): Y:=t->1-cos(t): # cycloid

> p1:=plot([X, Y, a..b], scaling=constrained): %;

> int(sqrt(diff(X(t), t)^2+diff(Y(t), t)^2), t=a..b); L1:=evalf(%);

> Mx:=evalf(int(Y(t)*sqrt(diff(X(t),t)^2+diff(Y(t),t)^2), t=a..b)); My:=evalf(int(X(t)*sqrt(diff(X(t),t)^2+diff(Y(t),t)^2),t=a..b)); rc1:=[My/L1,Mx/L1];

> A:=i->[X(a+i*(b-a)/n), Y(a+i*(b-a)/n)]: n:=10:

> p2:=plot([seq(A(i), i=0..n)], scaling=constrained): %;

> plots[display]([p1, p2], scaling=constrained);

> L2:=evalf(sum(linalg[norm](A(i+1)-A(i), 2), i=0..n-1)); d:=abs(L1-L2)/L1;

> rc2:=evalm(sum(0.5*(A(i+1)+A(i))*linalg[norm](A(i+1)-A(i), 2),i=0..n-1)/L2); evalm(rc1-rc2);

Curves in polar coordinates: (g)-(h).

> a:=0: b:=2*Pi: rho:=t -> 0.2*t: n:=12: # Archimedes' spiral

> p1:= plots[polarplot](rho, a..b, scaling=constrained): %;

> L1:=evalf(int(sqrt(rho(t)^2+diff(rho(t),t)^2), t=a..b));

> T:=i->a+i*(b-a)/n: p2:=plot([seq([rho(T(i)), T(i)], i=0..n)], coords=polar):

> plots[display]([p1, p2], scaling=constrained);

> L2:=evalf(sum(sqrt(rho(T(i+1))^2+rho(T(i))^2-2*rho(T(i+1))*rho(T(i))*cos((b-a)/n)),i=0..n-1)); d:=(L2-L1)/L1;

Parametrized space curves (i)-(k).

> restart:

> a:=0: b:=2*Pi: X:=t->cos(t)^2: Y:=t->cos(t)*sin(t): Z:=t -> sin(t): #curve of Viviani

> p1:=plots[spacecurve]([X, Y, Z, a..b], axes=boxed): %;

> int(sqrt(simplify(diff(X(t),t)^2+diff(Y(t),t)^2+diff(Z(t),t)^2)), t=a..b); L1:=evalf(%);

> Myz:=evalf(int(X(t)*sqrt(simplify(diff(X(t),t)^2+diff(Y(t), t)^2+diff(Z(t), t)^2, trig)), t=a..b)); Mxz:=evalf(int(Y(t)*sqrt(simplify(diff(X(t), t)^2+diff(Y(t), t)^2+diff(Z(t), t)^2, trig)), t=a..b)); Mxy:=evalf(int(Z(t)*sqrt(simplify(diff(X(t), t)^2+diff(Y(t), t)^2+diff(Z(t), t)^2, trig)), t=a..b));

> rc1:=[Myz, Mxz, Mxy]/L1;

> T:=i->a+i*(b-a)/n; A:=i->[X(T(i)), Y(T(i)), Z(T(i))]: n:=8:

> p2:=plots[pointplot3d]([seq(A(i), i=0..n)], style=line):

> plots[display]([p1, p2], scaling=constrained);

> L2:=evalf(sum(linalg[norm](A(i+1)-A(i), 2), i=0..n-1)); d:=(L2-L1)/L1;

> rc2:=evalm(sum(0.5*(A(i+1)+A(i))*linalg[norm](A(i+1)-A(i), 2), i=0..n-1)/L2); evalm(rc1-rc2);

>