Re: Compiling SingularValueDecomposition
- To: mathgroup at smc.vnet.net
- Subject: [mg89855] Re: Compiling SingularValueDecomposition
- From: Jean-Marc Gulliet <jeanmarc.gulliet at gmail.com>
- Date: Mon, 23 Jun 2008 02:45:32 -0400 (EDT)
- Organization: The Open University, Milton Keynes, UK
- References: <g3ihod$fc9$1@smc.vnet.net> <g3kunv$7e2$1@smc.vnet.net>
Jean-Marc Gulliet wrote:
> Frank Hu wrote:
>> Hi, group,
>>
>> I'm trying to speed up a piece of code that uses
>> SingularValueDecomposition by compiling it. But I couldn't get it to
>> compile. For a simple demonstration, try
>>
>> fx=Compile[{{x, _Real, 2}}, SingularValueDecomposition[x],
>> {{SingularValueDecomposition, _Real, 3}}]
>>
>> and fx[[4]] is
>> {{1, 5}, {54, Function[{x}, SingularValueDecomposition[x]], 3, 2, 0,
>> 3, 2, 1}, {2}}
>>
>> The "Function" there tells fx wasn't compiled successfully. Calling
>> fx[{{1., 2.}, {3., 4.}}] will generate the following warnings
>>
>> CompiledFunction::cfte: Compiled expression {<<1>>} should be a rank
>> 2 tensor of machine-size real numbers. >>
>> CompiledFunction::cfex: Could not complete external evaluation at
>> instruction 2; proceeding with uncompiled evaluation. >>
>
> Hi Frank,
>
> You must write the intermediate functions as full expressions, that is
> SingularValueDecomposition[_] rather than just
> SingularValueDecomposition. The following works as expected:
>
>
> In[1]:= fx = Compile[{{x, _Real, 2}}, SingularValueDecomposition[x],
> {{SingularValueDecomposition[_], _Real, 3}}]
> fx[[4]]
> fx[{{1., 2.}, {3., 4.}}]
>
> Out[1]= CompiledFunction[]
>
> Out[2]= {{1, 5}, {54, Function[{x}, SingularValueDecomposition[x]], 3,
> 2, 0, 3,
>
> 3, 1}, {2}}
>
> Out[3]= {{{-0.404554, 0.914514}, {-0.914514, -0.404554}},
>
> {{5.46499, 0.}, {0., 0.365966}},
>
> {{-0.576048, -0.817416}, {-0.817416, 0.576048}}}
Note that the code above will failed to compiled when the input is not a
*square* matrix because Mathematica expects the returned value to be a
*full* array. (In a full array all parts at a particular level must be
lists of the same length.) When you feed the function fx with a
rectangular (non square) matrix, the result returned by
SingularValueDecomposition is a list of unequal-length elements.
In[1]:= fx = Compile[{{x, _Real, 2}},
SingularValueDecomposition[x],
{{SingularValueDecomposition[_], _Real, 3}}]
f1 = fx[{{1., 2.}, {3., 4.}}]
f2 = fx[{{1., 2.}, {3., 4.}, {5., 6.}}]
f3 = fx[{{1., 2., 3.}, {4., 5., 6.}}]
({ArrayDepth[#1], ArrayQ[#1]} & ) /@ {f1, f2, f3}
Out[1]= CompiledFunction[]
Out[2]= {{{-0.404554, 0.914514}, {-0.914514, -0.404554}},
{{5.46499, 0.}, {0., 0.365966}},
{{-0.576048, -0.817416}, {-0.817416, 0.576048}}}
During evaluation of In[1]:= CompiledFunction::cfte:Compiled expression
{<<1>>} should be a rank 3 tensor \
of machine-size real numbers. >>
During evaluation of In[1]:= CompiledFunction::cfex:Could not complete
external evaluation at instruction \
2; proceeding with uncompiled evaluation. >>
Out[3]= {{{-0.229848, 0.883461, 0.408248}, {-0.524745, 0.240782,
-0.816497},
{-0.819642, -0.401896, 0.408248}},
{{9.52552, 0.}, {0., 0.514301}, {0., 0.}},
{{-0.619629, -0.784894}, {-0.784894, 0.619629}}}
During evaluation of In[1]:= CompiledFunction::cfte:Compiled expression \
{{{-0.386318,-0.922366},{-0.922366,0.386318}},{{<<1>>},<<1>>},<<1>>}
should \
be a rank 3 tensor of machine-size real numbers. >>
During evaluation of In[1]:= CompiledFunction::cfex:Could not complete
external evaluation at instruction \
2; proceeding with uncompiled evaluation. >>
Out[4]= {{{-0.386318, -0.922366}, {-0.922366, 0.386318}},
{{9.50803, 0., 0.}, {0., 0.77287, 0.}},
{{-0.428667, 0.805964, 0.408248}, {-0.566307, 0.112382, -0.816497},
{-0.703947, -0.581199, 0.408248}}}
Out[5]= {{3, True}, {1, False}, {1, False}}
On the other hand, if we use a function such as SingularValueList, which
returns a one-dimensional vector, the compilation works with rectangular
inputs.
In[6]:= gx = Compile[{{x, _Real, 2}},
SingularValueList[x],
{{SingularValueList[_], _Real, 1}}]
g1 = gx[{{1., 2.}, {3., 4.}}]
g2 = gx[{{1., 2.}, {3., 4.}, {5., 6.}}]
g3 = gx[{{1., 2., 3.}, {4., 5., 6.}}]
({ArrayDepth[#1], ArrayQ[#1]} & ) /@ {g1, g2, g3}
Out[6]= CompiledFunction[]
Out[7]= {5.46499, 0.365966}
Out[8]= {9.52552, 0.514301}
Out[9]= {9.50803, 0.77287}
Out[10]= {{1, True}, {1, True}, {1, True}}
HTH,
-- Jean-Marc