赞
踩
- SUBROUTINE INIT_MATRIX(A, m, n, lda)
- DOUBLE PRECISION A(*)
- CALL SRAND(2024)
- DO i=1, m
- DO j=1, n
- A(i + lda*(j-1)) = RAND() + RAND()
- C WRITE(*, '(F8.4)') A(i)
- END DO
- END DO
- END
-
-
- SUBROUTINE PRINT_MATRIX(A, m, n, lda)
- DOUBLE PRECISION A(*)
- DO i=1, m
- DO 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 DO
- PRINT *,''
- END DO
- END
-
-
- PROGRAM MATRIX_EX
- INTEGER m, n, lda
- DOUBLE PRECISION A(20)
-
- lda = 4
- m = 4
- n = 5
- CALL INIT_MATRIX(A, m, n, lda)
-
- PRINT *, "A ="
- CALL PRINT_MATRIX(A, m, n, lda)
-
- PRINT *, "A(11) ="
- PRINT *, A(11)
-
- END
- matrix_ex: matrix_ex.f
- gfortran -g $< -o $@
-
-
- .PHONY: clean
- clean:
- -rm matrix_ex
- make
-
- ./matrix_ex
效果:
- * 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, N
- PARAMETER ( M = 6, N = 5 )
- INTEGER LDA, LDU, LDVT
- PARAMETER ( LDA = M, LDU = M, LDVT = N )
- INTEGER LWMAX
- PARAMETER ( 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 DGESVD
- EXTERNAL PRINT_MATRIX
- *
- * .. Intrinsic Functions ..
- INTRINSIC INT, MIN
- *
- * .. Executable Statements ..
- WRITE(*,*)'DGESVD Example Program Results'
- *
- * Query the optimal workspace.
- *
- LWORK = -1
- CALL 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 ) THEN
- WRITE(*,*)'The algorithm computing SVD failed to converge.'
- STOP
- END 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 )
- STOP
- END
- *
- * End of DGESVD Example.
- *
- * =============================================================================
- *
- * Auxiliary routine: printing a matrix.
- *
- SUBROUTINE PRINT_MATRIX( DESC, M, N, A, LDA )
- CHARACTER*(*) DESC
- INTEGER M, N, LDA
- DOUBLE PRECISION A( LDA, * )
- *
- INTEGER I, J
- *
- WRITE(*,*)
- WRITE(*,*) DESC
- DO I = 1, M
- WRITE(*,9998) ( A( I, J ), J = 1, N )
- END DO
- *
- 9998 FORMAT( 11(:,1X,F6.2) )
- RETURN
- END
Makefile
- svd_dgesve_ex: svd_dgesve_ex.f
- gfortran -g $< ../lapack-3.11/liblapack.a ../lapack-3.11/librefblas.a -o $@
运行
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。