MathGroup Archive 1997

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Q: About constants inside Distribute

  • To: mathgroup at smc.vnet.net
  • Subject: [mg9258] Re: Q: About constants inside Distribute
  • From: Robert Villegas <villegas>
  • Date: Fri, 24 Oct 1997 01:01:13 -0400
  • Organization: Wolfram Research
  • Sender: owner-wri-mathgroup at wolfram.com

> The rule:
> 
>    rule = Dot[s___,(k:(a[i_]))m___, r___]:> k Dot[s,m,r]
> 
> works fine separating constants (a[i_]) in the distribution of
> a Dot product:
> 
>    f = O1*a[1] + O2*a[2]
>    g = P1*a[3] 
>    h = O2*a[3] 
> 
> In[44]:= Distribute[f.g.h]//.rule
> 
>                   2                         2
> Out[44]= a[1] a[3]  O1 . P1 . O2 + a[2] a[3]  O2 . P1 . O2
> 
> 
> However, the rule does not work if product of constants (a[i_]*a[j_] )
> are included. For example, taking:
> 
>    g = P1*a[3]*a[1]
> 
> 
> In[45]:= g = P1*a[3]*a[1]
> 
> Out[45]= P1 a[1] a[3]
> 
> In[46]:= Distribute[f.g.h]//.rule
> 
>              2
> Out[46]= a[1]  a[3] O1 . P1 . a[3] . O2 + a[1] a[2] a[3] O2 . P1 . a[3] . O2


In case you're impatient for a working rule, here's a solution:

   ruleNew = Dot[s___,(k:(a[i_])) m_, r___]:> k Dot[s, m, r]

The only change is to make m a single object _, not a sequence __,
taking advantage of a special property of Times.  It's very
instructive, though, to see why your m__ didn't work and how it can be
made to work, so here is some debugging and explanation below.


   Since "//. rule" repeatedly performs "/. rule" until the result is
the same for two consecutive stages, let's apply "/. rule" one step at
a time and examine each result with FullForm to see exactly where this
goes awry. We'll confine ourselves to the first term of
Distribute[f.g.h] since the problem reveals itself there.

First, the definitions:

In[1]:= {f, g, h} = {O1*a[1], P1*a[1]*a[3], O2*a[3]} 

Out[1]= {O1 a[1], P1 a[1] a[3], O2 a[3]}

In[2]:= dot = Dot[f, g, h] 

Out[2]= (O1 a[1]) . (P1 a[1] a[3]) . (O2 a[3])

In[3]:= rule = Dot[s___,(k:(a[i_]))m___, r___]:> k Dot[s,m,r]

Out[3]= (s___) . ((k:a[i_]) (m___)) . (r___) :> k s . m . r


(* Before rules are applied: *)

In[4]:= FullForm[dot]

Out[4]//FullForm= Dot[Times[O1, a[1]], Times[P1, a[1], a[3]], Times[O2,
a[3]]]


(* Application 1 of rule *)

In[5]:= FullForm[dot /. rule]

Out[5]//FullForm= Times[a[1], Dot[O1, Times[P1, a[1], a[3]], Times[O2,
a[3]]]]

The first Times[...] subexpression has dissolved, its a[1] taken outside
and its O1 becoming an element of Dot.


(* Application 2 of rule *)

In[6]:= FullForm[% /. rule]

Out[6]//FullForm= Times[Power[a[1], 2], Dot[O1, P1, a[3], Times[O2,
a[3]]]]

The second Times[...] subexpression dissolves, its a[1] taken outside,
and its other two elements, P1 and a[3], becoming elements of Dot. And
that's the problem!  Since P1 and a[3] were elements of a Times, they
should stay that way.


(* Application 3 of rule *)

In[7]:= FullForm[% /. rule]

Out[7]//FullForm= Times[Power[a[1], 2], a[3], Dot[O1, P1, a[3], O2]]

(* Application 4 of rule *)

In[8]:= FullForm[% /. rule]

Out[8]//FullForm= Times[Power[a[1], 2], a[3], Dot[O1, P1, a[3], O2]]

Step 3 finishes the substitutions (step 4 effects no change) by
dissolving the last Times[...] subexpression, pulling the a[3] outside
and leaving the O2 as an element of Dot.  No surprises at this stage
and nothing wrong.


So we've discovered that the error is in the way this term from the
left-hand side

   (k:(a[i_]))m___

is used on the right-hand side.  We forgot that m might represent two or
more elements of Times that must retain the head Times on the
right-hand, otherwise they'll just dissolve into the whole Dot
expression.  When m represents only one element, it's not a problem to
wrap it in Times on the right-hand side, because it will evaluate and
lose the head Times.


Here is the step-by-step procedure done with the new rule, which works
fine:

In[9]:= ruleNew = Dot[s___,(k:(a[i_]))m___, r___]:> k Dot[s, Times[m],
r]

Out[9]= (s___) . ((k:a[i_]) (m___)) . (r___) :> k s . Times[m] . r

In[10]:= FullForm[dot /. ruleNew]

Out[10]//FullForm= 
   Times[a[1], Dot[O1, Times[P1, a[1], a[3]], Times[O2, a[3]]]]

In[11]:= FullForm[% /. ruleNew]

Out[11]//FullForm= 
   Times[Power[a[1], 2], Dot[O1, Times[P1, a[3]], Times[O2, a[3]]]]

In[12]:= FullForm[% /. ruleNew]

Out[12]//FullForm= Times[Power[a[1], 2], a[3], Dot[O1, P1, Times[O2,
a[3]]]]

In[13]:= FullForm[% /. ruleNew]

Out[13]//FullForm= Times[Power[a[1], 2], Power[a[3], 2], Dot[O1, P1,
O2]]

In[14]:= FullForm[% /. ruleNew]

Out[14]//FullForm= Times[Power[a[1], 2], Power[a[3], 2], Dot[O1, P1,
O2]]


   This example is instructive because it exhibits a common oversight
when using the sequence patterns

  p__        (* BlankSequence *)
  p___       (* BlankNullSequence *)
  p..        (* Repeated *)
  p...       (* RepeatedNull *)

The oversight is forgetting that the value of p on the right-hand side
will be a headless pile of elements that dissolve into the surrounding
head -- it is not a separate subexpression any more!  Any time you use
one of these, think about how you use it on the right-hand side: is it
okay if the sequence's elements get separated and become arguments of
the enclosing expression?  Or should you enclose them in a separate
head to keep them together?  That is, think

    f[p__]  :>  h[a, p, b]    (* or: *)

    f[p__]  :>  h[a, g[p], b]


   By the way, in your case we happen to be dealing with Times, which is
a Flat head so the pattern-matcher automatically groups its elements if
needed to make a match.  If we avoid sequence patterns and use m_, then
m will consume as many factors as needed under an extra head of Times,
so its value on the right side will be something like Times[a[1], O1,
a[3], ...]

In[21]:= rule3 = Dot[s___,(k:(a[i_])) m_, r___]:> k Dot[s, m, r]

Out[21]= (s___) . ((k:a[i_]) (m_)) . (r___) :> k s . m . r

In[22]:= dot //. rule3

             2     2
Out[22]= a[1]  a[3]  O1 . P1 . O2


Robby Villegas


  • Prev by Date: Re: trig expansion
  • Next by Date: tessalation of a sphere
  • Previous by thread: Re: Q: About constants inside Distribute
  • Next by thread: Help - Mathematica bombs out on loading <<Calculus`LaplaceTransform`