|\^/|     Maple 15 (IBM INTEL LINUX)
._|\|   |/|_. Copyright (c) Maplesoft, a division of Waterloo Maple Inc. 2011
 \  MAPLE  /  All rights reserved. Maple is a trademark of
 <____ ____>  Waterloo Maple Inc.
      |       Type ? for help.
# ======================================================================== #
#
# Copyright 2012, Gabor Hegedus, Josef Schicho, and Hans-Peter Schroecker
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
# General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
#
# ======================================================================== #
> 
> read("bonds.mpl"):
> 
# define a spherical four-bar linkage
# ===================================
> h1 := [0, 1, 0, 0]:
> h2 := [0, 0, 1, 0]:
> h3 := [0, 0, 0, 1]:
> h4 := [0, 3/5, 4/5, 0]: 
> 
# sanity checks of input data
> for i from 1 to 4 do
>   if Mul(h||i, Cj(h||i)) <> [1, 0, 0, 0] then
>     print("Error: Input quaternions are not normalized!"):
>     stop:
>   end if:
>   if h||i + Cj(h||i) <> [0, 0, 0, 0] then
>     print("Error: Input quaternions are not rotations quaternions!"):
>     stop:
>   end if:
> od:
> 
# Ideal of configuration curve
# ============================
> 
> for i from 1 to 4 do:
>   u||i := [t||i, 0, 0, 0] - h||i:
> od:
> eqs := RVec(Mulall(u1, u2, u3, u4)):
> 
> es := {op(eqs[2..-1]), eqs[1]*u-1}:
> vars := [t1, t2, u, t3, t4]:
> Gb1 := Groebner[Basis]([op(es)], tdeg(t1, t2, u, t3, t4)):
> egb := select(T->not (u in indets(T)), Gb1):
> egb := Groebner[Basis]({op(egb)}, plex(t4, t3, t2, t1)):
> configcurve := allvalues(solve(egb, [t1, t2, t3, t4])):
> simplify(configcurve[1]);
               2         4        2      1/2
       5 - 5 t4  + (25 t4  - 14 t4  + 25)
[[t1 = -------------------------------------,
                       6 t4

              2             4        2      1/2
         -5 t4  - 5 + (25 t4  - 14 t4  + 25)
    t2 = --------------------------------------,
                          8 t4

                       2           4        2      1/2
                  25 t4    5 (25 t4  - 14 t4  + 25)
    t3 = - 7/24 + ------ - ---------------------------, t4 = t4]]
                    24                 24

> simplify(configcurve[2]);
                  2         4        2      1/2
         -5 + 5 t4  + (25 t4  - 14 t4  + 25)
[[t1 = - --------------------------------------,
                          6 t4

               2             4        2      1/2
           5 t4  + 5 + (25 t4  - 14 t4  + 25)
    t2 = - -------------------------------------,
                           8 t4

                       2           4        2      1/2
                  25 t4    5 (25 t4  - 14 t4  + 25)
    t3 = - 7/24 + ------ + ---------------------------, t4 = t4]]
                    24                 24

> 
# Compute bonds
# =============
> vars:=[t1, t2, t3, t4]:
> eqss := [op(egb), eqs[1]]:
> bonds := solve(eqss, vars):
> rr := simplify(convert(bonds, radical)):
> for r in rr do
>   print(r):
> od:
                     [t1 = -3 I, t2 = -I, t3 = -3, t4 = I]

                    [t1 = -1/3 I, t2 = I, t3 = 1/3, t4 = I]

                  [t1 = -I, t2 = -1, t3 = I, t4 = 4/5 + 3/5 I]

                 [t1 = -I, t2 = 1, t3 = -I, t4 = -4/5 + 3/5 I]

> stop:
memory used=12.1MB, alloc=10.1MB, time=0.35