Skip to content

Commit a22ee97

Browse files
martin-frbgmartin-frbg
authored
Merge branch 'OpenMathLib:develop' into lapack1286
2 parents 100988f + dd3ba97 commit a22ee97

32 files changed

Lines changed: 5802 additions & 5269 deletions

‎Makefile.system‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1648,7 +1648,11 @@ ifneq ($(C_COMPILER), SUN)
16481648
CCOMMON_OPT += -UASMNAME -UASMFNAME -UNAME -UCNAME -UCHAR_NAME -UCHAR_CNAME
16491649
endif
16501650
endif
1651+
ifeq ( ,$(findstring pizfix,$(CEXTRALIB)))
16511652
CCOMMON_OPT += -DASMNAME=$(FU)$(*F) -DASMFNAME=$(FU)$(*F)$(BU) -DNAME=$(*F)$(BU) -DCNAME=$(*F) -DCHAR_NAME=\"$(*F)$(BU)\" -DCHAR_CNAME=\"$(*F)\"
1653+
else
1654+
CCOMMON_OPT += -DASMNAME=$(FU)pizlonated_$(*F) -DASMFNAME=$(FU)pizlonated_$(*F)$(BU) -DNAME=$(*F)$(BU) -DCNAME=$(*F) -DCHAR_NAME=\"$(*F)$(BU)\" -DCHAR_CNAME=\"$(*F)\"
1655+
endif
16521656

16531657
ifeq ($(CORE), PPC440)
16541658
CCOMMON_OPT += -DALLOC_QALLOC

