1，源码
SUBROUTINE INIT_MATRIX(A, m, n, lda)DOUBLE PRECISION A(*)CALL SRAND(2024)DO i=1, mDO j=1, nA(i + lda*(j-1)) = RAND() + RAND()
C WRITE(*, '(F8.4)') A(i)END DOEND DOENDSUBROUTINE PRINT_MATRIX(A, m, n, lda)DOUBLE PRECISION A(*)DO i=1, mDO j=1, n
C WRITE(*, '(F8.4 , $)') A(i + lda*(j-1))WRITE(*, 20) A(i + lda*(j-1))
20 FORMAT (1X,F8.4, $)END DOPRINT *,''END DOENDPROGRAM MATRIX_EXINTEGER m, n, ldaDOUBLE PRECISION A(20)lda = 4m = 4n = 5CALL INIT_MATRIX(A, m, n, lda)PRINT *, "A ="CALL PRINT_MATRIX(A, m, n, lda)PRINT *, "A(11) ="PRINT *, A(11)ENDMakefile
matrix_ex: matrix_ex.fgfortran -g $< -o $@.PHONY: clean
clean:-rm matrix_ex2，执行
make./matrix_ex效果:
3, dgesvd example
* Copyright (C) 2009-2015 Intel Corporation. All Rights Reserved.
* The information and material ("Material") provided below is owned by Intel
* Corporation or its suppliers or licensors, and title to such Material remains
* with Intel Corporation or its suppliers or licensors. The Material contains
* proprietary information of Intel or its suppliers and licensors. The Material
* is protected by worldwide copyright laws and treaty provisions. No part of
* the Material may be copied, reproduced, published, uploaded, posted,
* transmitted, or distributed in any way without Intel's prior express written
* permission. No license under any patent, copyright or other intellectual
* property rights in the Material is granted to or conferred upon you, either
* expressly, by implication, inducement, estoppel or otherwise. Any license
* under such intellectual property rights must be express and approved by Intel
* in writing.
* =============================================================================
*
* DGESVD Example.
* ==============
*
* Program computes the singular value decomposition of a general
* rectangular matrix A:
*
* 8.79 9.93 9.83 5.45 3.16
* 6.11 6.91 5.04 -0.27 7.98
* -9.15 -7.93 4.86 4.85 3.01
* 9.57 1.64 8.83 0.74 5.80
* -3.49 4.02 9.80 10.00 4.27
* 9.84 0.15 -8.99 -6.02 -5.31
*
* Description.
* ============
*
* The routine computes the singular value decomposition (SVD) of a real
* m-by-n matrix A, optionally computing the left and/or right singular
* vectors. The SVD is written as
*
* A = U*SIGMA*VT
*
* where SIGMA is an m-by-n matrix which is zero except for its min(m,n)
* diagonal elements, U is an m-by-m orthogonal matrix and VT (V transposed)
* is an n-by-n orthogonal matrix. The diagonal elements of SIGMA
* are the singular values of A; they are real and non-negative, and are
* returned in descending order. The first min(m, n) columns of U and V are
* the left and right singular vectors of A.
*
* Note that the routine returns VT, not V.
*
* Example Program Results.
* ========================
*
* DGESVD Example Program Results
*
* Singular values
* 27.47 22.64 8.56 5.99 2.01
*
* Left singular vectors (stored columnwise)
* -0.59 0.26 0.36 0.31 0.23
* -0.40 0.24 -0.22 -0.75 -0.36
* -0.03 -0.60 -0.45 0.23 -0.31
* -0.43 0.24 -0.69 0.33 0.16
* -0.47 -0.35 0.39 0.16 -0.52
* 0.29 0.58 -0.02 0.38 -0.65
*
* Right singular vectors (stored rowwise)
* -0.25 -0.40 -0.69 -0.37 -0.41
* 0.81 0.36 -0.25 -0.37 -0.10
* -0.26 0.70 -0.22 0.39 -0.49
* 0.40 -0.45 0.25 0.43 -0.62
* -0.22 0.14 0.59 -0.63 -0.44
* =============================================================================
*
* .. Parameters ..INTEGER M, NPARAMETER ( M = 6, N = 5 )INTEGER LDA, LDU, LDVTPARAMETER ( LDA = M, LDU = M, LDVT = N )INTEGER LWMAXPARAMETER ( LWMAX = 1000 )
*
* .. Local Scalars ..INTEGER INFO, LWORK
*
* .. Local Arrays ..DOUBLE PRECISION A( LDA, N ), U( LDU, M ), VT( LDVT, N ), S( N ),$ WORK( LWMAX )DATA A/$ 8.79, 6.11,-9.15, 9.57,-3.49, 9.84,$ 9.93, 6.91,-7.93, 1.64, 4.02, 0.15,$ 9.83, 5.04, 4.86, 8.83, 9.80,-8.99,$ 5.45,-0.27, 4.85, 0.74,10.00,-6.02,$ 3.16, 7.98, 3.01, 5.80, 4.27,-5.31$ /
*
* .. External Subroutines ..EXTERNAL DGESVDEXTERNAL PRINT_MATRIX
*
* .. Intrinsic Functions ..INTRINSIC INT, MIN
*
* .. Executable Statements ..WRITE(*,*)'DGESVD Example Program Results'
*
* Query the optimal workspace.
*LWORK = -1CALL DGESVD( 'All', 'All', M, N, A, LDA, S, U, LDU, VT, LDVT,$ WORK, LWORK, INFO )LWORK = MIN( LWMAX, INT( WORK( 1 ) ) )
*
* Compute SVD.
*CALL DGESVD( 'All', 'All', M, N, A, LDA, S, U, LDU, VT, LDVT,$ WORK, LWORK, INFO )
*
* Check for convergence.
*IF( INFO.GT.0 ) THENWRITE(*,*)'The algorithm computing SVD failed to converge.'STOPEND IF
*
* Print singular values.
*CALL PRINT_MATRIX( 'Singular values', 1, N, S, 1 )
*
* Print left singular vectors.
*CALL PRINT_MATRIX( 'Left singular vectors (stored columnwise)',$ M, N, U, LDU )
*
* Print right singular vectors.
*CALL PRINT_MATRIX( 'Right singular vectors (stored rowwise)',$ N, N, VT, LDVT )STOPEND
*
* End of DGESVD Example.
*
* =============================================================================
*
* Auxiliary routine: printing a matrix.
*SUBROUTINE PRINT_MATRIX( DESC, M, N, A, LDA )CHARACTER*(*) DESCINTEGER M, N, LDADOUBLE PRECISION A( LDA, * )
*INTEGER I, J
*WRITE(*,*)WRITE(*,*) DESCDO I = 1, MWRITE(*,9998) ( A( I, J ), J = 1, N )END DO
*9998 FORMAT( 11(:,1X,F6.2) )RETURNENDMakefile
svd_dgesve_ex: svd_dgesve_ex.fgfortran -g $< ../lapack-3.11/liblapack.a ../lapack-3.11/librefblas.a -o $@运行
