PWC 248

PWC 248

I use PDL with Inline:Pdlpp for both challenges this week. I found these docs helpful: the perldocs on PDL::PP and PDL::API, the chapter on PP in the PDL book, and the Practical Magic guide to using PP.

Challenge 2 (submatrix sum)

We are given a matrix, say A, with M columns and N rows. We are asked to find a matrix, say B, with M-1 columns and N-1 rows, such that the entry for column i and row j of B, or B(i,j) using column-first notation is given by:

B(i,j) = A(i,j)+A(i,j+1)+A(i+1,j+1)+A(i+1,j)

With PP, we write XS code which uses a traditional pair of nested loops to create B. Here is the PP definition:

pp_def('submatrix_sum',
    Pars => 'a(m,n); [o]b(p,q);',
    RedoDimsCode => '$SIZE(p)=$SIZE(m)-1; $SIZE(q)=$SIZE(n)-1;',
    Code => q{
        PDL_Indx i, j, m_size, n_size;
        m_size=$SIZE(m)-1;
        n_size=$SIZE(n)-1;
        for (i=0; i < m_size; i++) {
            for (j=0; j < n_size; j++) {
                $b(p=>i, q=>j)=
                $a(m=>i, n=>j)+$a(m=>i+1, n=>j+1)+
                $a(m=>i, n=>j+1)+$a(m=>i+1, n=>j);
            }
        }
    },
);

The Pars section specifies the function signature. The submatrix_sum PDL method operates on a 2D PDL labelled a, with the two dimensions labelled m and n respectively. The output is also a 2D PDL, labelled b, with the two dimensions labelled p and q respectively.

We need to tell the compiler the size of B along dimensions p and q. We do this in the RedoDimsCode section.

Finally comes the main Code section. This is pretty much vanilla C with some additional custom preprocessor directives.

Challenge 1 (Shortest distance)

We are given a string say str and a single character say chr. We first find all the positions in the string where the character is chr. We then output an array showing for every position in the string, how far it is from the nearest occurrence of chr. For example, if the string is 'abaaaba' and the character is 'b', the output would be: [1,0,1,2,1,0,1]. The intermediate array showing the positions matching 'b' is [1,5]. (This is Perl, so arrays are zero-indexed).

My original plan was to set up two pp_def entries, one for an inner loop function that computes as output the sub-array showing the nearest elements, given as input a pair of consecutive items matching chr. For example, following the example above, the input could be the pair [1,5], and the output the sub-array 0,1,2,1,0 at positions 1 through 5 of the final output array. 

I wanted to call this from another outer-loop function also defined in PP which calls the inner-loop function. This turned out to be very tricky and I was not able to get it done in time. I ended up just writing an inner-loop function and calling it from a PDL (Perl 5) outer loop, which is not very efficient. I also however wrote a single PP function that includes both the outer and inner loop, in a different script labelled ch-1a.pl.

Here is the pp_def for the original standalone inner-loop function.

pp_def(
    'between',
    Pars => 'indx a1(); indx a2(); int [o]b(n);',
    RedoDimsCode => '$SIZE(n)=$a2()-$a1()+1;',
    Code => q{
        PDL_Indx i,c_mid, n_size=$SIZE(n);
        if (0==(n_size % 2)){
           c_mid=n_size/2;
        }
        else {
            c_mid=(n_size+1)/2;
        }
       for (i=0; i < c_mid; i++) {
            $b(n=>i)=i;             
       }
       for (i=c_mid; i < n_size; i++) {
            $b(n=>i)=n_size-i-1;
       }
    },
);

The logic reflects that if the gap between the two input indices leads to an odd-number-sized sub-array of distances, then the sub-array of distances would be of the pattern:

0,1,2,...x,x-1,...2,1,0

For example, for input-pair [1,5], we have [0,1,2,1,0].

If the gap leads to an even-number-sized sub-array of distances, then the sub-array of distances would be of the pattern:

0,1,2,...x-1,x,x,x-1,...2,1,0

For example, for input-pair [1,6], we have [0,1,2,2,1,0].

The XS function 'between' is called in PDL as a method, for example: pdl(1)->between(6).

I call it from a routine written in PDL (Perl 5):