‎common_linux.h‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -70,7 +70,7 @@ extern long int syscall (long int __sysno, ...);
7070
static inline int my_mbind(void *addr, unsigned long len, int mode,
7171
unsigned long *nodemask, unsigned long maxnode,
7272
unsigned flags) {
73-
#if defined (__LSB_VERSION__) || defined(ARCH_ZARCH)
73+
#if defined (__LSB_VERSION__) || defined(ARCH_ZARCH) || defined(__PIZLONATOR_WAS_HERE__)
7474
// So far, LSB (Linux Standard Base) don't support syscall().
7575
// https://lsbbugs.linuxfoundation.org/show_bug.cgi?id=3482
7676
return 0;

‎common_x86_64.h‎

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -83,11 +83,15 @@ static __inline void blas_lock(volatile BLASULONG *address){
8383
while (*address) {YIELDING;}
8484

8585
#ifndef C_MSVC
86+
#ifdef __PIZLONATOR_WAS_HERE__
87+
ret=__atomic_exchange_n(address, 1, __ATOMIC_SEQ_CST);
88+
#else
8689
__asm__ __volatile__(
8790
"xchgl %0, %1\n"
8891
: "=r"(ret), "=m"(*address)
8992
: "0"(1), "m"(*address)
9093
: "memory");
94+
#endif
9195
#else
9296
ret=InterlockedExchange64((volatile LONG64 *)(address), 1);
9397
#endif
@@ -237,7 +241,7 @@ static __inline unsigned int blas_quickdivide(unsigned int x, unsigned int y){
237241

238242
y = blas_quick_divide_table[y];
239243

240-
__asm__ __volatile__ ("mull %0" :"=d" (result), "+a"(x) : "0" (y));
244+
__asm__ __volatile__ ("mull %0" :"=d" (result), "+a"(x) : "0" (y) : "cc");
241245
return result;
242246
}
243247
#endif

‎lapack-netlib/SRC/cbbcsd.f‎

Lines changed: 63 additions & 30 deletions
Original file line numberDiff line numberDiff line change
@@ -5,15 +5,13 @@
55
* Online html documentation available at
66
* http://www.netlib.org/lapack/explore-html/
77
*
8-
*> \htmlonly
98
*> Download CBBCSD + dependencies
109
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cbbcsd.f">
1110
*> [TGZ]</a>
1211
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cbbcsd.f">
1312
*> [ZIP]</a>
1413
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cbbcsd.f">
1514
*> [TXT]</a>
16-
*> \endhtmlonly
1715
*
1816
* Definition:
1917
* ===========
@@ -322,13 +320,15 @@
322320
*> \author Univ. of Colorado Denver
323321
*> \author NAG Ltd.
324322
*
325-
*> \ingroup complexOTHERcomputational
323+
*> \ingroup bbcsd
326324
*
327325
* =====================================================================
328-
SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
326+
SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
327+
$ Q,
329328
$ THETA, PHI, U1, LDU1, U2, LDU2, V1T, LDV1T,
330329
$ V2T, LDV2T, B11D, B11E, B12D, B12E, B21D, B21E,
331330
$ B22D, B22E, RWORK, LRWORK, INFO )
331+
IMPLICIT NONE
332332
*
333333
* -- LAPACK computational routine --
334334
* -- LAPACK is a software package provided by Univ. of Tennessee, --
@@ -372,7 +372,8 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
372372
$ UNFL, X1, X2, Y1, Y2
373373
*
374374
* .. External Subroutines ..
375-
EXTERNAL CLASR, CSCAL, CSWAP, SLARTGP, SLARTGS, SLAS2,
375+
EXTERNAL CLASR, CSCAL, CSWAP, SLARTGP, SLARTGS,
376+
$ SLAS2,
376377
$ XERBLA
377378
* ..
378379
* .. External Functions ..
@@ -417,7 +418,7 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
417418
*
418419
IF( INFO .EQ. 0 .AND. Q .EQ. 0 ) THEN
419420
LRWORKMIN = 1
420-
RWORK(1) = LRWORKMIN
421+
RWORK(1) = REAL( LRWORKMIN )
421422
RETURN
422423
END IF
423424
*
@@ -434,7 +435,7 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
434435
IV2TSN = IV2TCS + Q
435436
LRWORKOPT = IV2TSN + Q - 1
436437
LRWORKMIN = LRWORKOPT
437-
RWORK(1) = LRWORKOPT
438+
RWORK(1) = REAL( LRWORKOPT )
438439
IF( LRWORK .LT. LRWORKMIN .AND. .NOT. LQUERY ) THEN
439440
INFO = -28
440441
END IF
@@ -453,7 +454,7 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
453454
UNFL = SLAMCH( 'Safe minimum' )
454455
TOLMUL = MAX( TEN, MIN( HUNDRED, EPS**MEIGHTH ) )
455456
TOL = TOLMUL*EPS
456-
THRESH = MAX( TOL, MAXITR*Q*Q*UNFL )
457+
THRESH = MAX( TOL, REAL( MAXITR*Q*Q )*UNFL )
457458
*
458459
* Test for negligible sines or cosines
459460
*
@@ -559,9 +560,11 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
559560
*
560561
* Compute shifts for B11 and B21 and use the lesser
561562
*
562-
CALL SLAS2( B11D(IMAX-1), B11E(IMAX-1), B11D(IMAX), SIGMA11,
563+
CALL SLAS2( B11D(IMAX-1), B11E(IMAX-1), B11D(IMAX),
564+
$ SIGMA11,
563565
$ DUMMY )
564-
CALL SLAS2( B21D(IMAX-1), B21E(IMAX-1), B21D(IMAX), SIGMA21,
566+
CALL SLAS2( B21D(IMAX-1), B21E(IMAX-1), B21D(IMAX),
567+
$ SIGMA21,
565568
$ DUMMY )
566569
*
567570
IF( SIGMA11 .LE. SIGMA21 ) THEN
@@ -613,7 +616,9 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
613616
*
614617
* Chase the bulges in B11(IMIN+1,IMIN) and B21(IMIN+1,IMIN)
615618
*
616-
IF( B11D(IMIN)**2+B11BULGE**2 .GT. THRESH**2 ) THEN
619+
IF( B11D(IMIN)**2+B11BULGE**2 .GT.
620+
$ (THRESH*MAX( ABS(B11D(IMIN)),
621+
$ ABS(B11D(IMIN+1)), UNFL ))**2 ) THEN
617622
CALL SLARTGP( B11BULGE, B11D(IMIN), RWORK(IU1SN+IMIN-1),
618623
$ RWORK(IU1CS+IMIN-1), R )
619624
ELSE IF( MU .LE. NU ) THEN
@@ -623,7 +628,9 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
623628
CALL SLARTGS( B12D( IMIN ), B12E( IMIN ), NU,
624629
$ RWORK(IU1CS+IMIN-1), RWORK(IU1SN+IMIN-1) )
625630
END IF
626-
IF( B21D(IMIN)**2+B21BULGE**2 .GT. THRESH**2 ) THEN
631+
IF( B21D(IMIN)**2+B21BULGE**2 .GT.
632+
$ (THRESH*MAX( ABS(B21D(IMIN)),
633+
$ ABS(B21D(IMIN+1)), UNFL ))**2 ) THEN
627634
CALL SLARTGP( B21BULGE, B21D(IMIN), RWORK(IU2SN+IMIN-1),
628635
$ RWORK(IU2CS+IMIN-1), R )
629636
ELSE IF( NU .LT. MU ) THEN
@@ -687,10 +694,18 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
687694
* Determine if there are bulges to chase or if a new direct
688695
* summand has been reached
689696
*
690-
RESTART11 = B11E(I-1)**2 + B11BULGE**2 .LE. THRESH**2
691-
RESTART21 = B21E(I-1)**2 + B21BULGE**2 .LE. THRESH**2
692-
RESTART12 = B12D(I-1)**2 + B12BULGE**2 .LE. THRESH**2
693-
RESTART22 = B22D(I-1)**2 + B22BULGE**2 .LE. THRESH**2
697+
RESTART11 = B11E(I-1)**2 + B11BULGE**2 .LE.
698+
$ (THRESH*MAX( ABS(B11D(I-1)), ABS(B11D(I)),
699+
$ UNFL ))**2
700+
RESTART21 = B21E(I-1)**2 + B21BULGE**2 .LE.
701+
$ (THRESH*MAX( ABS(B21D(I-1)), ABS(B21D(I)),
702+
$ UNFL ))**2
703+
RESTART12 = B12D(I-1)**2 + B12BULGE**2 .LE.
704+
$ (THRESH*MAX( ABS(B12E(I-1)), ABS(B12D(I)),
705+
$ UNFL ))**2
706+
RESTART22 = B22D(I-1)**2 + B22BULGE**2 .LE.
707+
$ (THRESH*MAX( ABS(B22E(I-1)), ABS(B22D(I)),
708+
$ UNFL ))**2
694709
*
695710
* If possible, chase bulges from B11(I-1,I+1), B12(I-1,I),
696711
* B21(I-1,I+1), and B22(I-1,I). If necessary, restart bulge-
@@ -718,10 +733,12 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
718733
CALL SLARTGP( Y2, Y1, RWORK(IV2TSN+I-1-1),
719734
$ RWORK(IV2TCS+I-1-1), R )
720735
ELSE IF( .NOT. RESTART12 .AND. RESTART22 ) THEN
721-
CALL SLARTGP( B12BULGE, B12D(I-1), RWORK(IV2TSN+I-1-1),
736+
CALL SLARTGP( B12BULGE, B12D(I-1),
737+
$ RWORK(IV2TSN+I-1-1),
722738
$ RWORK(IV2TCS+I-1-1), R )
723739
ELSE IF( RESTART12 .AND. .NOT. RESTART22 ) THEN
724-
CALL SLARTGP( B22BULGE, B22D(I-1), RWORK(IV2TSN+I-1-1),
740+
CALL SLARTGP( B22BULGE, B22D(I-1),
741+
$ RWORK(IV2TSN+I-1-1),
725742
$ RWORK(IV2TCS+I-1-1), R )
726743
ELSE IF( NU .LT. MU ) THEN
727744
CALL SLARTGS( B12E(I-1), B12D(I), NU,
@@ -770,17 +787,26 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
770787
* Determine if there are bulges to chase or if a new direct
771788
* summand has been reached
772789
*
773-
RESTART11 = B11D(I)**2 + B11BULGE**2 .LE. THRESH**2
774-
RESTART12 = B12E(I-1)**2 + B12BULGE**2 .LE. THRESH**2
775-
RESTART21 = B21D(I)**2 + B21BULGE**2 .LE. THRESH**2
776-
RESTART22 = B22E(I-1)**2 + B22BULGE**2 .LE. THRESH**2
790+
RESTART11 = B11D(I)**2 + B11BULGE**2 .LE.
791+
$ (THRESH*MAX( ABS(B11E(I)), ABS(B11D(I+1)),
792+
$ UNFL ))**2
793+
RESTART12 = B12E(I-1)**2 + B12BULGE**2 .LE.
794+
$ (THRESH*MAX( ABS(B12D(I)), ABS(B12E(I)),
795+
$ UNFL ))**2
796+
RESTART21 = B21D(I)**2 + B21BULGE**2 .LE.
797+
$ (THRESH*MAX( ABS(B21E(I)), ABS(B21D(I+1)),
798+
$ UNFL ))**2
799+
RESTART22 = B22E(I-1)**2 + B22BULGE**2 .LE.
800+
$ (THRESH*MAX( ABS(B22D(I)), ABS(B22E(I)),
801+
$ UNFL ))**2
777802
*
778803
* If possible, chase bulges from B11(I+1,I), B12(I+1,I-1),
779804
* B21(I+1,I), and B22(I+1,I-1). If necessary, restart bulge-
780805
* chasing by applying the original shift again.
781806
*
782807
IF( .NOT. RESTART11 .AND. .NOT. RESTART12 ) THEN
783-
CALL SLARTGP( X2, X1, RWORK(IU1SN+I-1), RWORK(IU1CS+I-1),
808+
CALL SLARTGP( X2, X1, RWORK(IU1SN+I-1),
809+
$ RWORK(IU1CS+I-1),
784810
$ R )
785811
ELSE IF( .NOT. RESTART11 .AND. RESTART12 ) THEN
786812
CALL SLARTGP( B11BULGE, B11D(I), RWORK(IU1SN+I-1),
@@ -789,14 +815,16 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
789815
CALL SLARTGP( B12BULGE, B12E(I-1), RWORK(IU1SN+I-1),
790816
$ RWORK(IU1CS+I-1), R )
791817
ELSE IF( MU .LE. NU ) THEN
792-
CALL SLARTGS( B11E(I), B11D(I+1), MU, RWORK(IU1CS+I-1),
818+
CALL SLARTGS( B11E(I), B11D(I+1), MU,
819+
$ RWORK(IU1CS+I-1),
793820
$ RWORK(IU1SN+I-1) )
794821
ELSE
795822
CALL SLARTGS( B12D(I), B12E(I), NU, RWORK(IU1CS+I-1),
796823
$ RWORK(IU1SN+I-1) )
797824
END IF
798825
IF( .NOT. RESTART21 .AND. .NOT. RESTART22 ) THEN
799-
CALL SLARTGP( Y2, Y1, RWORK(IU2SN+I-1), RWORK(IU2CS+I-1),
826+
CALL SLARTGP( Y2, Y1, RWORK(IU2SN+I-1),
827+
$ RWORK(IU2CS+I-1),
800828
$ R )
801829
ELSE IF( .NOT. RESTART21 .AND. RESTART22 ) THEN
802830
CALL SLARTGP( B21BULGE, B21D(I), RWORK(IU2SN+I-1),
@@ -805,7 +833,8 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
805833
CALL SLARTGP( B22BULGE, B22E(I-1), RWORK(IU2SN+I-1),
806834
$ RWORK(IU2CS+I-1), R )
807835
ELSE IF( NU .LT. MU ) THEN
808-
CALL SLARTGS( B21E(I), B21D(I+1), NU, RWORK(IU2CS+I-1),
836+
CALL SLARTGS( B21E(I), B21D(I+1), NU,
837+
$ RWORK(IU2CS+I-1),
809838
$ RWORK(IU2SN+I-1) )
810839
ELSE
811840
CALL SLARTGS( B22D(I), B22E(I), MU, RWORK(IU2CS+I-1),
@@ -857,8 +886,10 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
857886
*
858887
* Chase bulges from B12(IMAX-1,IMAX) and B22(IMAX-1,IMAX)
859888
*
860-
RESTART12 = B12D(IMAX-1)**2 + B12BULGE**2 .LE. THRESH**2
861-
RESTART22 = B22D(IMAX-1)**2 + B22BULGE**2 .LE. THRESH**2
889+
RESTART12 = B12D(IMAX-1)**2 + B12BULGE**2 .LE.
890+
$ (THRESH*MAX( ABS(B12E(IMAX-1)), UNFL ))**2
891+
RESTART22 = B22D(IMAX-1)**2 + B22BULGE**2 .LE.
892+
$ (THRESH*MAX( ABS(B22E(IMAX-1)), UNFL ))**2
862893
*
863894
IF( .NOT. RESTART12 .AND. .NOT. RESTART22 ) THEN
864895
CALL SLARTGP( Y2, Y1, RWORK(IV2TSN+IMAX-1-1),
@@ -991,7 +1022,8 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
9911022
IF( B12D(IMAX)+B22D(IMAX) .LT. 0 ) THEN
9921023
IF( WANTV2T ) THEN
9931024
IF( COLMAJOR ) THEN
994-
CALL CSCAL( M-Q, NEGONECOMPLEX, V2T(IMAX,1), LDV2T )
1025+
CALL CSCAL( M-Q, NEGONECOMPLEX, V2T(IMAX,1),
1026+
$ LDV2T )
9951027
ELSE
9961028
CALL CSCAL( M-Q, NEGONECOMPLEX, V2T(1,IMAX), 1 )
9971029
END IF
@@ -1058,7 +1090,8 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P, Q,
10581090
IF( WANTU2 )
10591091
$ CALL CSWAP( M-P, U2(1,I), 1, U2(1,MINI), 1 )
10601092
IF( WANTV1T )
1061-
$ CALL CSWAP( Q, V1T(I,1), LDV1T, V1T(MINI,1), LDV1T )
1093+
$ CALL CSWAP( Q, V1T(I,1), LDV1T, V1T(MINI,1),
1094+
$ LDV1T )
10621095
IF( WANTV2T )
10631096
$ CALL CSWAP( M-Q, V2T(I,1), LDV2T, V2T(MINI,1),
10641097
$ LDV2T )

0 commit comments

Comments
 (0)