diff --git a/scikits/learn/src/cblas/ATL_srefcopy.c b/scikits/learn/src/cblas/ATL_srefcopy.c new file mode 100644 index 0000000000000000000000000000000000000000..fec830e4fc9ea128c6cd48cfd79b3638c71d9acc --- /dev/null +++ b/scikits/learn/src/cblas/ATL_srefcopy.c @@ -0,0 +1,148 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflevel1.h" + +void ATL_srefcopy +( + const int N, + const float * X, + const int INCX, + float * Y, + const int INCY +) +{ +/* + * Purpose + * ======= + * + * ATL_srefcopy copies the entries of an n-vector x into an n-vector y. + * + * Arguments + * ========= + * + * N (input) const int + * On entry, N specifies the length of the vector x. N must be + * at least zero. Unchanged on exit. + * + * X (input) const float * + * On entry, X points to the first entry to be accessed of an + * incremented array of size equal to or greater than + * ( 1 + ( n - 1 ) * abs( INCX ) ) * sizeof( float ), + * that contains the vector x. Unchanged on exit. + * + * INCX (input) const int + * On entry, INCX specifies the increment for the elements of X. + * INCX must not be zero. Unchanged on exit. + * + * Y (input/output) float * + * On entry, Y points to the first entry to be accessed of an + * incremented array of size equal to or greater than + * ( 1 + ( n - 1 ) * abs( INCY ) ) * sizeof( float ), + * that contains the vector y. On exit, the entries of the in- + * cremented array X are copied into the entries of the incre- + * mented array Y. + * + * INCY (input) const int + * On entry, INCY specifies the increment for the elements of Y. + * INCY must not be zero. Unchanged on exit. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float x0, x1, x2, x3, x4, x5, x6, x7; + float * StX; + register int i; + int nu; + const int incX2 = 2 * INCX, incY2 = 2 * INCY, + incX3 = 3 * INCX, incY3 = 3 * INCY, + incX4 = 4 * INCX, incY4 = 4 * INCY, + incX5 = 5 * INCX, incY5 = 5 * INCY, + incX6 = 6 * INCX, incY6 = 6 * INCY, + incX7 = 7 * INCX, incY7 = 7 * INCY, + incX8 = 8 * INCX, incY8 = 8 * INCY; +/* .. + * .. Executable Statements .. + * + */ + if( N > 0 ) + { + if( ( nu = ( N >> 3 ) << 3 ) != 0 ) + { + StX = (float *)X + nu * INCX; + + do + { + x0 = (*X); x4 = X[incX4]; x1 = X[INCX ]; x5 = X[incX5]; + x2 = X[incX2]; x6 = X[incX6]; x3 = X[incX3]; x7 = X[incX7]; + + *Y = x0; Y[incY4] = x4; Y[INCY ] = x1; Y[incY5] = x5; + Y[incY2] = x2; Y[incY6] = x6; Y[incY3] = x3; Y[incY7] = x7; + + X += incX8; + Y += incY8; + + } while( X != StX ); + } + + for( i = N - nu; i != 0; i-- ) + { + x0 = (*X); + *Y = x0; + + X += INCX; + Y += INCY; + } + } +/* + * End of ATL_srefcopy + */ +} diff --git a/scikits/learn/src/cblas/ATL_srefrot.c b/scikits/learn/src/cblas/ATL_srefrot.c new file mode 100644 index 0000000000000000000000000000000000000000..e1a1d948fff5f1bc34f242e8eb98883eea9567ea --- /dev/null +++ b/scikits/learn/src/cblas/ATL_srefrot.c @@ -0,0 +1,170 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflevel1.h" + +void ATL_srefrot +( + const int N, + float * X, + const int INCX, + float * Y, + const int INCY, + const float C, + const float S +) +{ +/* + * Purpose + * ======= + * + * ATL_srefrot applies a plane rotation to the two n-vectors x and y. + * This routine computes: + * + * [ x_i ] [ c s ] [ x_i ] + * [ y_i ] = [ -s c ] [ y_i ] for all i = 1 .. n. + * + * If n <= 0 or if c = 1 and s = 0, this subroutine returns immediately. + * + * Arguments + * ========= + * + * N (input) const int + * On entry, N specifies the length of the vector x. N must be + * at least zero. Unchanged on exit. + * + * X (input/output) float * + * On entry, X points to the first entry to be accessed of an + * incremented array of size equal to or greater than + * ( 1 + ( n - 1 ) * abs( INCX ) ) * sizeof( float ), + * that contains the vector x. On exit, the entries of the in- + * cremented array X are rotated with the entries of the incre- + * mented array Y. + * + * INCX (input) const int + * On entry, INCX specifies the increment for the elements of X. + * INCX must not be zero. Unchanged on exit. + * + * Y (input/output) float * + * On entry, Y points to the first entry to be accessed of an + * incremented array of size equal to or greater than + * ( 1 + ( n - 1 ) * abs( INCY ) ) * sizeof( float ), + * that contains the vector y. On exit, the entries of the in- + * cremented array Y are rotated with the entries of the incre- + * mented array X. + * + * INCY (input) const int + * On entry, INCY specifies the increment for the elements of Y. + * INCY must not be zero. Unchanged on exit. + * + * C (input) const float + * On entry, C specifies the scalar c definiting the plane rota- + * tion. Unchanged on exit. + * + * S (input) const float + * On entry, S specifies the scalar s definiting the plane rota- + * tion. Unchanged on exit. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float x0, x1, x2, x3, y0, y1, y2, y3; + register const float co = C, si = S; + float * StX; + register int i; + int nu; + const int incX2 = 2 * INCX, incY2 = 2 * INCY, + incX3 = 3 * INCX, incY3 = 3 * INCY, + incX4 = 4 * INCX, incY4 = 4 * INCY; +/* .. + * .. Executable Statements .. + * + */ + if( ( N > 0 ) && !Msone( co, si ) ) + { + if( ( nu = ( N >> 2 ) << 2 ) != 0 ) + { + StX = (float *)X + nu * INCX; + + do + { + x0 = (*X); y0 = (*Y); + x1 = X[INCX ]; y1 = Y[INCY ]; + x2 = X[incX2]; y2 = Y[incY2]; + x3 = X[incX3]; y3 = Y[incY3]; + + *X = co * x0 + si * y0; *Y = co * y0 - si * x0; + X[INCX ] = co * x1 + si * y1; Y[INCY ] = co * y1 - si * x1; + X[incX2] = co * x2 + si * y2; Y[incY2] = co * y2 - si * x2; + X[incX3] = co * x3 + si * y3; Y[incY3] = co * y3 - si * x3; + + X += incX4; + Y += incY4; + + } while( X != StX ); + } + + for( i = N - nu; i != 0; i-- ) + { + x0 = (*X); + y0 = (*Y); + + *X = co * x0 + si * y0; + *Y = co * y0 - si * x0; + + X += INCX; + Y += INCY; + } + } +/* + * End of ATL_srefrot + */ +} diff --git a/scikits/learn/src/cblas/ATL_srefrotg.c b/scikits/learn/src/cblas/ATL_srefrotg.c new file mode 100644 index 0000000000000000000000000000000000000000..4692a03754f5425f057ecb4f0bc8591caef24d4d --- /dev/null +++ b/scikits/learn/src/cblas/ATL_srefrotg.c @@ -0,0 +1,146 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflevel1.h" + +void ATL_srefrotg +( + float * A, + float * B, + float * C, + float * S +) +{ +/* + * Purpose + * ======= + * + * ATL_srefrotg constructs a Givens plane rotation. Given the scalars a + * and b, this routine computes the following quantities: + * + * sigma = sgn(a) if |a| > |b|, sgn(b) otherwise; + * r = sigma * sqrt( a^2 + b^2 ); + * c = a / r if r <> 0, 1 otherwise; + * s = b / r if r <> 0, 0 otherwise. + * + * The numbers c, s and r then satisfy the matrix equation: + * + * [ c s ] [ a ] [ r ] + * [ -s c ] [ b ] = [ 0 ]. + * + * The introduction of sigma is not essential to the computation of a + * Givens rotation matrix, but it permits later stable reconstruction of + * c and s from just one stored number. For this purpose, this routine + * also computes + * + * s if |a| > |b|, + * z = 1 / c if |b| >= |a| and c <> 0, + * 1 if c = 0. + * + * This subroutine returns r overwriting a, and z overwriting b, as well + * as returning c and s. If one later wishes to reconstruct c and s from + * z, it can be done as follows: + * + * if z = 1, set c = 0 and s = 1, + * if |z| < 1, set c = sqrt(1 - z^2) and s = z, + * if |z| > 1, set c = 1 / z and s = sqrt(1 - c^2). + * + * See ``Basic Linear Algebra Subprograms for Fortran Usage'' by C. Law- + * son, R. Hanson, D. Kincaid and F. Krogh, ACM Transactions on Mathema- + * tical Software, 1979, 5(3) pp 308-323, for further information. + * + * Arguments + * ========= + * + * A (input/output) float * + * On entry, A specifies the scalar a. On exit, A is overwritten + * by the scalar r defined above. + * + * B (input/output) float * + * On entry, B specifies the scalar b. On exit, B is overwritten + * by the scalar z defined above. + * + * C (output) float * + * On exit, C specifies the scalar c defined above. + * + * S (output) float * + * On exit, S specifies the scalar s defined above. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float absa, absb, roe, scale, r, tmpa, tmpb, z; +/* .. + * .. Executable Statements .. + * + */ + absa = Msabs( *A ); absb = Msabs( *B ); + roe = ( absa > absb ? (*A) : (*B) ); scale = absa + absb; + + if( scale != ATL_sZERO ) + { + tmpa = (*A) / scale; tmpb = (*B) / scale; + if( roe < ATL_sZERO ) + { r = - scale * sqrt( ( tmpa * tmpa ) + ( tmpb * tmpb ) ); } + else + { r = scale * sqrt( ( tmpa * tmpa ) + ( tmpb * tmpb ) ); } + *C = (*A) / r; *S = (*B) / r; z = ATL_sONE; + if( absa > absb ) { z = *S; } + if( ( absb >= absa ) && ( (*C) != ATL_sZERO ) ) + { z = ATL_sONE / (*C); } + } + else { *C = ATL_sONE; *S = r = z = ATL_sZERO; } + + *A = r; *B = z; +/* + * End of ATL_srefrotg + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsv.c b/scikits/learn/src/cblas/ATL_sreftrsv.c new file mode 100644 index 0000000000000000000000000000000000000000..f77d179215a7c03d0e4aa5c00329bea0c83d1ecd --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsv.c @@ -0,0 +1,213 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsv +( + const enum ATLAS_UPLO UPLO, + const enum ATLAS_TRANS TRANS, + const enum ATLAS_DIAG DIAG, + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsv solves one of the systems of equations + * + * A * x = b, or A' * x = b, + * + * where b and x are n-element vectors and A is an n by n unit, or non- + * unit, upper or lower triangular matrix. + * + * No test for singularity or near-singularity is included in this + * routine. Such tests must be performed before calling this routine. + * + * Arguments + * ========= + * + * UPLO (input) const enum ATLAS_UPLO + * On entry, UPLO specifies whether the matrix is an upper or + * lower triangular matrix as follows: + * + * UPLO = AtlasUpper A is an upper triangular matrix. + * + * UPLO = AtlasLower A is a lower triangular matrix. + * + * Unchanged on exit. + * + * TRANS (input) const enum ATLAS_TRANS + * On entry, TRANS specifies the equations to be solved as fol- + * lows: + * + * TRANS = AtlasNoTrans A * x = b, + * + * TRANS = AtlasConj A * x = b, + * + * TRANS = AtlasTrans A' * x = b, + * + * TRANS = AtlasTrans A' * x = b. + * + * Unchanged on exit. + * + * DIAG (input) const enum ATLAS_DIAG + * On entry, DIAG specifies whether or not A is unit triangu- + * lar as follows: + * + * DIAG = AtlasUnit A is assumed to be unit triangular, + * + * DIAG = AtlasNonUnit A is not assumed to be unit trian- + * gular. + * + * Unchanged on exit. + * + * N (input) const int + * On entry, N specifies the order of the matrix A. N must be at + * least zero. Unchanged on exit. + * + * A (input) const float * + * On entry, A points to an array of size equal to or greater + * than LDA * n * sizeof( float ). Before entry with + * UPLO = AtlasUpper, the leading n by n upper triangular part + * of the array A must contain the upper triangular matrix and + * the strictly lower triangular part of A is not referenced. + * Before entry with UPLO = AtlasLower, the leading n by n lower + * triangular part of the array A must contain the lower trian- + * gular matrix and the strictly upper triangular part of A is + * not referenced. Unchanged on exit. + * + * Note that when DIAG = AtlasUnit, the diagonal elements of A + * are not referenced either, but are assumed to be unity. + * + * LDA (input) const int + * On entry, LDA specifies the leading dimension of A as decla- + * red in the calling (sub) program. LDA must be at least + * MAX( 1, n ). Unchanged on exit. + * + * X (input/output) float * + * On entry, X points to the first entry to be accessed of an + * incremented array of size equal to or greater than + * ( 1 + ( n - 1 ) * abs( INCX ) ) * sizeof( float ), + * that contains the vector x. Before entry, the incremented ar- + * ray X must contain the n element right-hand side vector b. On + * exit, X is overwritten with the solution vector x. + * + * INCX (input) const int + * On entry, INCX specifies the increment for the elements of X. + * INCX must not be zero. Unchanged on exit. + * + * --------------------------------------------------------------------- + */ +/* .. + * .. Executable Statements .. + * + */ + if( N == 0 ) return; + + if( UPLO == AtlasUpper ) + { + if( ( TRANS == AtlasNoTrans ) || ( TRANS == AtlasConj ) ) + { + if( DIAG == AtlasNonUnit ) + { + ATL_sreftrsvUNN( N, A, LDA, X, INCX ); + } + else + { + ATL_sreftrsvUNU( N, A, LDA, X, INCX ); + } + } + else + { + if( DIAG == AtlasNonUnit ) + { + ATL_sreftrsvUTN( N, A, LDA, X, INCX ); + } + else + { + ATL_sreftrsvUTU( N, A, LDA, X, INCX ); + } + } + } + else + { + if( ( TRANS == AtlasNoTrans ) || ( TRANS == AtlasConj ) ) + { + if( DIAG == AtlasNonUnit ) + { + ATL_sreftrsvLNN( N, A, LDA, X, INCX ); + } + else + { + ATL_sreftrsvLNU( N, A, LDA, X, INCX ); + } + } + else + { + if( DIAG == AtlasNonUnit ) + { + ATL_sreftrsvLTN( N, A, LDA, X, INCX ); + } + else + { + ATL_sreftrsvLTU( N, A, LDA, X, INCX ); + } + } + } +/* + * End of ATL_sreftrsv + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvLNN.c b/scikits/learn/src/cblas/ATL_sreftrsvLNN.c new file mode 100644 index 0000000000000000000000000000000000000000..109a68884845bc635115e506d56e7cbc7912e8bf --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvLNN.c @@ -0,0 +1,94 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvLNN +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvLNN( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasLower, AtlasNoTrans, AtlasNonUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx, ldap1 = LDA + 1; +/* .. + * .. Executable Statements .. + * + */ + for( j = 0, jaj = 0, jx = 0; j < N; j++, jaj += ldap1, jx += INCX ) + { + X[jx] /= A[jaj]; t0 = X[jx]; + for( i = j+1, iaij = jaj+1, ix = jx + INCX; + i < N; i++, iaij += 1, ix += INCX ) { X[ix] -= t0 * A[iaij]; } + } +/* + * End of ATL_sreftrsvLNN + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvLNU.c b/scikits/learn/src/cblas/ATL_sreftrsvLNU.c new file mode 100644 index 0000000000000000000000000000000000000000..f97e51be13966fdb163960297b5cc527656fa92f --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvLNU.c @@ -0,0 +1,94 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvLNU +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvLNU( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasLower, AtlasNoTrans, AtlasUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx, ldap1 = LDA + 1; +/* .. + * .. Executable Statements .. + * + */ + for( j = 0, jaj = 0, jx = 0; j < N; j++, jaj += ldap1, jx += INCX ) + { + t0 = X[jx]; + for( i = j+1, iaij = jaj+1, ix = jx + INCX; + i < N; i++, iaij += 1, ix += INCX ) { X[ix] -= t0 * A[iaij]; } + } +/* + * End of ATL_sreftrsvLNU + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvLTN.c b/scikits/learn/src/cblas/ATL_sreftrsvLTN.c new file mode 100644 index 0000000000000000000000000000000000000000..923f84e709f2e37d6bd309b9da3c9a61bda57797 --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvLTN.c @@ -0,0 +1,96 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvLTN +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvLTN( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasLower, AtlasTrans, AtlasNonUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx, ldap1 = LDA + 1; +/* .. + * .. Executable Statements .. + * + */ + for( j = N-1, jaj = (N-1)*(ldap1), jx = (N-1)*INCX; + j >= 0; j--, jaj -= ldap1, jx -= INCX ) + { + t0 = X[jx]; + for( i = j+1, iaij = 1+jaj, ix = jx + INCX; + i < N; i++, iaij += 1, ix += INCX ) { t0 -= A[iaij] * X[ix]; } + t0 /= A[jaj]; X[jx] = t0; + } +/* + * End of ATL_sreftrsvLTN + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvLTU.c b/scikits/learn/src/cblas/ATL_sreftrsvLTU.c new file mode 100644 index 0000000000000000000000000000000000000000..40288f92aa8cfb8f41b8e51a44287086b1daa35d --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvLTU.c @@ -0,0 +1,96 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvLTU +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvLTU( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasLower, AtlasTrans, AtlasUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx, ldap1 = LDA + 1; +/* .. + * .. Executable Statements .. + * + */ + for( j = N-1, jaj = (N-1)*(ldap1), jx = (N-1)*INCX; + j >= 0; j--, jaj -= ldap1, jx -= INCX ) + { + t0 = X[jx]; + for( i = j+1, iaij = 1+jaj, ix = jx + INCX; + i < N; i++, iaij += 1, ix += INCX ) { t0 -= A[iaij] * X[ix]; } + X[jx] = t0; + } +/* + * End of ATL_sreftrsvLTU + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvUNN.c b/scikits/learn/src/cblas/ATL_sreftrsvUNN.c new file mode 100644 index 0000000000000000000000000000000000000000..6e4b465d095a61e6b60482f06e481f68265c879c --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvUNN.c @@ -0,0 +1,95 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvUNN +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvUNN( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasUpper, AtlasNoTrans, AtlasNonUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx; +/* .. + * .. Executable Statements .. + * + */ + for( j = N-1, jaj = (N-1)*LDA, jx = (N-1)*INCX; + j >= 0; j--, jaj -= LDA, jx -= INCX ) + { + X[jx] /= A[j+jaj]; t0 = X[jx]; + for( i = 0, iaij = jaj, ix = 0; i < j; i++, iaij += 1, ix += INCX ) + { X[ix] -= t0 * A[iaij]; } + } +/* + * End of ATL_sreftrsvUNN + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvUNU.c b/scikits/learn/src/cblas/ATL_sreftrsvUNU.c new file mode 100644 index 0000000000000000000000000000000000000000..a3c2acaf764612b2d84f9630de822af2e96d3b1c --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvUNU.c @@ -0,0 +1,95 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvUNU +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvUNU( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasUpper, AtlasNoTrans, AtlasUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx; +/* .. + * .. Executable Statements .. + * + */ + for( j = N-1, jaj = (N-1)*LDA, jx = (N-1)*INCX; + j >= 0; j--, jaj -= LDA, jx -= INCX ) + { + t0 = X[jx]; + for( i = 0, iaij = jaj, ix = 0; i < j; i++, iaij += 1, ix += INCX ) + { X[ix] -= t0 * A[iaij]; } + } +/* + * End of ATL_sreftrsvUNU + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvUTN.c b/scikits/learn/src/cblas/ATL_sreftrsvUTN.c new file mode 100644 index 0000000000000000000000000000000000000000..6a689fa0ebcc2d4ae689294ba9092e8b26628650 --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvUTN.c @@ -0,0 +1,95 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvUTN +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvUTN( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasUpper, AtlasTrans, AtlasNonUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx; +/* .. + * .. Executable Statements .. + * + */ + for( j = 0, jaj = 0,jx = 0; j < N; j++, jaj += LDA, jx += INCX ) + { + t0 = X[jx]; + for( i = 0, iaij = jaj, ix = 0; i < j; i++, iaij += 1, ix += INCX ) + { t0 -= A[iaij] * X[ix]; } + t0 /= A[iaij]; X[jx] = t0; + } +/* + * End of ATL_sreftrsvUTN + */ +} diff --git a/scikits/learn/src/cblas/ATL_sreftrsvUTU.c b/scikits/learn/src/cblas/ATL_sreftrsvUTU.c new file mode 100644 index 0000000000000000000000000000000000000000..244700aa6e9c4621359ea92ec31d2c118716dd54 --- /dev/null +++ b/scikits/learn/src/cblas/ATL_sreftrsvUTU.c @@ -0,0 +1,95 @@ +/* --------------------------------------------------------------------- + * + * -- Automatically Tuned Linear Algebra Software (ATLAS) + * (C) Copyright 2000 All Rights Reserved + * + * -- ATLAS routine -- Version 3.9.24 -- December 25, 2000 + * + * Author : Antoine P. Petitet + * Originally developed at the University of Tennessee, + * Innovative Computing Laboratory, Knoxville TN, 37996-1301, USA. + * + * --------------------------------------------------------------------- + * + * -- Copyright notice and Licensing terms: + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in + * the documentation and/or other materials provided with the distri- + * bution. + * 3. The name of the University, the ATLAS group, or the names of its + * contributors may not be used to endorse or promote products deri- + * ved from this software without specific written permission. + * + * -- Disclaimer: + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT + * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR + * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY + * OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPE- + * CIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED + * TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, + * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEO- + * RY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (IN- + * CLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + * + * --------------------------------------------------------------------- + */ +/* + * Include files + */ +#include "atlas_refmisc.h" +#include "atlas_reflvl2.h" +#include "atlas_reflevel2.h" + +void ATL_sreftrsvUTU +( + const int N, + const float * A, + const int LDA, + float * X, + const int INCX +) +{ +/* + * Purpose + * ======= + * + * ATL_sreftrsvUTU( ... ) + * + * <=> + * + * ATL_sreftrsv( AtlasUpper, AtlasTrans, AtlasUnit, ... ) + * + * See ATL_sreftrsv for details. + * + * --------------------------------------------------------------------- + */ +/* + * .. Local Variables .. + */ + register float t0; + int i, iaij, ix, j, jaj, jx; +/* .. + * .. Executable Statements .. + * + */ + for( j = 0, jaj = 0, jx = 0; j < N; j++, jaj += LDA, jx += INCX ) + { + t0 = X[jx]; + for( i = 0, iaij = jaj, ix = 0; i < j; i++, iaij += 1, ix += INCX ) + { t0 -= A[iaij] * X[ix]; } + X[jx] = t0; + } +/* + * End of ATL_sreftrsvUTU + */ +} diff --git a/scikits/learn/src/cblas/atlas_ssysinfo.h b/scikits/learn/src/cblas/atlas_ssysinfo.h new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/scikits/learn/src/cblas/cblas_scopy.c b/scikits/learn/src/cblas/cblas_scopy.c new file mode 100644 index 0000000000000000000000000000000000000000..38feca50a67db4a78e5c4c6a875b8e1f38f75b06 --- /dev/null +++ b/scikits/learn/src/cblas/cblas_scopy.c @@ -0,0 +1,52 @@ +/* + * Automatically Tuned Linear Algebra Software v3.9.25 + * (C) Copyright 1999 R. Clint Whaley + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in the + * documentation and/or other materials provided with the distribution. + * 3. The name of the ATLAS group or the names of its contributers may + * not be used to endorse or promote products derived from this + * software without specific written permission. + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED + * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR + * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE ATLAS GROUP OR ITS CONTRIBUTORS + * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR + * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF + * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS + * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN + * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) + * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE + * POSSIBILITY OF SUCH DAMAGE. + * + */ + +#define SREAL +#include "atlas_misc.h" +#ifdef ATL_USEPTHREADS + #include "atlas_ptalias1.h" +#endif +#include "atlas_level1.h" +#include "cblas.h" + +void cblas_scopy(const int N, const float *X, const int incX, + float *Y, const int incY) +{ + if (N > 0) + { + if (incX < 0) + { + if (incY < 0) ATL_scopy(N, X, -incX, Y, -incY); + else ATL_scopy(N, X+(1-N)*incX, incX, Y, incY); + } + else if (incY < 0) ATL_scopy(N, X+(N-1)*incX, -incX, Y, -incY); + else ATL_scopy(N, X, incX, Y, incY); + } +} diff --git a/scikits/learn/src/cblas/cblas_srot.c b/scikits/learn/src/cblas/cblas_srot.c new file mode 100644 index 0000000000000000000000000000000000000000..b2365d3866605657aab44ffea3c6e88981ccb660 --- /dev/null +++ b/scikits/learn/src/cblas/cblas_srot.c @@ -0,0 +1,60 @@ +/* + * Automatically Tuned Linear Algebra Software v3.9.25 + * (C) Copyright 1999 R. Clint Whaley + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in the + * documentation and/or other materials provided with the distribution. + * 3. The name of the ATLAS group or the names of its contributers may + * not be used to endorse or promote products derived from this + * software without specific written permission. + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED + * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR + * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE ATLAS GROUP OR ITS CONTRIBUTORS + * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR + * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF + * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS + * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN + * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) + * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE + * POSSIBILITY OF SUCH DAMAGE. + * + */ + +#define SREAL +#include "atlas_misc.h" +#ifdef ATL_USEPTHREADS + #include "atlas_ptalias1.h" +#endif +#include "atlas_level1.h" +#include "cblas.h" + +void cblas_srot(const int N, float *X, const int incX, + float *Y, const int incY, const float c, const float s) +{ + float *x = X, *y = Y; + int incx = incX, incy = incY; + + if (N > 0) + { + if (incX < 0) + { + if (incY < 0) { incx = -incx; incy = -incy; } + else x += -incX * ((N-1)); + } + else if (incY < 0) + { + incy = -incy; + incx = -incx; + x += (N-1)*(incX); + } + ATL_srot(N, x, incx, y, incy, c, s); + } +} diff --git a/scikits/learn/src/cblas/cblas_srotg.c b/scikits/learn/src/cblas/cblas_srotg.c new file mode 100644 index 0000000000000000000000000000000000000000..56bb01757c66afb79f2df9a6dc37d4d887ab4d9f --- /dev/null +++ b/scikits/learn/src/cblas/cblas_srotg.c @@ -0,0 +1,42 @@ +/* + * Automatically Tuned Linear Algebra Software v3.9.25 + * (C) Copyright 1999 R. Clint Whaley + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in the + * documentation and/or other materials provided with the distribution. + * 3. The name of the ATLAS group or the names of its contributers may + * not be used to endorse or promote products derived from this + * software without specific written permission. + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED + * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR + * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE ATLAS GROUP OR ITS CONTRIBUTORS + * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR + * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF + * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS + * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN + * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) + * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE + * POSSIBILITY OF SUCH DAMAGE. + * + */ + +#define SREAL +#include "atlas_misc.h" +#ifdef ATL_USEPTHREADS + #include "atlas_ptalias1.h" +#endif +#include "atlas_level1.h" +#include "cblas.h" + +void cblas_srotg(float *a, float *b, float *c, float *s) +{ + Mjoin(PATL,rotg)(a, b, c, s); +} diff --git a/scikits/learn/src/cblas/cblas_strsv.c b/scikits/learn/src/cblas/cblas_strsv.c new file mode 100644 index 0000000000000000000000000000000000000000..367ffded02f14ae36821e1bf831cdf319a586f3a --- /dev/null +++ b/scikits/learn/src/cblas/cblas_strsv.c @@ -0,0 +1,87 @@ +/* + * Automatically Tuned Linear Algebra Software v3.9.25 + * (C) Copyright 1999 R. Clint Whaley + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions, and the following disclaimer in the + * documentation and/or other materials provided with the distribution. + * 3. The name of the ATLAS group or the names of its contributers may + * not be used to endorse or promote products derived from this + * software without specific written permission. + * + * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS + * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED + * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR + * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE ATLAS GROUP OR ITS CONTRIBUTORS + * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR + * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF + * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS + * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN + * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) + * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE + * POSSIBILITY OF SUCH DAMAGE. + * + */ + +#define SREAL +#include "atlas_misc.h" +#include "cblas.h" +#ifdef ATL_USEPTHREADS + #include "atlas_ptalias2.h" +#endif +#include "atlas_level2.h" + +void cblas_strsv(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo, + const enum CBLAS_TRANSPOSE TA, const enum CBLAS_DIAG Diag, + const int N, const float *A, const int lda, + float *X, const int incX) +{ + int info = 2000; + enum CBLAS_UPLO uplo; + enum CBLAS_TRANSPOSE ta; + #define x X + +#ifndef NoCblasErrorChecks + if (Order != CblasColMajor && Order != CblasRowMajor) + info = cblas_errprn(1, info, "Order must be %d or %d, but is set to %d", + CblasRowMajor, CblasColMajor, Order); + if (Uplo != CblasUpper && Uplo != CblasLower) + info = cblas_errprn(2, info, "UPLO must be %d or %d, but is set to %d", + CblasUpper, CblasLower, Uplo); + if (TA != CblasNoTrans && TA != CblasTrans && TA != CblasConjTrans) + info = cblas_errprn(3, info, + "TransA must be %d, %d or %d, but is set to %d", + CblasNoTrans, CblasTrans, CblasConjTrans, TA); + if (Diag != CblasUnit && Diag != CblasNonUnit) + info = cblas_errprn(4, info, "DIAG must be %d or %d, but is set to %d", + CblasUnit, CblasNonUnit, Diag); + + if (N < 0) info = cblas_errprn(5, info, + "N cannot be less than zero; is set to %d.", N); + if (lda < N || lda < 1) + info = cblas_errprn(7, info, "lda must be >= MAX(N,1): lda=%d N=%d", + lda, N); + if (!incX) info = cblas_errprn(9, info, + "incX cannot be zero; is set to %d.", incX); + if (info != 2000) + { + cblas_xerbla(info, "cblas_strsv", ""); + return; + } +#endif + if (incX < 0) x += (1-N)*incX; + if (Order == CblasColMajor) + ATL_strsv(Uplo, TA, Diag, N, A, lda, x, incX); + else + { + uplo = ( (Uplo == CblasUpper) ? CblasLower : CblasUpper ); + if (TA == CblasNoTrans) ta = CblasTrans; + else ta = CblasNoTrans; + ATL_strsv(uplo, ta, Diag, N, A, lda, x, incX); + } +} diff --git a/scikits/learn/utils/arrayfuncs.c b/scikits/learn/utils/arrayfuncs.c index 6e65702eda574223aae81bb454e19844f11c518a..068f615421efecd38cd4881cb1a93f975e87ee0e 100644 --- a/scikits/learn/utils/arrayfuncs.c +++ b/scikits/learn/utils/arrayfuncs.c @@ -1,4 +1,4 @@ -/* Generated by Cython 0.12.1 on Fri Nov 5 15:03:01 2010 */ +/* Generated by Cython 0.12.1 on Mon Nov 22 14:49:19 2010 */ #define PY_SSIZE_T_CLEAN #include "Python.h" @@ -434,34 +434,6 @@ static void __Pyx_RaiseArgtupleInvalid(const char* func_name, int exact, static int __Pyx_ParseOptionalKeywords(PyObject *kwds, PyObject **argnames[], PyObject *kwds2, PyObject *values[], Py_ssize_t num_pos_args, const char* function_name); /*proto*/ -/* Run-time type information about structs used with buffers */ -struct __Pyx_StructField_; - -typedef struct { - const char* name; /* for error messages only */ - struct __Pyx_StructField_* fields; - size_t size; /* sizeof(type) */ - char typegroup; /* _R_eal, _C_omplex, Signed _I_nt, _U_nsigned int, _S_truct, _P_ointer, _O_bject */ -} __Pyx_TypeInfo; - -typedef struct __Pyx_StructField_ { - __Pyx_TypeInfo* type; - const char* name; - size_t offset; -} __Pyx_StructField; - -typedef struct { - __Pyx_StructField* field; - size_t parent_offset; -} __Pyx_BufFmt_StackElem; - - -static CYTHON_INLINE void __Pyx_SafeReleaseBuffer(Py_buffer* info); -static int __Pyx_GetBufferAndValidate(Py_buffer* buf, PyObject* obj, __Pyx_TypeInfo* dtype, int flags, int nd, int cast, __Pyx_BufFmt_StackElem* stack); - -static CYTHON_INLINE void __Pyx_ErrRestore(PyObject *type, PyObject *value, PyObject *tb); /*proto*/ -static CYTHON_INLINE void __Pyx_ErrFetch(PyObject **type, PyObject **value, PyObject **tb); /*proto*/ - static CYTHON_INLINE void __Pyx_RaiseNeedMoreValuesError(Py_ssize_t index); static CYTHON_INLINE void __Pyx_RaiseTooManyValuesError(void); @@ -477,21 +449,14 @@ static void __Pyx_UnpackTupleError(PyObject *, Py_ssize_t index); /*proto*/ static int __Pyx_ArgTypeTest(PyObject *obj, PyTypeObject *type, int none_allowed, const char *name, int exact); /*proto*/ -#if PY_MAJOR_VERSION < 3 -static int __Pyx_GetBuffer(PyObject *obj, Py_buffer *view, int flags); -static void __Pyx_ReleaseBuffer(Py_buffer *view); -#else -#define __Pyx_GetBuffer PyObject_GetBuffer -#define __Pyx_ReleaseBuffer PyBuffer_Release -#endif - -Py_ssize_t __Pyx_zeros[] = {0, 0}; -Py_ssize_t __Pyx_minusones[] = {-1, -1}; static PyObject *__Pyx_Import(PyObject *name, PyObject *from_list); /*proto*/ static PyObject *__Pyx_GetName(PyObject *dict, PyObject *name); /*proto*/ +static CYTHON_INLINE void __Pyx_ErrRestore(PyObject *type, PyObject *value, PyObject *tb); /*proto*/ +static CYTHON_INLINE void __Pyx_ErrFetch(PyObject **type, PyObject **value, PyObject **tb); /*proto*/ + static void __Pyx_Raise(PyObject *type, PyObject *value, PyObject *tb); /*proto*/ #if CYTHON_CCOMPLEX @@ -643,7 +608,6 @@ static CYTHON_INLINE PyObject *__pyx_f_5numpy_get_array_base(PyArrayObject *); / static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float *, Py_ssize_t); /*proto*/ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(double *, Py_ssize_t); /*proto*/ -static __Pyx_TypeInfo __Pyx_TypeInfo_nn___pyx_t_7scikits_5learn_5utils_10arrayfuncs_DOUBLE = { "scikits.learn.utils.arrayfuncs.DOUBLE", NULL, sizeof(__pyx_t_7scikits_5learn_5utils_10arrayfuncs_DOUBLE), 'R' }; #define __Pyx_MODULE_NAME "scikits.learn.utils.arrayfuncs" int __pyx_module_is_main_scikits__learn__utils__arrayfuncs = 0; @@ -652,15 +616,16 @@ static PyObject *__pyx_builtin_ValueError; static PyObject *__pyx_builtin_range; static PyObject *__pyx_builtin_RuntimeError; static char __pyx_k_1[] = "Unsupported dtype for array X"; -static char __pyx_k_2[] = "ndarray is not C contiguous"; -static char __pyx_k_3[] = "ndarray is not Fortran contiguous"; -static char __pyx_k_4[] = "Non-native byte order not supported"; -static char __pyx_k_5[] = "unknown dtype code in numpy.pxd (%d)"; -static char __pyx_k_6[] = "Format string allocated too short, see comment in numpy.pxd"; -static char __pyx_k_7[] = "Format string allocated too short."; -static char __pyx_k_8[] = "\nSmall collection of auxiliary functions that operate on arrays\n\n"; -static char __pyx_k_9[] = "min_pos (line 42)"; -static char __pyx_k_10[] = "solve_triangular (line 73)"; +static char __pyx_k_2[] = "Unsupported or inconsistent dtype in arrays X, y"; +static char __pyx_k_3[] = "ndarray is not C contiguous"; +static char __pyx_k_4[] = "ndarray is not Fortran contiguous"; +static char __pyx_k_5[] = "Non-native byte order not supported"; +static char __pyx_k_6[] = "unknown dtype code in numpy.pxd (%d)"; +static char __pyx_k_7[] = "Format string allocated too short, see comment in numpy.pxd"; +static char __pyx_k_8[] = "Format string allocated too short."; +static char __pyx_k_9[] = "\nSmall collection of auxiliary functions that operate on arrays\n\n"; +static char __pyx_k_10[] = "min_pos (line 47)"; +static char __pyx_k_11[] = "solve_triangular (line 78)"; static char __pyx_k__B[] = "B"; static char __pyx_k__H[] = "H"; static char __pyx_k__I[] = "I"; @@ -713,13 +678,14 @@ static char __pyx_k__RuntimeError[] = "RuntimeError"; static char __pyx_k__solve_triangular[] = "solve_triangular"; static PyObject *__pyx_kp_s_1; static PyObject *__pyx_kp_u_10; -static PyObject *__pyx_kp_u_2; +static PyObject *__pyx_kp_u_11; +static PyObject *__pyx_kp_s_2; static PyObject *__pyx_kp_u_3; static PyObject *__pyx_kp_u_4; static PyObject *__pyx_kp_u_5; static PyObject *__pyx_kp_u_6; static PyObject *__pyx_kp_u_7; -static PyObject *__pyx_kp_u_9; +static PyObject *__pyx_kp_u_8; static PyObject *__pyx_n_s__L; static PyObject *__pyx_n_s__RuntimeError; static PyObject *__pyx_n_s__ValueError; @@ -756,7 +722,7 @@ static PyObject *__pyx_n_s__type_num; static PyObject *__pyx_n_s__y; static PyObject *__pyx_int_15; -/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":42 +/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":47 * * * def min_pos(np.ndarray X): # <<<<<<<<<<<<<< @@ -765,7 +731,7 @@ static PyObject *__pyx_int_15; */ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject *__pyx_self, PyObject *__pyx_v_X); /*proto*/ -static char __pyx_doc_7scikits_5learn_5utils_10arrayfuncs_min_pos[] = "\n Find the minimum value of an array over positivie values\n "; +static char __pyx_doc_7scikits_5learn_5utils_10arrayfuncs_min_pos[] = "\n Find the minimum value of an array over positivie values\n\n Returns a huge value if none of the values are positive\n "; static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject *__pyx_self, PyObject *__pyx_v_X) { PyObject *__pyx_r = NULL; PyObject *__pyx_t_1 = NULL; @@ -775,28 +741,28 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject * __Pyx_RefNannySetupContext("min_pos"); __pyx_self = __pyx_self; __Pyx_INCREF((PyObject *)__pyx_v_X); - if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_X), __pyx_ptype_5numpy_ndarray, 1, "X", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 42; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_X), __pyx_ptype_5numpy_ndarray, 1, "X", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 47; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":46 - * Find the minimum value of an array over positivie values + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":53 + * Returns a huge value if none of the values are positive * """ * if X.dtype.name == 'float32': # <<<<<<<<<<<<<< * return _float_min_pos(<float *> X.data, X.size) * elif X.dtype.name == 'float64': */ - __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 46; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 53; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); - __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 46; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 53; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_2); __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; - __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float32), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 46; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float32), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 53; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; - __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 46; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 53; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; if (__pyx_t_3) { - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":47 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":54 * """ * if X.dtype.name == 'float32': * return _float_min_pos(<float *> X.data, X.size) # <<<<<<<<<<<<<< @@ -804,11 +770,11 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject * * return _double_min_pos(<double *> X.data, X.size) */ __Pyx_XDECREF(__pyx_r); - __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__size); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 47; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__size); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 54; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); - __pyx_t_4 = __Pyx_PyIndex_AsSsize_t(__pyx_t_1); if (unlikely((__pyx_t_4 == (Py_ssize_t)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 47; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_4 = __Pyx_PyIndex_AsSsize_t(__pyx_t_1); if (unlikely((__pyx_t_4 == (Py_ssize_t)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 54; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; - __pyx_t_1 = PyFloat_FromDouble(__pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(((float *)((PyArrayObject *)__pyx_v_X)->data), __pyx_t_4)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 47; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyFloat_FromDouble(__pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(((float *)((PyArrayObject *)__pyx_v_X)->data), __pyx_t_4)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 54; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); __pyx_r = __pyx_t_1; __pyx_t_1 = 0; @@ -816,26 +782,26 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject * goto __pyx_L5; } - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":48 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":55 * if X.dtype.name == 'float32': * return _float_min_pos(<float *> X.data, X.size) * elif X.dtype.name == 'float64': # <<<<<<<<<<<<<< * return _double_min_pos(<double *> X.data, X.size) * else: */ - __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 48; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 55; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); - __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 48; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 55; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_2); __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; - __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float64), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 48; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float64), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 55; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; - __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 48; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 55; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; if (__pyx_t_3) { - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":49 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":56 * return _float_min_pos(<float *> X.data, X.size) * elif X.dtype.name == 'float64': * return _double_min_pos(<double *> X.data, X.size) # <<<<<<<<<<<<<< @@ -843,11 +809,11 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject * * raise ValueError('Unsupported dtype for array X') */ __Pyx_XDECREF(__pyx_r); - __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__size); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 49; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyObject_GetAttr(__pyx_v_X, __pyx_n_s__size); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 56; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); - __pyx_t_4 = __Pyx_PyIndex_AsSsize_t(__pyx_t_1); if (unlikely((__pyx_t_4 == (Py_ssize_t)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 49; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_4 = __Pyx_PyIndex_AsSsize_t(__pyx_t_1); if (unlikely((__pyx_t_4 == (Py_ssize_t)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 56; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; - __pyx_t_1 = PyFloat_FromDouble(__pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(((double *)((PyArrayObject *)__pyx_v_X)->data), __pyx_t_4)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 49; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyFloat_FromDouble(__pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(((double *)((PyArrayObject *)__pyx_v_X)->data), __pyx_t_4)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 56; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); __pyx_r = __pyx_t_1; __pyx_t_1 = 0; @@ -856,24 +822,24 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject * } /*else*/ { - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":51 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":58 * return _double_min_pos(<double *> X.data, X.size) * else: * raise ValueError('Unsupported dtype for array X') # <<<<<<<<<<<<<< * * */ - __pyx_t_1 = PyTuple_New(1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 51; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_1 = PyTuple_New(1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 58; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_1); __Pyx_INCREF(((PyObject *)__pyx_kp_s_1)); PyTuple_SET_ITEM(__pyx_t_1, 0, ((PyObject *)__pyx_kp_s_1)); __Pyx_GIVEREF(((PyObject *)__pyx_kp_s_1)); - __pyx_t_2 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_1, NULL); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 51; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_2 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_1, NULL); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 58; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_2); __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; __Pyx_Raise(__pyx_t_2, 0, 0); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; - {__pyx_filename = __pyx_f[0]; __pyx_lineno = 51; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + {__pyx_filename = __pyx_f[0]; __pyx_lineno = 58; __pyx_clineno = __LINE__; goto __pyx_L1_error;} } __pyx_L5:; @@ -891,7 +857,7 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_min_pos(PyObject * return __pyx_r; } -/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":54 +/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":61 * * * cdef float _float_min_pos (float *X, Py_ssize_t size): # <<<<<<<<<<<<<< @@ -910,7 +876,7 @@ static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float * int __pyx_t_5; __Pyx_RefNannySetupContext("_float_min_pos"); - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":56 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":63 * cdef float _float_min_pos (float *X, Py_ssize_t size): * cdef Py_ssize_t i * cdef float min_val = DBL_MAX # <<<<<<<<<<<<<< @@ -919,7 +885,7 @@ static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float * */ __pyx_v_min_val = DBL_MAX; - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":57 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":64 * cdef Py_ssize_t i * cdef float min_val = DBL_MAX * for i in range(size): # <<<<<<<<<<<<<< @@ -930,7 +896,7 @@ static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float * for (__pyx_t_2 = 0; __pyx_t_2 < __pyx_t_1; __pyx_t_2+=1) { __pyx_v_i = __pyx_t_2; - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":58 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":65 * cdef float min_val = DBL_MAX * for i in range(size): * if X[i] > 0. and X[i] < min_val: # <<<<<<<<<<<<<< @@ -946,7 +912,7 @@ static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float * } if (__pyx_t_5) { - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":59 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":66 * for i in range(size): * if X[i] > 0. and X[i] < min_val: * min_val = X[i] # <<<<<<<<<<<<<< @@ -959,7 +925,7 @@ static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float * __pyx_L5:; } - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":60 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":67 * if X[i] > 0. and X[i] < min_val: * min_val = X[i] * return min_val # <<<<<<<<<<<<<< @@ -975,7 +941,7 @@ static float __pyx_f_7scikits_5learn_5utils_10arrayfuncs__float_min_pos(float * return __pyx_r; } -/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":63 +/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":70 * * * cdef double _double_min_pos (double *X, Py_ssize_t size): # <<<<<<<<<<<<<< @@ -994,7 +960,7 @@ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(doubl int __pyx_t_5; __Pyx_RefNannySetupContext("_double_min_pos"); - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":65 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":72 * cdef double _double_min_pos (double *X, Py_ssize_t size): * cdef Py_ssize_t i * cdef np.float64_t min_val = FLT_MAX # <<<<<<<<<<<<<< @@ -1003,7 +969,7 @@ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(doubl */ __pyx_v_min_val = FLT_MAX; - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":66 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":73 * cdef Py_ssize_t i * cdef np.float64_t min_val = FLT_MAX * for i in range(size): # <<<<<<<<<<<<<< @@ -1014,7 +980,7 @@ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(doubl for (__pyx_t_2 = 0; __pyx_t_2 < __pyx_t_1; __pyx_t_2+=1) { __pyx_v_i = __pyx_t_2; - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":67 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":74 * cdef np.float64_t min_val = FLT_MAX * for i in range(size): * if X[i] > 0. and X[i] < min_val: # <<<<<<<<<<<<<< @@ -1030,7 +996,7 @@ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(doubl } if (__pyx_t_5) { - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":68 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":75 * for i in range(size): * if X[i] > 0. and X[i] < min_val: * min_val = X[i] # <<<<<<<<<<<<<< @@ -1043,12 +1009,12 @@ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(doubl __pyx_L5:; } - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":69 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":76 * if X[i] > 0. and X[i] < min_val: * min_val = X[i] * return min_val # <<<<<<<<<<<<<< * - * @cython.boundscheck(False) + * def solve_triangular (np.ndarray X, np.ndarray y): */ __pyx_r = __pyx_v_min_val; goto __pyx_L0; @@ -1059,31 +1025,28 @@ static double __pyx_f_7scikits_5learn_5utils_10arrayfuncs__double_min_pos(doubl return __pyx_r; } -/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":73 - * @cython.boundscheck(False) - * @cython.wraparound(False) - * def solve_triangular ( # <<<<<<<<<<<<<< - * np.ndarray[DOUBLE, ndim=2] X, - * np.ndarray[DOUBLE, ndim=1] y +/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":78 + * return min_val + * + * def solve_triangular (np.ndarray X, np.ndarray y): # <<<<<<<<<<<<<< + * """ + * Solves a triangular system (overwrites y) */ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_solve_triangular(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/ -static char __pyx_doc_7scikits_5learn_5utils_10arrayfuncs_solve_triangular[] = "\n Solves a triangular system (overwrites y)\n\n TODO: option for upper, lower, etc.\n "; +static char __pyx_doc_7scikits_5learn_5utils_10arrayfuncs_solve_triangular[] = "\n Solves a triangular system (overwrites y)\n\n Note: The lapack function to solve triangular systems was added to\n scipy v0.9. Remove this when we stop supporting earlier versions.\n "; static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_solve_triangular(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) { PyArrayObject *__pyx_v_X = 0; PyArrayObject *__pyx_v_y = 0; int __pyx_v_lda; - Py_buffer __pyx_bstruct_y; - Py_ssize_t __pyx_bstride_0_y = 0; - Py_ssize_t __pyx_bshape_0_y = 0; - Py_buffer __pyx_bstruct_X; - Py_ssize_t __pyx_bstride_0_X = 0; - Py_ssize_t __pyx_bstride_1_X = 0; - Py_ssize_t __pyx_bshape_0_X = 0; - Py_ssize_t __pyx_bshape_1_X = 0; PyObject *__pyx_r = NULL; - int __pyx_t_1; - size_t __pyx_t_2; + PyObject *__pyx_t_1 = NULL; + PyObject *__pyx_t_2 = NULL; + int __pyx_t_3; + int __pyx_t_4; + int __pyx_t_5; + int __pyx_t_6; + size_t __pyx_t_7; static PyObject **__pyx_pyargnames[] = {&__pyx_n_s__X,&__pyx_n_s__y,0}; __Pyx_RefNannySetupContext("solve_triangular"); __pyx_self = __pyx_self; @@ -1105,11 +1068,11 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_solve_triangular(P values[1] = PyDict_GetItem(__pyx_kwds, __pyx_n_s__y); if (likely(values[1])) kw_args--; else { - __Pyx_RaiseArgtupleInvalid("solve_triangular", 1, 2, 2, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 73; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + __Pyx_RaiseArgtupleInvalid("solve_triangular", 1, 2, 2, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 78; __pyx_clineno = __LINE__; goto __pyx_L3_error;} } } if (unlikely(kw_args > 0)) { - if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, PyTuple_GET_SIZE(__pyx_args), "solve_triangular") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 73; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, PyTuple_GET_SIZE(__pyx_args), "solve_triangular") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 78; __pyx_clineno = __LINE__; goto __pyx_L3_error;} } __pyx_v_X = ((PyArrayObject *)values[0]); __pyx_v_y = ((PyArrayObject *)values[1]); @@ -1121,78 +1084,179 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_solve_triangular(P } goto __pyx_L4_argument_unpacking_done; __pyx_L5_argtuple_error:; - __Pyx_RaiseArgtupleInvalid("solve_triangular", 1, 2, 2, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 73; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + __Pyx_RaiseArgtupleInvalid("solve_triangular", 1, 2, 2, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 78; __pyx_clineno = __LINE__; goto __pyx_L3_error;} __pyx_L3_error:; __Pyx_AddTraceback("scikits.learn.utils.arrayfuncs.solve_triangular"); return NULL; __pyx_L4_argument_unpacking_done:; - __pyx_bstruct_X.buf = NULL; - __pyx_bstruct_y.buf = NULL; - if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_X), __pyx_ptype_5numpy_ndarray, 1, "X", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 74; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_y), __pyx_ptype_5numpy_ndarray, 1, "y", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 75; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - { - __Pyx_BufFmt_StackElem __pyx_stack[1]; - if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_bstruct_X, (PyObject*)__pyx_v_X, &__Pyx_TypeInfo_nn___pyx_t_7scikits_5learn_5utils_10arrayfuncs_DOUBLE, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 73; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_INCREF((PyObject *)__pyx_v_X); + __Pyx_INCREF((PyObject *)__pyx_v_y); + if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_X), __pyx_ptype_5numpy_ndarray, 1, "X", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 78; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_y), __pyx_ptype_5numpy_ndarray, 1, "y", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 78; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":87 + * cdef int lda + * + * if X.dtype.name == 'float64' and y.dtype.name == 'float64': # <<<<<<<<<<<<<< + * lda = <int> X.strides[0] / sizeof(double) + * + */ + __pyx_t_1 = PyObject_GetAttr(((PyObject *)__pyx_v_X), __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float64), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + if (__pyx_t_3) { + __pyx_t_1 = PyObject_GetAttr(((PyObject *)__pyx_v_y), __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float64), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + __pyx_t_4 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_4 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 87; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_5 = __pyx_t_4; + } else { + __pyx_t_5 = __pyx_t_3; } - __pyx_bstride_0_X = __pyx_bstruct_X.strides[0]; __pyx_bstride_1_X = __pyx_bstruct_X.strides[1]; - __pyx_bshape_0_X = __pyx_bstruct_X.shape[0]; __pyx_bshape_1_X = __pyx_bstruct_X.shape[1]; - { - __Pyx_BufFmt_StackElem __pyx_stack[1]; - if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_bstruct_y, (PyObject*)__pyx_v_y, &__Pyx_TypeInfo_nn___pyx_t_7scikits_5learn_5utils_10arrayfuncs_DOUBLE, PyBUF_FORMAT| PyBUF_STRIDES, 1, 0, __pyx_stack) == -1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 73; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + if (__pyx_t_5) { + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":88 + * + * if X.dtype.name == 'float64' and y.dtype.name == 'float64': + * lda = <int> X.strides[0] / sizeof(double) # <<<<<<<<<<<<<< + * + * cblas_dtrsv (CblasRowMajor, CblasLower, CblasNoTrans, + */ + __pyx_t_6 = ((int)(__pyx_v_X->strides[0])); + __pyx_t_7 = (sizeof(double)); + if (unlikely(__pyx_t_7 == 0)) { + PyErr_Format(PyExc_ZeroDivisionError, "integer division or modulo by zero"); + {__pyx_filename = __pyx_f[0]; __pyx_lineno = 88; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + } + __pyx_v_lda = (__pyx_t_6 / __pyx_t_7); + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":92 + * cblas_dtrsv (CblasRowMajor, CblasLower, CblasNoTrans, + * CblasNonUnit, <int> X.shape[0], <double *> X.data, + * lda, <double *> y.data, 1); # <<<<<<<<<<<<<< + * + * elif X.dtype.name == 'float32' and y.dtype.name == 'float32': + */ + cblas_dtrsv(CblasRowMajor, CblasLower, CblasNoTrans, CblasNonUnit, ((int)(__pyx_v_X->dimensions[0])), ((double *)__pyx_v_X->data), __pyx_v_lda, ((double *)__pyx_v_y->data), 1); + goto __pyx_L6; } - __pyx_bstride_0_y = __pyx_bstruct_y.strides[0]; - __pyx_bshape_0_y = __pyx_bstruct_y.shape[0]; - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":82 - * TODO: option for upper, lower, etc. - * """ - * cdef int lda = <int> X.strides[0] / sizeof(double) # <<<<<<<<<<<<<< + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":94 + * lda, <double *> y.data, 1); + * + * elif X.dtype.name == 'float32' and y.dtype.name == 'float32': # <<<<<<<<<<<<<< + * lda = <int> X.strides[0] / sizeof(float) + * + */ + __pyx_t_1 = PyObject_GetAttr(((PyObject *)__pyx_v_X), __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float32), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + __pyx_t_5 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_5 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + if (__pyx_t_5) { + __pyx_t_1 = PyObject_GetAttr(((PyObject *)__pyx_v_y), __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float32), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_4 = __pyx_t_3; + } else { + __pyx_t_4 = __pyx_t_5; + } + if (__pyx_t_4) { + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":95 + * + * elif X.dtype.name == 'float32' and y.dtype.name == 'float32': + * lda = <int> X.strides[0] / sizeof(float) # <<<<<<<<<<<<<< * - * cblas_dtrsv (CblasRowMajor, CblasLower, CblasNoTrans, + * cblas_strsv (CblasRowMajor, CblasLower, CblasNoTrans, + */ + __pyx_t_6 = ((int)(__pyx_v_X->strides[0])); + __pyx_t_7 = (sizeof(float)); + if (unlikely(__pyx_t_7 == 0)) { + PyErr_Format(PyExc_ZeroDivisionError, "integer division or modulo by zero"); + {__pyx_filename = __pyx_f[0]; __pyx_lineno = 95; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + } + __pyx_v_lda = (__pyx_t_6 / __pyx_t_7); + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":99 + * cblas_strsv (CblasRowMajor, CblasLower, CblasNoTrans, + * CblasNonUnit, <int> X.shape[0], <float *> X.data, + * lda, <float *> y.data, 1); # <<<<<<<<<<<<<< + * else: + * raise ValueError ('Unsupported or inconsistent dtype in arrays X, y') */ - __pyx_t_1 = ((int)(__pyx_v_X->strides[0])); - __pyx_t_2 = (sizeof(double)); - if (unlikely(__pyx_t_2 == 0)) { - PyErr_Format(PyExc_ZeroDivisionError, "integer division or modulo by zero"); - {__pyx_filename = __pyx_f[0]; __pyx_lineno = 82; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + cblas_strsv(CblasRowMajor, CblasLower, CblasNoTrans, CblasNonUnit, ((int)(__pyx_v_X->dimensions[0])), ((float *)__pyx_v_X->data), __pyx_v_lda, ((float *)__pyx_v_y->data), 1); + goto __pyx_L6; } - __pyx_v_lda = (__pyx_t_1 / __pyx_t_2); + /*else*/ { - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":86 - * cblas_dtrsv (CblasRowMajor, CblasLower, CblasNoTrans, - * CblasNonUnit, <int> X.shape[0], <double *> X.data, - * lda, <double *> y.data, 1); # <<<<<<<<<<<<<< + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":101 + * lda, <float *> y.data, 1); + * else: + * raise ValueError ('Unsupported or inconsistent dtype in arrays X, y') # <<<<<<<<<<<<<< * * */ - cblas_dtrsv(CblasRowMajor, CblasLower, CblasNoTrans, CblasNonUnit, ((int)(__pyx_v_X->dimensions[0])), ((double *)__pyx_v_X->data), __pyx_v_lda, ((double *)__pyx_v_y->data), 1); + __pyx_t_1 = PyTuple_New(1); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 101; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_INCREF(((PyObject *)__pyx_kp_s_2)); + PyTuple_SET_ITEM(__pyx_t_1, 0, ((PyObject *)__pyx_kp_s_2)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_s_2)); + __pyx_t_2 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_1, NULL); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 101; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __Pyx_Raise(__pyx_t_2, 0, 0); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + {__pyx_filename = __pyx_f[0]; __pyx_lineno = 101; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + } + __pyx_L6:; __pyx_r = Py_None; __Pyx_INCREF(Py_None); goto __pyx_L0; __pyx_L1_error:; - { PyObject *__pyx_type, *__pyx_value, *__pyx_tb; - __Pyx_ErrFetch(&__pyx_type, &__pyx_value, &__pyx_tb); - __Pyx_SafeReleaseBuffer(&__pyx_bstruct_y); - __Pyx_SafeReleaseBuffer(&__pyx_bstruct_X); - __Pyx_ErrRestore(__pyx_type, __pyx_value, __pyx_tb);} + __Pyx_XDECREF(__pyx_t_1); + __Pyx_XDECREF(__pyx_t_2); __Pyx_AddTraceback("scikits.learn.utils.arrayfuncs.solve_triangular"); __pyx_r = NULL; - goto __pyx_L2; __pyx_L0:; - __Pyx_SafeReleaseBuffer(&__pyx_bstruct_y); - __Pyx_SafeReleaseBuffer(&__pyx_bstruct_X); - __pyx_L2:; + __Pyx_DECREF((PyObject *)__pyx_v_X); + __Pyx_DECREF((PyObject *)__pyx_v_y); __Pyx_XGIVEREF(__pyx_r); __Pyx_RefNannyFinishContext(); return __pyx_r; } -/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":89 +/* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":104 + * * + * def cholesky_delete (np.ndarray L, int go_out): # <<<<<<<<<<<<<< * - * def cholesky_delete ( # <<<<<<<<<<<<<< - * np.ndarray[DOUBLE, ndim=2] L, - * int go_out): + * cdef int n = <int> L.shape[0] */ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_cholesky_delete(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/ @@ -1201,14 +1265,12 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_cholesky_delete(Py int __pyx_v_go_out; int __pyx_v_n; int __pyx_v_m; - Py_buffer __pyx_bstruct_L; - Py_ssize_t __pyx_bstride_0_L = 0; - Py_ssize_t __pyx_bstride_1_L = 0; - Py_ssize_t __pyx_bshape_0_L = 0; - Py_ssize_t __pyx_bshape_1_L = 0; PyObject *__pyx_r = NULL; - int __pyx_t_1; - size_t __pyx_t_2; + PyObject *__pyx_t_1 = NULL; + PyObject *__pyx_t_2 = NULL; + int __pyx_t_3; + int __pyx_t_4; + size_t __pyx_t_5; static PyObject **__pyx_pyargnames[] = {&__pyx_n_s__L,&__pyx_n_s__go_out,0}; __Pyx_RefNannySetupContext("cholesky_delete"); __pyx_self = __pyx_self; @@ -1230,79 +1292,138 @@ static PyObject *__pyx_pf_7scikits_5learn_5utils_10arrayfuncs_cholesky_delete(Py values[1] = PyDict_GetItem(__pyx_kwds, __pyx_n_s__go_out); if (likely(values[1])) kw_args--; else { - __Pyx_RaiseArgtupleInvalid("cholesky_delete", 1, 2, 2, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 89; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + __Pyx_RaiseArgtupleInvalid("cholesky_delete", 1, 2, 2, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 104; __pyx_clineno = __LINE__; goto __pyx_L3_error;} } } if (unlikely(kw_args > 0)) { - if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, PyTuple_GET_SIZE(__pyx_args), "cholesky_delete") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 89; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, PyTuple_GET_SIZE(__pyx_args), "cholesky_delete") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 104; __pyx_clineno = __LINE__; goto __pyx_L3_error;} } __pyx_v_L = ((PyArrayObject *)values[0]); - __pyx_v_go_out = __Pyx_PyInt_AsInt(values[1]); if (unlikely((__pyx_v_go_out == (int)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 91; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + __pyx_v_go_out = __Pyx_PyInt_AsInt(values[1]); if (unlikely((__pyx_v_go_out == (int)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 104; __pyx_clineno = __LINE__; goto __pyx_L3_error;} } else if (PyTuple_GET_SIZE(__pyx_args) != 2) { goto __pyx_L5_argtuple_error; } else { __pyx_v_L = ((PyArrayObject *)PyTuple_GET_ITEM(__pyx_args, 0)); - __pyx_v_go_out = __Pyx_PyInt_AsInt(PyTuple_GET_ITEM(__pyx_args, 1)); if (unlikely((__pyx_v_go_out == (int)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 91; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + __pyx_v_go_out = __Pyx_PyInt_AsInt(PyTuple_GET_ITEM(__pyx_args, 1)); if (unlikely((__pyx_v_go_out == (int)-1) && PyErr_Occurred())) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 104; __pyx_clineno = __LINE__; goto __pyx_L3_error;} } goto __pyx_L4_argument_unpacking_done; __pyx_L5_argtuple_error:; - __Pyx_RaiseArgtupleInvalid("cholesky_delete", 1, 2, 2, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 89; __pyx_clineno = __LINE__; goto __pyx_L3_error;} + __Pyx_RaiseArgtupleInvalid("cholesky_delete", 1, 2, 2, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 104; __pyx_clineno = __LINE__; goto __pyx_L3_error;} __pyx_L3_error:; __Pyx_AddTraceback("scikits.learn.utils.arrayfuncs.cholesky_delete"); return NULL; __pyx_L4_argument_unpacking_done:; - __pyx_bstruct_L.buf = NULL; - if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_L), __pyx_ptype_5numpy_ndarray, 1, "L", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 90; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - { - __Pyx_BufFmt_StackElem __pyx_stack[1]; - if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_bstruct_L, (PyObject*)__pyx_v_L, &__Pyx_TypeInfo_nn___pyx_t_7scikits_5learn_5utils_10arrayfuncs_DOUBLE, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 89; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - } - __pyx_bstride_0_L = __pyx_bstruct_L.strides[0]; __pyx_bstride_1_L = __pyx_bstruct_L.strides[1]; - __pyx_bshape_0_L = __pyx_bstruct_L.shape[0]; __pyx_bshape_1_L = __pyx_bstruct_L.shape[1]; + __Pyx_INCREF((PyObject *)__pyx_v_L); + if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_L), __pyx_ptype_5numpy_ndarray, 1, "L", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 104; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":93 - * int go_out): + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":106 + * def cholesky_delete (np.ndarray L, int go_out): * * cdef int n = <int> L.shape[0] # <<<<<<<<<<<<<< - * cdef int m = <int> L.strides[0] / sizeof (double) - * c_cholesky_delete (m, n, <double *> L.data, go_out) + * cdef int m + * */ __pyx_v_n = ((int)(__pyx_v_L->dimensions[0])); - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":94 + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":109 + * cdef int m * - * cdef int n = <int> L.shape[0] - * cdef int m = <int> L.strides[0] / sizeof (double) # <<<<<<<<<<<<<< - * c_cholesky_delete (m, n, <double *> L.data, go_out) - */ - __pyx_t_1 = ((int)(__pyx_v_L->strides[0])); - __pyx_t_2 = (sizeof(double)); - if (unlikely(__pyx_t_2 == 0)) { - PyErr_Format(PyExc_ZeroDivisionError, "integer division or modulo by zero"); - {__pyx_filename = __pyx_f[0]; __pyx_lineno = 94; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + * if L.dtype.name == 'float64': # <<<<<<<<<<<<<< + * m = <int> L.strides[0] / sizeof (double) + * double_cholesky_delete (m, n, <double *> L.data, go_out) + */ + __pyx_t_1 = PyObject_GetAttr(((PyObject *)__pyx_v_L), __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 109; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 109; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float64), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 109; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 109; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + if (__pyx_t_3) { + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":110 + * + * if L.dtype.name == 'float64': + * m = <int> L.strides[0] / sizeof (double) # <<<<<<<<<<<<<< + * double_cholesky_delete (m, n, <double *> L.data, go_out) + * elif L.dtype.name == 'float32': + */ + __pyx_t_4 = ((int)(__pyx_v_L->strides[0])); + __pyx_t_5 = (sizeof(double)); + if (unlikely(__pyx_t_5 == 0)) { + PyErr_Format(PyExc_ZeroDivisionError, "integer division or modulo by zero"); + {__pyx_filename = __pyx_f[0]; __pyx_lineno = 110; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + } + __pyx_v_m = (__pyx_t_4 / __pyx_t_5); + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":111 + * if L.dtype.name == 'float64': + * m = <int> L.strides[0] / sizeof (double) + * double_cholesky_delete (m, n, <double *> L.data, go_out) # <<<<<<<<<<<<<< + * elif L.dtype.name == 'float32': + * m = <int> L.strides[0] / sizeof (float) + */ + double_cholesky_delete(__pyx_v_m, __pyx_v_n, ((double *)__pyx_v_L->data), __pyx_v_go_out); + goto __pyx_L6; } - __pyx_v_m = (__pyx_t_1 / __pyx_t_2); - /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":95 - * cdef int n = <int> L.shape[0] - * cdef int m = <int> L.strides[0] / sizeof (double) - * c_cholesky_delete (m, n, <double *> L.data, go_out) # <<<<<<<<<<<<<< + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":112 + * m = <int> L.strides[0] / sizeof (double) + * double_cholesky_delete (m, n, <double *> L.data, go_out) + * elif L.dtype.name == 'float32': # <<<<<<<<<<<<<< + * m = <int> L.strides[0] / sizeof (float) + * float_cholesky_delete (m, n, <float *> L.data, go_out) + */ + __pyx_t_1 = PyObject_GetAttr(((PyObject *)__pyx_v_L), __pyx_n_s__dtype); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 112; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __pyx_t_2 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__name); if (unlikely(!__pyx_t_2)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 112; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_2); + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + __pyx_t_1 = PyObject_RichCompare(__pyx_t_2, ((PyObject *)__pyx_n_s__float32), Py_EQ); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 112; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_GOTREF(__pyx_t_1); + __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; + __pyx_t_3 = __Pyx_PyObject_IsTrue(__pyx_t_1); if (unlikely(__pyx_t_3 < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 112; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0; + if (__pyx_t_3) { + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":113 + * double_cholesky_delete (m, n, <double *> L.data, go_out) + * elif L.dtype.name == 'float32': + * m = <int> L.strides[0] / sizeof (float) # <<<<<<<<<<<<<< + * float_cholesky_delete (m, n, <float *> L.data, go_out) + * */ - cholesky_delete(__pyx_v_m, __pyx_v_n, ((double *)__pyx_v_L->data), __pyx_v_go_out); + __pyx_t_4 = ((int)(__pyx_v_L->strides[0])); + __pyx_t_5 = (sizeof(float)); + if (unlikely(__pyx_t_5 == 0)) { + PyErr_Format(PyExc_ZeroDivisionError, "integer division or modulo by zero"); + {__pyx_filename = __pyx_f[0]; __pyx_lineno = 113; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + } + __pyx_v_m = (__pyx_t_4 / __pyx_t_5); + + /* "/home/fabian/dev/scikit-learn/scikits/learn/utils/arrayfuncs.pyx":114 + * elif L.dtype.name == 'float32': + * m = <int> L.strides[0] / sizeof (float) + * float_cholesky_delete (m, n, <float *> L.data, go_out) # <<<<<<<<<<<<<< + * + */ + float_cholesky_delete(__pyx_v_m, __pyx_v_n, ((float *)__pyx_v_L->data), __pyx_v_go_out); + goto __pyx_L6; + } + __pyx_L6:; __pyx_r = Py_None; __Pyx_INCREF(Py_None); goto __pyx_L0; __pyx_L1_error:; - { PyObject *__pyx_type, *__pyx_value, *__pyx_tb; - __Pyx_ErrFetch(&__pyx_type, &__pyx_value, &__pyx_tb); - __Pyx_SafeReleaseBuffer(&__pyx_bstruct_L); - __Pyx_ErrRestore(__pyx_type, __pyx_value, __pyx_tb);} + __Pyx_XDECREF(__pyx_t_1); + __Pyx_XDECREF(__pyx_t_2); __Pyx_AddTraceback("scikits.learn.utils.arrayfuncs.cholesky_delete"); __pyx_r = NULL; - goto __pyx_L2; __pyx_L0:; - __Pyx_SafeReleaseBuffer(&__pyx_bstruct_L); - __pyx_L2:; + __Pyx_DECREF((PyObject *)__pyx_v_L); __Pyx_XGIVEREF(__pyx_r); __Pyx_RefNannyFinishContext(); return __pyx_r; @@ -1437,9 +1558,9 @@ static int __pyx_pf_5numpy_7ndarray___getbuffer__(PyObject *__pyx_v_self, Py_buf */ __pyx_t_4 = PyTuple_New(1); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 205; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_4); - __Pyx_INCREF(((PyObject *)__pyx_kp_u_2)); - PyTuple_SET_ITEM(__pyx_t_4, 0, ((PyObject *)__pyx_kp_u_2)); - __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_2)); + __Pyx_INCREF(((PyObject *)__pyx_kp_u_3)); + PyTuple_SET_ITEM(__pyx_t_4, 0, ((PyObject *)__pyx_kp_u_3)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_3)); __pyx_t_5 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_4, NULL); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 205; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0; @@ -1483,9 +1604,9 @@ static int __pyx_pf_5numpy_7ndarray___getbuffer__(PyObject *__pyx_v_self, Py_buf */ __pyx_t_5 = PyTuple_New(1); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 209; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); - __Pyx_INCREF(((PyObject *)__pyx_kp_u_3)); - PyTuple_SET_ITEM(__pyx_t_5, 0, ((PyObject *)__pyx_kp_u_3)); - __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_3)); + __Pyx_INCREF(((PyObject *)__pyx_kp_u_4)); + PyTuple_SET_ITEM(__pyx_t_5, 0, ((PyObject *)__pyx_kp_u_4)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_4)); __pyx_t_4 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_5, NULL); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 209; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_4); __Pyx_DECREF(__pyx_t_5); __pyx_t_5 = 0; @@ -1760,9 +1881,9 @@ static int __pyx_pf_5numpy_7ndarray___getbuffer__(PyObject *__pyx_v_self, Py_buf */ __pyx_t_4 = PyTuple_New(1); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 247; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_4); - __Pyx_INCREF(((PyObject *)__pyx_kp_u_4)); - PyTuple_SET_ITEM(__pyx_t_4, 0, ((PyObject *)__pyx_kp_u_4)); - __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_4)); + __Pyx_INCREF(((PyObject *)__pyx_kp_u_5)); + PyTuple_SET_ITEM(__pyx_t_4, 0, ((PyObject *)__pyx_kp_u_5)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_5)); __pyx_t_5 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_4, NULL); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 247; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); __Pyx_DECREF(__pyx_t_4); __pyx_t_4 = 0; @@ -2004,7 +2125,7 @@ static int __pyx_pf_5numpy_7ndarray___getbuffer__(PyObject *__pyx_v_self, Py_buf */ __pyx_t_5 = PyInt_FromLong(__pyx_v_t); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 266; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); - __pyx_t_4 = PyNumber_Remainder(((PyObject *)__pyx_kp_u_5), __pyx_t_5); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 266; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_4 = PyNumber_Remainder(((PyObject *)__pyx_kp_u_6), __pyx_t_5); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 266; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_4); __Pyx_DECREF(__pyx_t_5); __pyx_t_5 = 0; __pyx_t_5 = PyTuple_New(1); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 266; __pyx_clineno = __LINE__; goto __pyx_L1_error;} @@ -2513,9 +2634,9 @@ static CYTHON_INLINE char *__pyx_f_5numpy__util_dtypestring(PyArray_Descr *__pyx */ __pyx_t_5 = PyTuple_New(1); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 786; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); - __Pyx_INCREF(((PyObject *)__pyx_kp_u_6)); - PyTuple_SET_ITEM(__pyx_t_5, 0, ((PyObject *)__pyx_kp_u_6)); - __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_6)); + __Pyx_INCREF(((PyObject *)__pyx_kp_u_7)); + PyTuple_SET_ITEM(__pyx_t_5, 0, ((PyObject *)__pyx_kp_u_7)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_7)); __pyx_t_3 = PyObject_Call(__pyx_builtin_RuntimeError, __pyx_t_5, NULL); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 786; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_3); __Pyx_DECREF(__pyx_t_5); __pyx_t_5 = 0; @@ -2570,9 +2691,9 @@ static CYTHON_INLINE char *__pyx_f_5numpy__util_dtypestring(PyArray_Descr *__pyx */ __pyx_t_3 = PyTuple_New(1); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 790; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_3); - __Pyx_INCREF(((PyObject *)__pyx_kp_u_4)); - PyTuple_SET_ITEM(__pyx_t_3, 0, ((PyObject *)__pyx_kp_u_4)); - __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_4)); + __Pyx_INCREF(((PyObject *)__pyx_kp_u_5)); + PyTuple_SET_ITEM(__pyx_t_3, 0, ((PyObject *)__pyx_kp_u_5)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_5)); __pyx_t_5 = PyObject_Call(__pyx_builtin_ValueError, __pyx_t_3, NULL); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 790; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; @@ -2679,9 +2800,9 @@ static CYTHON_INLINE char *__pyx_f_5numpy__util_dtypestring(PyArray_Descr *__pyx */ __pyx_t_3 = PyTuple_New(1); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 810; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_3); - __Pyx_INCREF(((PyObject *)__pyx_kp_u_7)); - PyTuple_SET_ITEM(__pyx_t_3, 0, ((PyObject *)__pyx_kp_u_7)); - __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_7)); + __Pyx_INCREF(((PyObject *)__pyx_kp_u_8)); + PyTuple_SET_ITEM(__pyx_t_3, 0, ((PyObject *)__pyx_kp_u_8)); + __Pyx_GIVEREF(((PyObject *)__pyx_kp_u_8)); __pyx_t_5 = PyObject_Call(__pyx_builtin_RuntimeError, __pyx_t_3, NULL); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 810; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; @@ -3029,7 +3150,7 @@ static CYTHON_INLINE char *__pyx_f_5numpy__util_dtypestring(PyArray_Descr *__pyx * f += 1 * else: */ - __pyx_t_3 = PyNumber_Remainder(((PyObject *)__pyx_kp_u_5), __pyx_v_t); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 831; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_t_3 = PyNumber_Remainder(((PyObject *)__pyx_kp_u_6), __pyx_v_t); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 831; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_3); __pyx_t_5 = PyTuple_New(1); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 831; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_5); @@ -3253,7 +3374,7 @@ static void __pyx_init_filenames(void); /*proto*/ static struct PyModuleDef __pyx_moduledef = { PyModuleDef_HEAD_INIT, __Pyx_NAMESTR("arrayfuncs"), - __Pyx_DOCSTR(__pyx_k_8), /* m_doc */ + __Pyx_DOCSTR(__pyx_k_9), /* m_doc */ -1, /* m_size */ __pyx_methods /* m_methods */, NULL, /* m_reload */ @@ -3266,13 +3387,14 @@ static struct PyModuleDef __pyx_moduledef = { static __Pyx_StringTabEntry __pyx_string_tab[] = { {&__pyx_kp_s_1, __pyx_k_1, sizeof(__pyx_k_1), 0, 0, 1, 0}, {&__pyx_kp_u_10, __pyx_k_10, sizeof(__pyx_k_10), 0, 1, 0, 0}, - {&__pyx_kp_u_2, __pyx_k_2, sizeof(__pyx_k_2), 0, 1, 0, 0}, + {&__pyx_kp_u_11, __pyx_k_11, sizeof(__pyx_k_11), 0, 1, 0, 0}, + {&__pyx_kp_s_2, __pyx_k_2, sizeof(__pyx_k_2), 0, 0, 1, 0}, {&__pyx_kp_u_3, __pyx_k_3, sizeof(__pyx_k_3), 0, 1, 0, 0}, {&__pyx_kp_u_4, __pyx_k_4, sizeof(__pyx_k_4), 0, 1, 0, 0}, {&__pyx_kp_u_5, __pyx_k_5, sizeof(__pyx_k_5), 0, 1, 0, 0}, {&__pyx_kp_u_6, __pyx_k_6, sizeof(__pyx_k_6), 0, 1, 0, 0}, {&__pyx_kp_u_7, __pyx_k_7, sizeof(__pyx_k_7), 0, 1, 0, 0}, - {&__pyx_kp_u_9, __pyx_k_9, sizeof(__pyx_k_9), 0, 1, 0, 0}, + {&__pyx_kp_u_8, __pyx_k_8, sizeof(__pyx_k_8), 0, 1, 0, 0}, {&__pyx_n_s__L, __pyx_k__L, sizeof(__pyx_k__L), 0, 0, 1, 1}, {&__pyx_n_s__RuntimeError, __pyx_k__RuntimeError, sizeof(__pyx_k__RuntimeError), 0, 0, 1, 1}, {&__pyx_n_s__ValueError, __pyx_k__ValueError, sizeof(__pyx_k__ValueError), 0, 0, 1, 1}, @@ -3310,8 +3432,8 @@ static __Pyx_StringTabEntry __pyx_string_tab[] = { {0, 0, 0, 0, 0, 0, 0} }; static int __Pyx_InitCachedBuiltins(void) { - __pyx_builtin_ValueError = __Pyx_GetName(__pyx_b, __pyx_n_s__ValueError); if (!__pyx_builtin_ValueError) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 51; __pyx_clineno = __LINE__; goto __pyx_L1_error;} - __pyx_builtin_range = __Pyx_GetName(__pyx_b, __pyx_n_s__range); if (!__pyx_builtin_range) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 57; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_builtin_ValueError = __Pyx_GetName(__pyx_b, __pyx_n_s__ValueError); if (!__pyx_builtin_ValueError) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 58; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + __pyx_builtin_range = __Pyx_GetName(__pyx_b, __pyx_n_s__range); if (!__pyx_builtin_range) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 64; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __pyx_builtin_RuntimeError = __Pyx_GetName(__pyx_b, __pyx_n_s__RuntimeError); if (!__pyx_builtin_RuntimeError) {__pyx_filename = __pyx_f[1]; __pyx_lineno = 786; __pyx_clineno = __LINE__; goto __pyx_L1_error;} return 0; __pyx_L1_error:; @@ -3364,7 +3486,7 @@ PyMODINIT_FUNC PyInit_arrayfuncs(void) #endif /*--- Module creation code ---*/ #if PY_MAJOR_VERSION < 3 - __pyx_m = Py_InitModule4(__Pyx_NAMESTR("arrayfuncs"), __pyx_methods, __Pyx_DOCSTR(__pyx_k_8), 0, PYTHON_API_VERSION); + __pyx_m = Py_InitModule4(__Pyx_NAMESTR("arrayfuncs"), __pyx_methods, __Pyx_DOCSTR(__pyx_k_9), 0, PYTHON_API_VERSION); #else __pyx_m = PyModule_Create(&__pyx_moduledef); #endif @@ -3418,14 +3540,14 @@ PyMODINIT_FUNC PyInit_arrayfuncs(void) __pyx_t_3 = __Pyx_GetAttrString(__pyx_t_2, "__doc__"); __Pyx_GOTREF(__pyx_t_3); __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; - if (PyDict_SetItem(__pyx_t_1, ((PyObject *)__pyx_kp_u_9), __pyx_t_3) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + if (PyDict_SetItem(__pyx_t_1, ((PyObject *)__pyx_kp_u_10), __pyx_t_3) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; __pyx_t_3 = PyObject_GetAttr(__pyx_m, __pyx_n_s__solve_triangular); if (unlikely(!__pyx_t_3)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_GOTREF(__pyx_t_3); __pyx_t_2 = __Pyx_GetAttrString(__pyx_t_3, "__doc__"); __Pyx_GOTREF(__pyx_t_2); __Pyx_DECREF(__pyx_t_3); __pyx_t_3 = 0; - if (PyDict_SetItem(__pyx_t_1, ((PyObject *)__pyx_kp_u_10), __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;} + if (PyDict_SetItem(__pyx_t_1, ((PyObject *)__pyx_kp_u_11), __pyx_t_2) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(__pyx_t_2); __pyx_t_2 = 0; if (PyObject_SetAttr(__pyx_m, __pyx_n_s____test__, ((PyObject *)__pyx_t_1)) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;} __Pyx_DECREF(((PyObject *)__pyx_t_1)); __pyx_t_1 = 0; @@ -3589,477 +3711,6 @@ bad: return -1; } -static CYTHON_INLINE int __Pyx_IsLittleEndian(void) { - unsigned int n = 1; - return *(unsigned char*)(&n) != 0; -} - -typedef struct { - __Pyx_StructField root; - __Pyx_BufFmt_StackElem* head; - size_t fmt_offset; - int new_count, enc_count; - int is_complex; - char enc_type; - char packmode; -} __Pyx_BufFmt_Context; - -static void __Pyx_BufFmt_Init(__Pyx_BufFmt_Context* ctx, - __Pyx_BufFmt_StackElem* stack, - __Pyx_TypeInfo* type) { - stack[0].field = &ctx->root; - stack[0].parent_offset = 0; - ctx->root.type = type; - ctx->root.name = "buffer dtype"; - ctx->root.offset = 0; - ctx->head = stack; - ctx->head->field = &ctx->root; - ctx->fmt_offset = 0; - ctx->head->parent_offset = 0; - ctx->packmode = '@'; - ctx->new_count = 1; - ctx->enc_count = 0; - ctx->enc_type = 0; - ctx->is_complex = 0; - while (type->typegroup == 'S') { - ++ctx->head; - ctx->head->field = type->fields; - ctx->head->parent_offset = 0; - type = type->fields->type; - } -} - -static int __Pyx_BufFmt_ParseNumber(const char** ts) { - int count; - const char* t = *ts; - if (*t < '0' || *t > '9') { - return -1; - } else { - count = *t++ - '0'; - while (*t >= '0' && *t < '9') { - count *= 10; - count += *t++ - '0'; - } - } - *ts = t; - return count; -} - -static void __Pyx_BufFmt_RaiseUnexpectedChar(char ch) { - char msg[] = {ch, 0}; - PyErr_Format(PyExc_ValueError, "Unexpected format string character: '%s'", msg); -} - -static const char* __Pyx_BufFmt_DescribeTypeChar(char ch, int is_complex) { - switch (ch) { - case 'b': return "'char'"; - case 'B': return "'unsigned char'"; - case 'h': return "'short'"; - case 'H': return "'unsigned short'"; - case 'i': return "'int'"; - case 'I': return "'unsigned int'"; - case 'l': return "'long'"; - case 'L': return "'unsigned long'"; - case 'q': return "'long long'"; - case 'Q': return "'unsigned long long'"; - case 'f': return (is_complex ? "'complex float'" : "'float'"); - case 'd': return (is_complex ? "'complex double'" : "'double'"); - case 'g': return (is_complex ? "'complex long double'" : "'long double'"); - case 'T': return "a struct"; - case 'O': return "Python object"; - case 'P': return "a pointer"; - case 0: return "end"; - default: return "unparseable format string"; - } -} - -static size_t __Pyx_BufFmt_TypeCharToStandardSize(char ch, int is_complex) { - switch (ch) { - case '?': case 'c': case 'b': case 'B': return 1; - case 'h': case 'H': return 2; - case 'i': case 'I': case 'l': case 'L': return 4; - case 'q': case 'Q': return 8; - case 'f': return (is_complex ? 8 : 4); - case 'd': return (is_complex ? 16 : 8); - case 'g': { - PyErr_SetString(PyExc_ValueError, "Python does not define a standard format string size for long double ('g').."); - return 0; - } - case 'O': case 'P': return sizeof(void*); - default: - __Pyx_BufFmt_RaiseUnexpectedChar(ch); - return 0; - } -} - -static size_t __Pyx_BufFmt_TypeCharToNativeSize(char ch, int is_complex) { - switch (ch) { - case 'c': case 'b': case 'B': return 1; - case 'h': case 'H': return sizeof(short); - case 'i': case 'I': return sizeof(int); - case 'l': case 'L': return sizeof(long); - #ifdef HAVE_LONG_LONG - case 'q': case 'Q': return sizeof(PY_LONG_LONG); - #endif - case 'f': return sizeof(float) * (is_complex ? 2 : 1); - case 'd': return sizeof(double) * (is_complex ? 2 : 1); - case 'g': return sizeof(long double) * (is_complex ? 2 : 1); - case 'O': case 'P': return sizeof(void*); - default: { - __Pyx_BufFmt_RaiseUnexpectedChar(ch); - return 0; - } - } -} - -typedef struct { char c; short x; } __Pyx_st_short; -typedef struct { char c; int x; } __Pyx_st_int; -typedef struct { char c; long x; } __Pyx_st_long; -typedef struct { char c; float x; } __Pyx_st_float; -typedef struct { char c; double x; } __Pyx_st_double; -typedef struct { char c; long double x; } __Pyx_st_longdouble; -typedef struct { char c; void *x; } __Pyx_st_void_p; -#ifdef HAVE_LONG_LONG -typedef struct { char c; PY_LONG_LONG x; } __Pyx_s_long_long; -#endif - -static size_t __Pyx_BufFmt_TypeCharToAlignment(char ch, int is_complex) { - switch (ch) { - case '?': case 'c': case 'b': case 'B': return 1; - case 'h': case 'H': return sizeof(__Pyx_st_short) - sizeof(short); - case 'i': case 'I': return sizeof(__Pyx_st_int) - sizeof(int); - case 'l': case 'L': return sizeof(__Pyx_st_long) - sizeof(long); -#ifdef HAVE_LONG_LONG - case 'q': case 'Q': return sizeof(__Pyx_s_long_long) - sizeof(PY_LONG_LONG); -#endif - case 'f': return sizeof(__Pyx_st_float) - sizeof(float); - case 'd': return sizeof(__Pyx_st_double) - sizeof(double); - case 'g': return sizeof(__Pyx_st_longdouble) - sizeof(long double); - case 'P': case 'O': return sizeof(__Pyx_st_void_p) - sizeof(void*); - default: - __Pyx_BufFmt_RaiseUnexpectedChar(ch); - return 0; - } -} - -static size_t __Pyx_BufFmt_TypeCharToGroup(char ch, int is_complex) { - switch (ch) { - case 'c': case 'b': case 'h': case 'i': case 'l': case 'q': return 'I'; - case 'B': case 'H': case 'I': case 'L': case 'Q': return 'U'; - case 'f': case 'd': case 'g': return (is_complex ? 'C' : 'R'); - case 'O': return 'O'; - case 'P': return 'P'; - default: { - __Pyx_BufFmt_RaiseUnexpectedChar(ch); - return 0; - } - } -} - -static void __Pyx_BufFmt_RaiseExpected(__Pyx_BufFmt_Context* ctx) { - if (ctx->head == NULL || ctx->head->field == &ctx->root) { - const char* expected; - const char* quote; - if (ctx->head == NULL) { - expected = "end"; - quote = ""; - } else { - expected = ctx->head->field->type->name; - quote = "'"; - } - PyErr_Format(PyExc_ValueError, - "Buffer dtype mismatch, expected %s%s%s but got %s", - quote, expected, quote, - __Pyx_BufFmt_DescribeTypeChar(ctx->enc_type, ctx->is_complex)); - } else { - __Pyx_StructField* field = ctx->head->field; - __Pyx_StructField* parent = (ctx->head - 1)->field; - PyErr_Format(PyExc_ValueError, - "Buffer dtype mismatch, expected '%s' but got %s in '%s.%s'", - field->type->name, __Pyx_BufFmt_DescribeTypeChar(ctx->enc_type, ctx->is_complex), - parent->type->name, field->name); - } -} - -static int __Pyx_BufFmt_ProcessTypeChunk(__Pyx_BufFmt_Context* ctx) { - char group; - size_t size, offset; - if (ctx->enc_type == 0) return 0; - group = __Pyx_BufFmt_TypeCharToGroup(ctx->enc_type, ctx->is_complex); - do { - __Pyx_StructField* field = ctx->head->field; - __Pyx_TypeInfo* type = field->type; - - if (ctx->packmode == '@' || ctx->packmode == '^') { - size = __Pyx_BufFmt_TypeCharToNativeSize(ctx->enc_type, ctx->is_complex); - } else { - size = __Pyx_BufFmt_TypeCharToStandardSize(ctx->enc_type, ctx->is_complex); - } - if (ctx->packmode == '@') { - int align_at = __Pyx_BufFmt_TypeCharToAlignment(ctx->enc_type, ctx->is_complex); - int align_mod_offset; - if (align_at == 0) return -1; - align_mod_offset = ctx->fmt_offset % align_at; - if (align_mod_offset > 0) ctx->fmt_offset += align_at - align_mod_offset; - } - - if (type->size != size || type->typegroup != group) { - if (type->typegroup == 'C' && type->fields != NULL) { - /* special case -- treat as struct rather than complex number */ - size_t parent_offset = ctx->head->parent_offset + field->offset; - ++ctx->head; - ctx->head->field = type->fields; - ctx->head->parent_offset = parent_offset; - continue; - } - - __Pyx_BufFmt_RaiseExpected(ctx); - return -1; - } - - offset = ctx->head->parent_offset + field->offset; - if (ctx->fmt_offset != offset) { - PyErr_Format(PyExc_ValueError, - "Buffer dtype mismatch; next field is at offset %"PY_FORMAT_SIZE_T"d " - "but %"PY_FORMAT_SIZE_T"d expected", ctx->fmt_offset, offset); - return -1; - } - - ctx->fmt_offset += size; - - --ctx->enc_count; /* Consume from buffer string */ - - /* Done checking, move to next field, pushing or popping struct stack if needed */ - while (1) { - if (field == &ctx->root) { - ctx->head = NULL; - if (ctx->enc_count != 0) { - __Pyx_BufFmt_RaiseExpected(ctx); - return -1; - } - break; /* breaks both loops as ctx->enc_count == 0 */ - } - ctx->head->field = ++field; - if (field->type == NULL) { - --ctx->head; - field = ctx->head->field; - continue; - } else if (field->type->typegroup == 'S') { - size_t parent_offset = ctx->head->parent_offset + field->offset; - if (field->type->fields->type == NULL) continue; /* empty struct */ - field = field->type->fields; - ++ctx->head; - ctx->head->field = field; - ctx->head->parent_offset = parent_offset; - break; - } else { - break; - } - } - } while (ctx->enc_count); - ctx->enc_type = 0; - ctx->is_complex = 0; - return 0; -} - -static int __Pyx_BufFmt_FirstPack(__Pyx_BufFmt_Context* ctx) { - if (ctx->enc_type != 0 || ctx->packmode != '@') { - PyErr_SetString(PyExc_ValueError, "Buffer packing mode currently only allowed at beginning of format string (this is a defect)"); - return -1; - } - return 0; -} - -static const char* __Pyx_BufFmt_CheckString(__Pyx_BufFmt_Context* ctx, const char* ts) { - int got_Z = 0; - while (1) { - switch(*ts) { - case 0: - if (ctx->enc_type != 0 && ctx->head == NULL) { - __Pyx_BufFmt_RaiseExpected(ctx); - return NULL; - } - if (__Pyx_BufFmt_ProcessTypeChunk(ctx) == -1) return NULL; - if (ctx->head != NULL) { - __Pyx_BufFmt_RaiseExpected(ctx); - return NULL; - } - return ts; - case ' ': - case 10: - case 13: - ++ts; - break; - case '<': - if (!__Pyx_IsLittleEndian()) { - PyErr_SetString(PyExc_ValueError, "Little-endian buffer not supported on big-endian compiler"); - return NULL; - } - if (__Pyx_BufFmt_FirstPack(ctx) == -1) return NULL; - ctx->packmode = '='; - ++ts; - break; - case '>': - case '!': - if (__Pyx_IsLittleEndian()) { - PyErr_SetString(PyExc_ValueError, "Big-endian buffer not supported on little-endian compiler"); - return NULL; - } - if (__Pyx_BufFmt_FirstPack(ctx) == -1) return NULL; - ctx->packmode = '='; - ++ts; - break; - case '=': - case '@': - case '^': - if (__Pyx_BufFmt_FirstPack(ctx) == -1) return NULL; - ctx->packmode = *ts++; - break; - case 'T': /* substruct */ - { - int i; - const char* ts_after_sub; - int struct_count = ctx->new_count; - ctx->new_count = 1; - ++ts; - if (*ts != '{') { - PyErr_SetString(PyExc_ValueError, "Buffer acquisition: Expected '{' after 'T'"); - return NULL; - } - ++ts; - ts_after_sub = ts; - for (i = 0; i != struct_count; ++i) { - ts_after_sub = __Pyx_BufFmt_CheckString(ctx, ts); - if (!ts_after_sub) return NULL; - } - ts = ts_after_sub; - } - break; - case '}': /* end of substruct; either repeat or move on */ - ++ts; - return ts; - case 'x': - if (__Pyx_BufFmt_ProcessTypeChunk(ctx) == -1) return NULL; - ctx->fmt_offset += ctx->new_count; - ctx->new_count = 1; - ctx->enc_count = 0; - ctx->enc_type = 0; - ++ts; - break; - case 'Z': - got_Z = 1; - ++ts; - if (*ts != 'f' && *ts != 'd' && *ts != 'g') { - __Pyx_BufFmt_RaiseUnexpectedChar('Z'); - return NULL; - } /* fall through */ - case 'c': case 'b': case 'B': case 'h': case 'H': case 'i': case 'I': - case 'l': case 'L': case 'q': case 'Q': - case 'f': case 'd': case 'g': - case 'O': - if (ctx->enc_type == *ts && got_Z == ctx->is_complex) { - /* Continue pooling same type */ - ctx->enc_count += ctx->new_count; - } else { - /* New type */ - if (__Pyx_BufFmt_ProcessTypeChunk(ctx) == -1) return NULL; - ctx->enc_count = ctx->new_count; - ctx->enc_type = *ts; - ctx->is_complex = got_Z; - } - ++ts; - ctx->new_count = 1; - got_Z = 0; - break; - default: - { - ctx->new_count = __Pyx_BufFmt_ParseNumber(&ts); - if (ctx->new_count == -1) { /* First char was not a digit */ - char msg[2] = { *ts, 0 }; - PyErr_Format(PyExc_ValueError, - "Does not understand character buffer dtype format string ('%s')", msg); - return NULL; - } - } - - } - } -} - -static CYTHON_INLINE void __Pyx_ZeroBuffer(Py_buffer* buf) { - buf->buf = NULL; - buf->obj = NULL; - buf->strides = __Pyx_zeros; - buf->shape = __Pyx_zeros; - buf->suboffsets = __Pyx_minusones; -} - -static int __Pyx_GetBufferAndValidate(Py_buffer* buf, PyObject* obj, __Pyx_TypeInfo* dtype, int flags, int nd, int cast, __Pyx_BufFmt_StackElem* stack) { - if (obj == Py_None) { - __Pyx_ZeroBuffer(buf); - return 0; - } - buf->buf = NULL; - if (__Pyx_GetBuffer(obj, buf, flags) == -1) goto fail; - if (buf->ndim != nd) { - PyErr_Format(PyExc_ValueError, - "Buffer has wrong number of dimensions (expected %d, got %d)", - nd, buf->ndim); - goto fail; - } - if (!cast) { - __Pyx_BufFmt_Context ctx; - __Pyx_BufFmt_Init(&ctx, stack, dtype); - if (!__Pyx_BufFmt_CheckString(&ctx, buf->format)) goto fail; - } - if ((unsigned)buf->itemsize != dtype->size) { - PyErr_Format(PyExc_ValueError, - "Item size of buffer (%"PY_FORMAT_SIZE_T"d byte%s) does not match size of '%s' (%"PY_FORMAT_SIZE_T"d byte%s)", - buf->itemsize, (buf->itemsize > 1) ? "s" : "", - dtype->name, - dtype->size, (dtype->size > 1) ? "s" : ""); - goto fail; - } - if (buf->suboffsets == NULL) buf->suboffsets = __Pyx_minusones; - return 0; -fail:; - __Pyx_ZeroBuffer(buf); - return -1; -} - -static CYTHON_INLINE void __Pyx_SafeReleaseBuffer(Py_buffer* info) { - if (info->buf == NULL) return; - if (info->suboffsets == __Pyx_minusones) info->suboffsets = NULL; - __Pyx_ReleaseBuffer(info); -} - -static CYTHON_INLINE void __Pyx_ErrRestore(PyObject *type, PyObject *value, PyObject *tb) { - PyObject *tmp_type, *tmp_value, *tmp_tb; - PyThreadState *tstate = PyThreadState_GET(); - - tmp_type = tstate->curexc_type; - tmp_value = tstate->curexc_value; - tmp_tb = tstate->curexc_traceback; - tstate->curexc_type = type; - tstate->curexc_value = value; - tstate->curexc_traceback = tb; - Py_XDECREF(tmp_type); - Py_XDECREF(tmp_value); - Py_XDECREF(tmp_tb); -} - -static CYTHON_INLINE void __Pyx_ErrFetch(PyObject **type, PyObject **value, PyObject **tb) { - PyThreadState *tstate = PyThreadState_GET(); - *type = tstate->curexc_type; - *value = tstate->curexc_value; - *tb = tstate->curexc_traceback; - - tstate->curexc_type = 0; - tstate->curexc_value = 0; - tstate->curexc_traceback = 0; -} - - static CYTHON_INLINE void __Pyx_RaiseNeedMoreValuesError(Py_ssize_t index) { PyErr_Format(PyExc_ValueError, #if PY_VERSION_HEX < 0x02050000 @@ -4143,30 +3794,6 @@ static int __Pyx_ArgTypeTest(PyObject *obj, PyTypeObject *type, int none_allowed return 0; } -#if PY_MAJOR_VERSION < 3 -static int __Pyx_GetBuffer(PyObject *obj, Py_buffer *view, int flags) { - #if PY_VERSION_HEX >= 0x02060000 - if (Py_TYPE(obj)->tp_flags & Py_TPFLAGS_HAVE_NEWBUFFER) - return PyObject_GetBuffer(obj, view, flags); - #endif - if (PyObject_TypeCheck(obj, __pyx_ptype_5numpy_ndarray)) return __pyx_pf_5numpy_7ndarray___getbuffer__(obj, view, flags); - else { - PyErr_Format(PyExc_TypeError, "'%100s' does not have the buffer interface", Py_TYPE(obj)->tp_name); - return -1; - } -} - -static void __Pyx_ReleaseBuffer(Py_buffer *view) { - PyObject* obj = view->obj; - if (obj) { -if (PyObject_TypeCheck(obj, __pyx_ptype_5numpy_ndarray)) __pyx_pf_5numpy_7ndarray___releasebuffer__(obj, view); - Py_DECREF(obj); - view->obj = NULL; - } -} - -#endif - static PyObject *__Pyx_Import(PyObject *name, PyObject *from_list) { PyObject *__import__ = 0; PyObject *empty_list = 0; @@ -4208,6 +3835,33 @@ static PyObject *__Pyx_GetName(PyObject *dict, PyObject *name) { return result; } +static CYTHON_INLINE void __Pyx_ErrRestore(PyObject *type, PyObject *value, PyObject *tb) { + PyObject *tmp_type, *tmp_value, *tmp_tb; + PyThreadState *tstate = PyThreadState_GET(); + + tmp_type = tstate->curexc_type; + tmp_value = tstate->curexc_value; + tmp_tb = tstate->curexc_traceback; + tstate->curexc_type = type; + tstate->curexc_value = value; + tstate->curexc_traceback = tb; + Py_XDECREF(tmp_type); + Py_XDECREF(tmp_value); + Py_XDECREF(tmp_tb); +} + +static CYTHON_INLINE void __Pyx_ErrFetch(PyObject **type, PyObject **value, PyObject **tb) { + PyThreadState *tstate = PyThreadState_GET(); + *type = tstate->curexc_type; + *value = tstate->curexc_value; + *tb = tstate->curexc_traceback; + + tstate->curexc_type = 0; + tstate->curexc_value = 0; + tstate->curexc_traceback = 0; +} + + #if PY_MAJOR_VERSION < 3 static void __Pyx_Raise(PyObject *type, PyObject *value, PyObject *tb) { Py_XINCREF(type); diff --git a/scikits/learn/utils/arrayfuncs.pyx b/scikits/learn/utils/arrayfuncs.pyx index eaffb35208a82769fcbab51642d14681fb569dd5..d75f6cc447d4d3016b4b7cc803e05d676ec2af93 100644 --- a/scikits/learn/utils/arrayfuncs.pyx +++ b/scikits/learn/utils/arrayfuncs.pyx @@ -28,13 +28,18 @@ cdef extern from "cblas.h": int N, double *A, int lda, double *X, int incX) + void cblas_strsv(CBLAS_ORDER Order, CBLAS_UPLO Uplo, + CBLAS_TRANSPOSE TransA, CBLAS_DIAG Diag, + int N, float *A, int lda, float *X, + int incX) cdef extern from "float.h": cdef double DBL_MAX cdef float FLT_MAX cdef extern from "src/cholesky_delete.c": - int c_cholesky_delete "cholesky_delete" (int m, int n, double *L, int go_out) + int double_cholesky_delete (int m, int n, double *L, int go_out) + int float_cholesky_delete (int m, int n, float *L, int go_out) ctypedef np.float64_t DOUBLE @@ -70,26 +75,41 @@ cdef double _double_min_pos (double *X, Py_ssize_t size): min_val = X[i] return min_val -def solve_triangular ( - np.ndarray[DOUBLE, ndim=2] X, - np.ndarray[DOUBLE, ndim=1] y - ): +def solve_triangular (np.ndarray X, np.ndarray y): """ Solves a triangular system (overwrites y) - TODO: option for upper, lower, etc. + Note: The lapack function to solve triangular systems was added to + scipy v0.9. Remove this when we stop supporting earlier versions. """ - cdef int lda = <int> X.strides[0] / sizeof(double) + cdef int lda + + if X.dtype.name == 'float64' and y.dtype.name == 'float64': + lda = <int> X.strides[0] / sizeof(double) + + cblas_dtrsv (CblasRowMajor, CblasLower, CblasNoTrans, + CblasNonUnit, <int> X.shape[0], <double *> X.data, + lda, <double *> y.data, 1); - cblas_dtrsv (CblasRowMajor, CblasLower, CblasNoTrans, - CblasNonUnit, <int> X.shape[0], <double *> X.data, - lda, <double *> y.data, 1); + elif X.dtype.name == 'float32' and y.dtype.name == 'float32': + lda = <int> X.strides[0] / sizeof(float) + cblas_strsv (CblasRowMajor, CblasLower, CblasNoTrans, + CblasNonUnit, <int> X.shape[0], <float *> X.data, + lda, <float *> y.data, 1); + else: + raise ValueError ('Unsupported or inconsistent dtype in arrays X, y') -def cholesky_delete ( - np.ndarray[DOUBLE, ndim=2] L, - int go_out): + +def cholesky_delete (np.ndarray L, int go_out): cdef int n = <int> L.shape[0] - cdef int m = <int> L.strides[0] / sizeof (double) - c_cholesky_delete (m, n, <double *> L.data, go_out) + cdef int m + + if L.dtype.name == 'float64': + m = <int> L.strides[0] / sizeof (double) + double_cholesky_delete (m, n, <double *> L.data, go_out) + elif L.dtype.name == 'float32': + m = <int> L.strides[0] / sizeof (float) + float_cholesky_delete (m, n, <float *> L.data, go_out) + diff --git a/scikits/learn/utils/src/cholesky_delete.c b/scikits/learn/utils/src/cholesky_delete.c index fe029de4b8b28ca79a15d5fecae4bc5845d51d25..7117ba0a5c14de7a04537a15f42bc174e7fbd994 100644 --- a/scikits/learn/utils/src/cholesky_delete.c +++ b/scikits/learn/utils/src/cholesky_delete.c @@ -13,7 +13,7 @@ * TODO: put transpose as an option * */ -int cholesky_delete (int m, int n, double *L, int go_out) { +int double_cholesky_delete (int m, int n, double *L, int go_out) { double c, s; @@ -45,3 +45,37 @@ int cholesky_delete (int m, int n, double *L, int go_out) { return 0; } + + +int float_cholesky_delete (int m, int n, float *L, int go_out) { + + float c, s; + + /* delete row go_out */ + float * _L = L + (go_out * m); + int i; + for (i = go_out; i < n - 1; ++i) { + cblas_scopy (i + 2, _L + m , 1, _L, 1); + _L += m; + } + + _L = L + (go_out * m); + for (i=go_out; i < n - 1; ++i) { + + cblas_srotg (_L + i, _L + i + 1, &c, &s); + if (_L[i] < 0) { + /* Diagonals cannot be negative */ + /* _L[i] = copysign(_L[i], 1.0); */ + c = -c; + s = -s; + } + _L[i+1] = 0.; /* just for cleanup */ + _L += m; + + cblas_srot (n - (i + 2), _L + i, m, _L + i + 1, + m, c, s); + + } + + return 0; +}