Jan Hoogenraad wrote:
> Victory.
Congratulations and happy Pldpp-ing!
Some suggestions below for the PDL testing...
> Attached the first lapack call I have compiled on an Ubunto linux machine.
> Dependencies:
> - f2c (for the f2c.h header file)
> - lapack (compiled package)
> - clapack.h (from internet, requested to be included in package)
> - C-compiler
> - PDL
>
> I'll try to convert a few more routines in the coming days.
>
> -----------------------------
> use blib; # when using inside the dist tree
> use PDL; # this must be called before (!) 'use Inline Pdlpp' calls
>
> # Example of call to lapack routine
>
> use Inline Pdlpp => Config =>
> INC => "-I$ENV{HOME}/include -I../../..",
> LIBS => "-L$ENV{HOME}/lib -llapack -lm",
> AUTO_INCLUDE => <<'EOINC',
> #include <math.h>
> #include "f2c.h" /* type definitions for a files */
> #include "clapack.h" /* Definitions of all routines */
>
> EOINC
> BOOT => '';
>
>
> use Inline Pdlpp; # the actual code is in the __Pdlpp__ block below
>
> test_dgsev();
>
>
>
> sub test_dgsev() {
>
> $n= pdl 4;
Perl scalar work fine here. e.g. $n = 4 should
be ok.
> $A = PDL->zeroes(double ,$n,$n);
>
> for $i (0.. $n-1 ) {
> for $j (0.. $n-1 ) {
> set $A,$i,$j , pow($i,$j)+$n*$j+1;
> }
> }
Looping through piddles element-by-element is
a slow process. For better performance, compactness,
and legibility (once you are familiar with the PDL
idioms) use the vectorized/threaded PDL ops. The
above can be calculated equivalently as:
$A .= $A->xvals**$A->yvals + 4*$A->yvals + 1;
The set() operation is mainly for debugging purposes
and has been largely replaced by slicing ops, especially
PDL::NiceSlice syntax.
> $oneone = pdl 1;
> $b = PDL->zeroes(double ,$n,$oneone);
> for $i (0.. $n-1 ) {
> set $b,$i,0, pow(2,$i+1);
> }
Can be calculated as:
$b .= 2**($b->xvals + 1);
> $ipiv = PDL->zeroes(long ,$n);
> $info = pdl -1;
> $qmat="\n\[\n\ \[\ 2\ \ 2\ \ 2\ \ 2\]\n\ \[\ 5\ \ 6\ \ 7\ \ 8\]\n\ \[\
> 9\ 10\ 13\ 18\]\n\ \[13\ 14\ 21\ 40\]\n\]\n";
This is easier to read if constructed with pdl:
$qmat = pdl [ [ 2, 2, 2, 2, ], [ 5, 6, 7 , 8, ], [ 9, 10, 13,18, ], [ 13, 14,
21, 40, ] ];
> $qmbt="\n\[\n\ \[\ 2\ \ 4\ \ 8\ 16\]\n\]\n";
> $qmipivt="\[0\ 0\ 0\ 0\]";
> $qma=("".$A);
> $qmb=("".$b);
> $qmipiv=("".$ipiv);
> print "$qma <> $qmat\n" if($qmat ne $qma);
Then you can just use $A != $qmat to test numerically
or approx to allow for floating point epsilon as in:
print "\$qmat == \$A\n" if all approx $qmat, $A;
since PDL overloads the comparison operators as
well as the arithmetic ones.
Cheers,
Chris
> print "$qmb <> $qmbt\n" if($qmbt ne $qmb);
> print "$qmipiv <> $qmipivt\n" if($qmipivt ne $qmipiv);
>
>
> PDL::dgsev($A,$ipiv,$b,$info);
>
>
> print "Error code returned: $info\n" if ($info);
>
> $qmbt="\n\[\n\ \[\ \ \ \ \-5\.3333333\ \ \ \ \ \ 1\.6666667\
> \-7\.4051008e\-17\ \ \ \ \ 0\.33333333\]\n\]\n";
> $qmipivt="\[1\ 4\ 4\ 4\]";
>
> $qmb=("".$b);
> $qmipiv=("".$ipiv);
>
> print "$qmb <> $qmbt\n" if($qmbt ne $qmb);
> print "$qmipiv <> $qmipivt\n" if($qmipivt ne $qmipiv);
>
> print "End of test\n"
>
> }
>
> __DATA__
>
> __Pdlpp__
>
> # Call dgsev
> pp_def('dgsev',
> Pars => '[o] a(lda,n); int [o] ipiv(n); [o] b(ldb,nrhs); int [o]
> info();',
> GenericTypes => [D],
> Code => 'int ldas = $SIZE(lda); int ns = $SIZE(n);
> int ldbs =$SIZE(ldb);int nrhss =$SIZE(nrhs);
> dgesv_(&ns, &nrhss, $P(a), &ldas, $P(ipiv), $P(b), &ldbs,
> $P(info) );',
> );
>
>
> -----------------------------
>
> Jan Hoogenraad wrote:
>> Please discard my previous posting. I now see I missed brackets in the
>> Pars section.
>>
>> I've now got the file to compile, and now get the error
>> Error in Pdgsev_:Wrong dims
>> at the moment of the call with the attached perl file.
>>
>> Apparently, I miss out on something with the dimensions...
>>
>
>
> _______________________________________________
> Perldl mailing list
> [email protected]
> http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
>
>
> ------------------------------------------------------------------------
>
>
> No virus found in this incoming message.
> Checked by AVG - www.avg.com
> Version: 8.5.339 / Virus Database: 270.12.48/2148 - Release Date: 06/01/09
> 06:09:00
>
_______________________________________________
Perldl mailing list
[email protected]
http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
