views:

893

answers:

2

I'm trying to figure out how to use Mathematica to solve systems of equations where some of the variables and coefficients are vectors. A simple example would be something like

A + Vt = Pt

where I know A, V, and the magnitude of P, and I have to solve for t and the direction of P. (Basically, given two rays A and B, where I know everything about A but only the origin and magnitude of B, figure out what the direction of B must be such that it intersects A.)

Now, I know how to solve this sort of thing by hand, but that's slow and error-prone, so I was hoping I could use Mathematica to speed things along and error-check me. However, I can't see how to get Mathematica to symbolically solve equations involving vectors like this.

I've looked in the VectorAnalysis package, without finding anything there that seems relevant; meanwhile the Linear Algebra package only seems to have a solver for linear systems (which this isn't, since I don't know t or P, just |P|).

I tried doing the simpleminded thing: expanding the vectors into their components (pretend they're 3D) and solving them as if I were trying to equate two parametric functions,

Solve[ 
      { Function[t, {Bx + Vx*t, By + Vy*t, Bz + Vz*t}][t] == 
          Function[t, {Px*t, Py*t, Pz*t}][t],
        Px^2 + Py^2 + Pz^2 == Q^2 } , 
      { t, Px, Py, Pz } 
     ]

but the "solution" that spits out is a huge mess of coefficients and congestion. It also forces me to expand out each of the dimensions I feed it.

What I want is a nice symbolic solution in terms of dot products, cross products, and norms:

alt text

But I can't see how to tell Solve that some of the coefficients are vectors instead of scalars.

Is this possible? Can Mathematica give me symbolic solutions on vectors? Or should I just stick with No.2 Pencil technology?

(Just to be clear, I'm not interested in the solution to the particular equation at top -- I'm asking if I can use Mathematica to solve computational geometry problems like that generally without my having to express everything as an explicit matrix of {Ax, Ay, Az}, etc.)

+4  A: 

I don't have a general solution for you by any means (MathForum may be the better way to go), but there are some tips that I can offer you. The first is to do the expansion of your vectors into components in a more systematic way. For instance, I would solve the equation you wrote as follows.

rawSol = With[{coords = {x, y, z}},
  Solve[
    Flatten[
     {A[#] + V[#] t == P[#] t & /@ coords,
     Total[P[#]^2 & /@ coords] == P^2}],
    Flatten[{t, P /@ coords}]]];

Then you can work with the rawSol variable more easily. Next, because you are referring the vector components in a uniform way (always matching the Mathematica pattern v_[x|y|z]), you can define rules that will aid in simplifying them. I played around a bit before coming up with the following rules:

vectorRules =
  {forms___ + vec_[x]^2 + vec_[y]^2 + vec_[z]^2 :> forms + vec^2,
   forms___ + c_. v1_[x]*v2_[x] + c_. v1_[y]*v2_[y] + c_. v1_[z]*v2_[z] :>
     forms + c v1\[CenterDot]v2};

These rules will simplify the relationships for vector norms and dot products (cross-products are left as a likely painful exercise for the reader). EDIT: rcollyer pointed out that you can make c optional in the rule for dot products, so you only need two rules for norms and dot products.

With these rules, I was immediately able to simplify the solution for t into a form very close to yours:

  In[3] := t /. rawSol //. vectorRules // Simplify // InputForm
  Out[3] = {(A \[CenterDot] V - Sqrt[A^2*(P^2 - V^2) + 
                                   (A \[CenterDot] V)^2])/(P^2 - V^2), 
            (A \[CenterDot] V + Sqrt[A^2*(P^2 - V^2) + 
                                   (A \[CenterDot] V)^2])/(P^2 - V^2)}

Like I said, it's not a complete way of solving these kinds of problems by any means, but if you're careful about casting the problem into terms that are easy to work with from a pattern-matching and rule-replacement standpoint, you can go pretty far.

Pillsy
The two dot-product rules can be merged into one by replacing the `c_` in the first rule with `c_.` or `c:1`, where the period tells Mathematica to use the standard default for multiplication.
rcollyer
+5  A: 

With Mathematica 7.0.1.0

Clear[A, V, P];
A = {1, 2, 3};
V = {4, 5, 6};
P = {P1, P2, P3};
Solve[A + V t == P, P]

outputs:

{{P1 -> 1 + 4 t, P2 -> 2 + 5 t, P3 -> 3 (1 + 2 t)}}

Typing out P = {P1, P2, P3} can be annoying if the array or matrix is large.

Clear[A, V, PP, P];
A = {1, 2, 3};
V = {4, 5, 6};
PP = Array[P, 3];
Solve[A + V t == PP, PP]

outputs:

{{P[1] -> 1 + 4 t, P[2] -> 2 + 5 t, P[3] -> 3 (1 + 2 t)}}

Matrix vector inner product:

Clear[A, xx, bb];
A = {{1, 5}, {6, 7}};
xx = Array[x, 2];
bb = Array[b, 2];
Solve[A.xx == bb, xx]

outputs:

{{x[1] -> 1/23 (-7 b[1] + 5 b[2]), x[2] -> 1/23 (6 b[1] - b[2])}}

Matrix multiplication:

Clear[A, BB, d];
A = {{1, 5}, {6, 7}};
BB = Array[B, {2, 2}];
d = {{6, 7}, {8, 9}};
Solve[A.BB == d]

outputs:

{{B[1, 1] -> -(2/23), B[2, 1] -> 28/23, B[1, 2] -> -(4/23), B[2, 2] -> 33/23}}

The dot product has an infix notation built in just use a period for the dot.

I do not think the cross product does however. This is how you use the Notation package to make one. "X" will become our infix form of Cross. I suggest coping the example from the Notation, Symbolize and InfixNotation tutorial. Also use the Notation Palette which helps abstract away some of the Box syntax.

Clear[X]
Needs["Notation`"]
Notation[x_ X y_\[DoubleLongLeftRightArrow]Cross[x_, y_]]
Notation[NotationTemplateTag[
  RowBox[{x_,  , X,  , y_,  }]] \[DoubleLongLeftRightArrow] 
  NotationTemplateTag[RowBox[{ , 
RowBox[{Cross, [, 
RowBox[{x_, ,, y_}], ]}]}]]]
{a, b, c} X {x, y, z}

outputs:

{-c y + b z, c x - a z, -b x + a y}

The above looks horrible but when using the Notation Palette it looks like:

Clear[X]
Needs["Notation`"]
Notation[x_ X y_\[DoubleLongLeftRightArrow]Cross[x_, y_]]
{a, b, c} X {x, y, z}

I have run into some quirks using the notation package in the past versions of mathematica so be careful.

Davorak
In `Solve[A + V t == P, P]` you've missed the `t` multiplying `P`.
rcollyer
You type cross product with `:esc: cross :esc:`.
KennyTM
@rcollyerMy mistake, should work fine with the addition.@KennyTMThanks for pointing that out.
Davorak