local *shortest_distance=sub {
    my ($str,$chr)=@_;
    ($str =~ /$chr/) || return 0;

    (1==length($chr)) || (die "Single character input needed:$!\n");

    my $retval=zeros(length($str)); #-- pdl to hold return value
    
    #-- get indices where the element equals chr
    my @grp=(grep {$chr eq substr($str,$_,1)} 0..length($str));
    
    #-- if the initial matching index is greater than 0 say x, fill 0:x
    ($grp[0] > 0) && ($retval->slice("0:$grp[0]") .= sequence($grp[0]+1)->slice('-1:0'));
    
    #-- fill in distances between elements of @grp
    map {
        $retval->slice("$grp[$_]:$grp[$_+1]") .=
        pdl($grp[$_])->between(pdl($grp[$_+1]));
    } 0 .. ($#grp-1);
    
    
    #-- if the last element of @grp is less than the length of $str, fill in the distance for the remaining elements.
    ($grp[-1] < (length($str)-1)) && ($retval->slice("$grp[-1]:") .= sequence(length($str)-$grp[-1]));
    
    return $retval;
 };

The array @grp contains the positions in str that match chr.

We loop through elements of @grp, calling the between operator at each stage. We also directly construct in Perl 5 the distances before the first index in @grp, as well as after the last index in @grp.

Potential improvements

If this was production code, I would tinker with it further, but given our one-week constraint, I will confine myself to reflecting on improvements in this blog.

Routines defined in PP should be capable of broadcasting to higher dimensions. My 'between' routine does mechanically broadcast, but not in a useful way. A better design could be to return a vector of the same length as the original string, with the distances inserted for the section between the input indices. For example, for input indices [1,6] with string length 10, the output would be [0,0,1,2,2,1,0,0,0,0,0], not just [0,1,2,2,1,0]. We could then say pass the routine more than one set of indices, say [[1,6],[6,8]] with the same  string length, say 10, to get an output like

[
[0,0,1,2,2,1,0,0,0,0,0],
[0,0,0,0,0,0,0,1,0,0,0] 
].
 
One can then just call ...->xchg(0,1)->sumover to get
 
[0,0,1,2,2,1,0,1,0,0,0].
 
Of course, this needs the pp_def to change the input definition (Pars) since we have to input the string length too. RedoDimsCode will also have to change.
 
A more efficient approach  

A more efficient approach is the one in ch-1a.pl, where I use a single pp_def function to take a list of indices and a string length as inputs, and return the shortest-distance array as output. Here is that single pp_def. 

pp_def(
    'shortest_distance_pp',
    Pars => 'indx a(m); indx lenstr(); [o]b(n);',
    RedoDimsCode => '$SIZE(n)=$lenstr();',
    GenericTypes => ['N'],
    Code => q{
        PDL_Indx i,j,size_m, size_n, ctr_stop;
        size_m=$SIZE(m);
        size_n=$SIZE(n);
        ctr_stop=$a(m=>0);
        if (0 < $a(m=>0)){
            for (i=0; i < ctr_stop; i++) {
                $b(n=>i)=ctr_stop-i;
            }
        }
        PDL_Indx a1,a2,c_mid,gap_size;
        for (j=0; j<size_m; j++) {
            a1=$a(m=>j);
            a2=$a(m=>j+1);
            gap_size=a2-a1+1;
            
            if (0==(gap_size % 2)){
                c_mid=gap_size/2;
            }
            else {
                c_mid=(gap_size+1)/2;
            }
            
            for (i=a1; i < a1+c_mid; i++) {
                $b(n=>i)=i-a1;             
            }
            for (i=a1+c_mid; i <= a2; i++) {
                $b(n=>i)=a2-i;
            }
            
        }       
        if (size_n > $a(m=>-1)+1){
            a1=$a(m=>size_m-1);
            for (i=a1; i < size_n; i++) {
                $b(n=>i)=i-a1;
            }
        }                        
    },
);

It follows the same logic as the Perl 5 outer loop + pp_def inner loop in ch-1.pl. But it would be faster and more scalable. Confining the input types to 'N' via the GenericTypes field ('N' refers to the PDL index type, usually a signed integer) reduces the work for the compiler which would otherwise create executable code unnecessarily for things like doubles etc.

A couple of very nice solutions from other contributors also using PDL: Jo S (blog post here) with a very compact almost one-liner solution that uses PDL even for the character handling; and W. Luis Mochán (blog post here) with a different very compact solution that is presented as both one-liner and short script.

Comments

Post a Comment

Popular posts from this blog

PWC 258

PWC 258 Last week (#257) I finished challenge 2 and my blog post for challenge 257 late last week, so had only submitted challenge 1. Here is my blog post for PWC 257 , including the full code for my challenge 2 solution. This week, I have done both challenges in perl 4.019 on DOSBOX. Both scripts also work fine on the latest Perl 5.38. I also did the challenges in Raku, where they are one-liners.  Challenge 1 (Count Even Digits Number) We are given a list of integers, and asked to count the number of elements whose base-10 representation has an even number of digits. In Raku , this is a one-liner. Here it is: raku -e 'say +@*ARGS.grep({$_.chars %% 2})' $@ It could be a one-liner in Perl 4 too, but since I run it on DOSBOX, my script would not be a bash shell script. Here is the short subroutine that does the job.  sub count_even_digits {  scalar( grep( ((length($_) % 2) == 0), @_) ); } The full code is here . Challenge 2 (Sum of values) We are given an integer k  an...
Read more

PWC 253

PWC 253 I do both challenges this week using old software on DOSBOX: perl 4.019, and for comparison, also python 1.4 beta. My solutions also work on the current versions of Perl 5 (5.38) and Python 2 (2.7.18) respectively. Challenge 1 (Split Strings) Links to code: perl4.019 python1.4b We are given a list of strings, and a single character say sep . We are asked to return the list, with the strings split into component substrings using sep as a separator. We should not include any empty strings in the output. In perl 4, I first stick a backslash in front of the separator character, so that it is escaped if entered in a regex. This is because the perl split function takes a regex as its first argument indicating the separator. If we have a character like '$' or ',', we want to escape it so that it is read literally when parsing the regex. Here is the key subroutine:   sub split_strings {  local ($separator, @words)=@_;  local (@retval)=(); #array to be returned  # ...
Read more

PWC 249

PWC 249 Challenge 1 (Equal pairs) We are given a list of integers with an even number of elements. We are asked to check if the list can be completely divided into sub-lists of size 2 with each sub-list containing the same integer repeated, like say (2,2). If so, we are asked to return the list of such sub-lists. This is straightforward in Perl 5. I use modern conveniences this time, though it could be done using only primitive Perl 4 syntax. Here is my subroutine. It returns a list of array references. local *equal_pairs = sub {     #-- return () if input is of odd length     ((scalar(@_) % 2) > 0) && (return ());          #-- count number of occurrences of each element using a hash     local %_ = ();     map {$_{$_}++} @_;          #-- if any element occurs an odd number of times, return ()     (grep {($_{$_} % 2) > 0} keys %_) &&...
Read more