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]> 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]
> https://lists.sourceforge.net/lists/listinfo/pdl-general
_______________________________________________
pdl-general mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/pdl-general
