> restart;
# The following examples is Jan Verschelde's examples:
# The purpose of this worksheet is to generate instances of multivariate
# polynomials with complex coefficients. Executing this worksheet
# generates the files exf1, exf2, exf3, and exf4.
#
# The Jacobian for a Stewart-Gough Platform
# We derive the Jacobian for a platform, first symbolically.
# 1. The Symbolic Jacobian Matrix.
# A rotation is represented by a quaternion q, as an element in
# projective 3-space:
> q := array(0..3):
> Q := sum(q['i']^2,'i'=0..3):
> Rhat := Matrix(3,3):
> R := Rhat/Q:
> Rhat[1,1] := q[0]^2 + q[1]^2 - q[2]^2 - q[3]^2:
> Rhat[2,2] := q[0]^2 - q[1]^2 + q[2]^2 - q[3]^2:
> Rhat[3,3] := q[0]^2 - q[1]^2 - q[2]^2 + q[3]^2:
> Rhat[1,2] := 2*(q[1]*q[2] + q[0]*q[3]):
> Rhat[2,1] := 2*(q[1]*q[2] - q[0]*q[3]):
> Rhat[1,3] := 2*(q[1]*q[3] - q[0]*q[2]):#
> Rhat[3,1] := 2*(q[1]*q[3] + q[0]*q[2]):
> Rhat[2,3] := 2*(q[2]*q[3] + q[0]*q[1]):
> Rhat[3,2] := 2*(q[2]*q[3] - q[0]*q[1]):
> R;
[ 2 2 2 2
[q[0] + q[1] - q[2] - q[3] , 2 q[1] q[2] + 2 q[0] q[3] ,
]
2 q[1] q[3] - 2 q[0] q[2]]
[ 2 2 2 2
[2 q[1] q[2] - 2 q[0] q[3] , q[0] - q[1] + q[2] - q[3] ,
]
2 q[2] q[3] + 2 q[0] q[1]]
[
[2 q[1] q[3] + 2 q[0] q[2] , 2 q[2] q[3] - 2 q[0] q[1] ,
2 2 2 2] /
q[0] - q[1] - q[2] + q[3] ] / (
/
2 2 2 2
q[0] + q[1] + q[2] + q[3] )
# Next we define the position vector pt:
> pt :=
:
# and we define the locations at endplate and the baseplate joints:
>
> for i from 1 to 6 do
> aa||i := :
> bb||i := :
> end do:
# Now we are ready to define the Jacobian matrix:
>
> jp := Matrix(6,6):
> for i from 1 to 6 do
> jpA||i :=
> LinearAlgebra:-ScalarMultiply(pt-bb||i,Q)+LinearAlgebra:-Multiply(Rhat
> ,aa||i);
> for j from 1 to 3 do
> jp[i,j] := jpA||i[j];
> end do;
> jpB||i :=
> LinearAlgebra:-CrossProduct(LinearAlgebra:-Multiply(Rhat,aa||i),pt-bb|
> |i);
> for j from 1 to 3 do
> jp[i,j+3] := jpB||i[j];
> end do;
> end do:
# We better do not ask here to compute the determinant...
>
# 2. Numerical Instances
# We have three vectors to specify: the position vector and the
# locations at endplate and the baseplate joints.
# For the purposes of applying homotopy continuation to factor the
# polynomials, we generate random complex number of modulus one, using
# 20 digit
>
# precision
> Digits:=20;
Digits := 20
> randcomp := proc()
> description `returns a random complex number of modulus one`;
> local angle;
> angle := stats[random,uniform[0,2*Pi]](1):
> evalc(cos(angle)+sin(angle)*I):
> end proc:
> randcomp();
-0.89780540910193948634 + 0.44039237889329905448 I
# Example 1
# In this first instance, we pick the locations and position vector at
# random:
> randomize:
> for i from 1 to 3 do
> p[i] := randcomp():
> end do:
> for i from 1 to 6 do
> for j from 1 to 3 do
> a||i[j] := randcomp();
> b||i[j] := randcomp();
> end do:
> end do:
> f1 := LinearAlgebra:-Determinant(jp):
> save f1, "verschelde/svw_exf1";
# Example 2
# Here we set the third element of all locations equal to zero:
> Digits:=20;
Digits := 20
> for i from 1 to 6 do
> a||i[3] := 0:
> b||i[3] := 0:
> end do:
> f2 := LinearAlgebra:-Determinant(jp):
> save f2,"verschelde/svw_exf2";
# Example 3
# In addition we set q[1] and q[2] equal to zero:
> Digits:=20;
Digits := 20
> q[1] := 0: q[2] := 0:
> f3 := LinearAlgebra:-Determinant(jp):
> f3:
> ff:=f3:
> save ff, "verschelde/svw_exf3";
# Example 4
# Finally, we let the point pt have arbitrary coordinates:
> Digits:=20;
Digits := 20
> for i from 1 to 3 do
> p[i] := evaln(p[i]);
> end do:
> pt := ;
[p[1]]
[ ]
pt := [p[2]]
[ ]
[p[3]]
> jp := Matrix(6,6);
> for i from 1 to 6 do
> jpA||i :=
> LinearAlgebra:-ScalarMultiply(pt-bb||i,Q)+LinearAlgebra:-Multiply(Rhat
> ,aa||i);
> for j from 1 to 3 do
> jp[i,j] := jpA||i[j];
> end do;
> jpB||i :=
> LinearAlgebra:-CrossProduct(LinearAlgebra:-Multiply(Rhat,aa||i),pt-bb|
> |i);
> for j from 1 to 3 do
> jp[i,j+3] := jpB||i[j];
> end do;
> end do:
[0 0 0 0 0 0]
[ ]
[0 0 0 0 0 0]
[ ]
[0 0 0 0 0 0]
jp := [ ]
[0 0 0 0 0 0]
[ ]
[0 0 0 0 0 0]
[ ]
[0 0 0 0 0 0]
> f4 :=LinearAlgebra:- Determinant(jp):
> save f4, "verschelde/svw_exf4";
>
>
>