Craig
thanks, I somehow missed the RedoDimsCode key (it is not mentioned in
the Book) when browsing docs. That's the missing step.

Ingo

On 27.05.19 18:30, Craig DeForest wrote:
Ingo —

I had a look at your code.  I see a couple minor things there (e.g.
you don’t need to declare “PDL_Indx m,l” in that context since they’re
predeclared by the Pars field) — but I don’t see anything obvious that
would suck up memory.

I looked up Hankel matrices on Wikipedia and found a basic definition.
Here’s a computed inline PP function that generates hankel matrices.
 It accepts a source linear array and Hankelizes it into a square
matrix, rounding the input size down to the next odd number.  It uses
the RedoDims mechanism to set the size of the output array, based on
the parameters in the input array.  (Note that NONE of the $SIZE
macros are set yet when RedoDimsCode gets called, so you can’t say
“$SIZE(w)=($SIZE(n)+1)/2” — you have to dive into the input PDL
structure as shown).

Does that help?

Cheers,
Craig

—— snip snip --

=head2 hankel

=for ref

Takes an input array and turns it into a square Hankel matrix of the
appropriate size for the array.
The input should have an odd number of values.  If it has an even
number, then it is rounded down to the
next odd value below.

=cut

*hankel = \&PDL::hankel;

no PDL::NiceSlice;
use Inline Pdlpp=><<'EOF';

pp_def('hankel',
Pars=>'source(n); [o]hnk(w,h)',
RedoDimsCode => <<'EORDC',
$SIZE(w) = ($PDL(source)->dims[0]+1)/2;
$SIZE(h) = ($PDL(source)->dims[0]+1)/2;
EORDC
Code=><<'EOC',
loop(h) %{
 loop(w) %{
    $hnk() = $source(n=>h+w);
                  %}
%}
EOC
);
EOF


On May 27, 2019, at 3:16 AM, Ingo Schmid <[email protected]
<mailto:[email protected]>> wrote:

Hi,

I've tried to create a PP function to calculate the Hankel matrix
from a vector. I used Inline pdlpp. Whenever I run the function, I
get an out of memory error. I tried explicit for loops as well. There
is no good example in the book or the docs where output dimensions
are calculated based on input. In this case l depends on n and m.

pp_def ('fid_to_hankel2',
        Pars=>'v(n); [o]hankel(l,m)',
        OtherPars=>'int ms => m',
        Code => pp_line_numbers(__LINE__, q{
                PDL_Indx ns = $SIZE(n);
                PDL_Indx m,l;
                l=ns-m+1;
                printf ("size n %d , l %d , m %d\n",ns,l,m);
                loop(n) %{
                loop(l) %{
                        loop(m) %{
$hankel()=$v(n=>l+m);
                                printf ("n %d = l %d + m %d\n",l+m,l,m);
                        %}
                %}
                %}
        }), # line_numbers
        );

I know that the Hankel matrix can be created efficiently using xvals,
yvals and index.

sub hankel {
        my $fid=shift;
        my $rank=shift;
        my $size=($rank,$fid->dim(0)-$rank+1);
        my $ind=xvals($size)+yvals($size);
        $fid->index($ind);
}

Nonetheless I would like to understand what's going on. Btw., that's
with pdl 2.019 on linux 64 bit. Other PP functions work alright.

Best wishes

Ingo


_______________________________________________
pdl-general mailing list
[email protected]
<mailto:[email protected]>
https://lists.sourceforge.net/lists/listinfo/pdl-general


_______________________________________________
pdl-general mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/pdl-general

Reply via email to