Télécharger dgemm.eso

Retour à la liste

Numérotation des lignes :

  1. C DGEMM SOURCE BP208322 20/09/18 21:15:51 10718
  2. *> \brief \b DGEMM
  3. *
  4. * =========== DOCUMENTATION ===========
  5. *
  6. * Online html documentation available at
  7. * http://www.netlib.org/lapack/explore-html/
  8. *
  9. * Definition:
  10. * ===========
  11. *
  12. * SUBROUTINE DGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
  13. *
  14. * .. Scalar Arguments ..
  15. * REAL ALPHA,BETA
  16. * INTEGER K,LDA,LDB,LDC,M,N
  17. * CHARACTER TRANSA,TRANSB
  18. * ..
  19. * .. Array Arguments ..
  20. * REAL A(LDA,*),B(LDB,*),C(LDC,*)
  21. * ..
  22. *
  23. *
  24. *> \par Purpose:
  25. * =============
  26. *>
  27. *> \verbatim
  28. *>
  29. *> DGEMM performs one of the matrix-matrix operations
  30. *>
  31. *> C := alpha*op( A )*op( B ) + beta*C,
  32. *>
  33. *> where op( X ) is one of
  34. *>
  35. *> op( X ) = X or op( X ) = X**T,
  36. *>
  37. *> alpha and beta are scalars, and A, B and C are matrices, with op( A )
  38. *> an m by k matrix, op( B ) a k by n matrix and C an m by n matrix.
  39. *> \endverbatim
  40. *
  41. * Arguments:
  42. * ==========
  43. *
  44. *> \param[in] TRANSA
  45. *> \verbatim
  46. *> TRANSA is CHARACTER*1
  47. *> On entry, TRANSA specifies the form of op( A ) to be used in
  48. *> the matrix multiplication as follows:
  49. *>
  50. *> TRANSA = 'N' or 'n', op( A ) = A.
  51. *>
  52. *> TRANSA = 'T' or 't', op( A ) = A**T.
  53. *>
  54. *> TRANSA = 'C' or 'c', op( A ) = A**T.
  55. *> \endverbatim
  56. *>
  57. *> \param[in] TRANSB
  58. *> \verbatim
  59. *> TRANSB is CHARACTER*1
  60. *> On entry, TRANSB specifies the form of op( B ) to be used in
  61. *> the matrix multiplication as follows:
  62. *>
  63. *> TRANSB = 'N' or 'n', op( B ) = B.
  64. *>
  65. *> TRANSB = 'T' or 't', op( B ) = B**T.
  66. *>
  67. *> TRANSB = 'C' or 'c', op( B ) = B**T.
  68. *> \endverbatim
  69. *>
  70. *> \param[in] M
  71. *> \verbatim
  72. *> M is INTEGER
  73. *> On entry, M specifies the number of rows of the matrix
  74. *> op( A ) and of the matrix C. M must be at least zero.
  75. *> \endverbatim
  76. *>
  77. *> \param[in] N
  78. *> \verbatim
  79. *> N is INTEGER
  80. *> On entry, N specifies the number of columns of the matrix
  81. *> op( B ) and the number of columns of the matrix C. N must be
  82. *> at least zero.
  83. *> \endverbatim
  84. *>
  85. *> \param[in] K
  86. *> \verbatim
  87. *> K is INTEGER
  88. *> On entry, K specifies the number of columns of the matrix
  89. *> op( A ) and the number of rows of the matrix op( B ). K must
  90. *> be at least zero.
  91. *> \endverbatim
  92. *>
  93. *> \param[in] ALPHA
  94. *> \verbatim
  95. *> ALPHA is REAL
  96. *> On entry, ALPHA specifies the scalar alpha.
  97. *> \endverbatim
  98. *>
  99. *> \param[in] A
  100. *> \verbatim
  101. *> A is REAL array, dimension ( LDA, ka ), where ka is
  102. *> k when TRANSA = 'N' or 'n', and is m otherwise.
  103. *> Before entry with TRANSA = 'N' or 'n', the leading m by k
  104. *> part of the array A must contain the matrix A, otherwise
  105. *> the leading k by m part of the array A must contain the
  106. *> matrix A.
  107. *> \endverbatim
  108. *>
  109. *> \param[in] LDA
  110. *> \verbatim
  111. *> LDA is INTEGER
  112. *> On entry, LDA specifies the first dimension of A as declared
  113. *> in the calling (sub) program. When TRANSA = 'N' or 'n' then
  114. *> LDA must be at least max( 1, m ), otherwise LDA must be at
  115. *> least max( 1, k ).
  116. *> \endverbatim
  117. *>
  118. *> \param[in] B
  119. *> \verbatim
  120. *> B is REAL array, dimension ( LDB, kb ), where kb is
  121. *> n when TRANSB = 'N' or 'n', and is k otherwise.
  122. *> Before entry with TRANSB = 'N' or 'n', the leading k by n
  123. *> part of the array B must contain the matrix B, otherwise
  124. *> the leading n by k part of the array B must contain the
  125. *> matrix B.
  126. *> \endverbatim
  127. *>
  128. *> \param[in] LDB
  129. *> \verbatim
  130. *> LDB is INTEGER
  131. *> On entry, LDB specifies the first dimension of B as declared
  132. *> in the calling (sub) program. When TRANSB = 'N' or 'n' then
  133. *> LDB must be at least max( 1, k ), otherwise LDB must be at
  134. *> least max( 1, n ).
  135. *> \endverbatim
  136. *>
  137. *> \param[in] BETA
  138. *> \verbatim
  139. *> BETA is REAL
  140. *> On entry, BETA specifies the scalar beta. When BETA is
  141. *> supplied as zero then C need not be set on input.
  142. *> \endverbatim
  143. *>
  144. *> \param[in,out] C
  145. *> \verbatim
  146. *> C is REAL array, dimension ( LDC, N )
  147. *> Before entry, the leading m by n part of the array C must
  148. *> contain the matrix C, except when beta is zero, in which
  149. *> case C need not be set on entry.
  150. *> On exit, the array C is overwritten by the m by n matrix
  151. *> ( alpha*op( A )*op( B ) + beta*C ).
  152. *> \endverbatim
  153. *>
  154. *> \param[in] LDC
  155. *> \verbatim
  156. *> LDC is INTEGER
  157. *> On entry, LDC specifies the first dimension of C as declared
  158. *> in the calling (sub) program. LDC must be at least
  159. *> max( 1, m ).
  160. *> \endverbatim
  161. *
  162. * Authors:
  163. * ========
  164. *
  165. *> \author Univ. of Tennessee
  166. *> \author Univ. of California Berkeley
  167. *> \author Univ. of Colorado Denver
  168. *> \author NAG Ltd.
  169. *
  170. *> \date December 2016
  171. *
  172. *> \ingroup single_blas_level3
  173. *
  174. *> \par Further Details:
  175. * =====================
  176. *>
  177. *> \verbatim
  178. *>
  179. *> Level 3 Blas routine.
  180. *>
  181. *> -- Written on 8-February-1989.
  182. *> Jack Dongarra, Argonne National Laboratory.
  183. *> Iain Duff, AERE Harwell.
  184. *> Jeremy Du Croz, Numerical Algorithms Group Ltd.
  185. *> Sven Hammarling, Numerical Algorithms Group Ltd.
  186. *> \endverbatim
  187. *>
  188. * =====================================================================
  189. SUBROUTINE DGEMM(TRANSA,TRANSB,M,N,K,ALPHA,A,LDA,B,LDB,BETA,C,LDC)
  190. *
  191. * -- Reference BLAS level3 routine (version 3.7.0) --
  192. * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
  193. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  194. * December 2016
  195.  
  196. IMPLICIT INTEGER(I-N)
  197. IMPLICIT REAL*8(A-H,O-Z)
  198. *
  199. * .. Scalar Arguments ..
  200. REAL*8 ALPHA,BETA
  201. INTEGER K,LDA,LDB,LDC,M,N
  202. CHARACTER TRANSA,TRANSB
  203. * ..
  204. * .. Array Arguments ..
  205. REAL*8 A(LDA,*),B(LDB,*),C(LDC,*)
  206. * ..
  207. *
  208. * =====================================================================
  209. * ..
  210. * .. External Subroutines ..
  211. * EXTERNAL XERBLA
  212. * .. Parameters ..
  213. REAL*8 ONE, ZERO
  214. PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
  215. * .. Local Scalars ..
  216. REAL*8 TEMP
  217. INTEGER I,INFO,J,L,NCOLA,NROWA,NROWB
  218. LOGICAL NOTA,NOTB
  219. *
  220. * Set NOTA and NOTB as true if A and B respectively are not
  221. * transposed and set NROWA, NCOLA and NROWB as the number of rows
  222. * and columns of A and the number of rows of B respectively.
  223. *
  224. NOTA = (TRANSA.EQ.'N')
  225. NOTB = (TRANSB.EQ.'N')
  226. IF (NOTA) THEN
  227. NROWA = M
  228. NCOLA = K
  229. ELSE
  230. NROWA = K
  231. NCOLA = M
  232. END IF
  233. IF (NOTB) THEN
  234. NROWB = K
  235. ELSE
  236. NROWB = N
  237. END IF
  238. *
  239. * Test the input parameters.
  240. *
  241. INFO = 0
  242. IF ((.NOT.NOTA) .AND. (.NOT.(TRANSA.EQ.'C')) .AND.
  243. + (.NOT.(TRANSA.EQ.'T'))) THEN
  244. INFO = 1
  245. ELSE IF ((.NOT.NOTB) .AND. (.NOT.(TRANSB.EQ.'C')) .AND.
  246. + (.NOT.(TRANSB.EQ.'T'))) THEN
  247. INFO = 2
  248. ELSE IF (M.LT.0) THEN
  249. INFO = 3
  250. ELSE IF (N.LT.0) THEN
  251. INFO = 4
  252. ELSE IF (K.LT.0) THEN
  253. INFO = 5
  254. ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
  255. INFO = 8
  256. ELSE IF (LDB.LT.MAX(1,NROWB)) THEN
  257. INFO = 10
  258. ELSE IF (LDC.LT.MAX(1,M)) THEN
  259. INFO = 13
  260. END IF
  261. IF (INFO.NE.0) THEN
  262. CALL XERBLA('SGEMM ',INFO)
  263. RETURN
  264. END IF
  265. *
  266. * Quick return if possible.
  267. *
  268. IF ((M.EQ.0) .OR. (N.EQ.0) .OR.
  269. + (((ALPHA.EQ.ZERO).OR. (K.EQ.0)).AND. (BETA.EQ.ONE))) RETURN
  270. *
  271. * And if alpha.eq.zero.
  272. *
  273. IF (ALPHA.EQ.ZERO) THEN
  274. IF (BETA.EQ.ZERO) THEN
  275. DO 20 J = 1,N
  276. DO 10 I = 1,M
  277. C(I,J) = ZERO
  278. 10 CONTINUE
  279. 20 CONTINUE
  280. ELSE
  281. DO 40 J = 1,N
  282. DO 30 I = 1,M
  283. C(I,J) = BETA*C(I,J)
  284. 30 CONTINUE
  285. 40 CONTINUE
  286. END IF
  287. RETURN
  288. END IF
  289. *
  290. * Start the operations.
  291. *
  292. IF (NOTB) THEN
  293. IF (NOTA) THEN
  294. *
  295. * Form C := alpha*A*B + beta*C.
  296. *
  297. DO 90 J = 1,N
  298. IF (BETA.EQ.ZERO) THEN
  299. DO 50 I = 1,M
  300. C(I,J) = ZERO
  301. 50 CONTINUE
  302. ELSE IF (BETA.NE.ONE) THEN
  303. DO 60 I = 1,M
  304. C(I,J) = BETA*C(I,J)
  305. 60 CONTINUE
  306. END IF
  307. DO 80 L = 1,K
  308. TEMP = ALPHA*B(L,J)
  309. DO 70 I = 1,M
  310. C(I,J) = C(I,J) + TEMP*A(I,L)
  311. 70 CONTINUE
  312. 80 CONTINUE
  313. 90 CONTINUE
  314. ELSE
  315. *
  316. * Form C := alpha*A**T*B + beta*C
  317. *
  318. DO 120 J = 1,N
  319. DO 110 I = 1,M
  320. TEMP = ZERO
  321. DO 100 L = 1,K
  322. TEMP = TEMP + A(L,I)*B(L,J)
  323. 100 CONTINUE
  324. IF (BETA.EQ.ZERO) THEN
  325. C(I,J) = ALPHA*TEMP
  326. ELSE
  327. C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  328. END IF
  329. 110 CONTINUE
  330. 120 CONTINUE
  331. END IF
  332. ELSE
  333. IF (NOTA) THEN
  334. *
  335. * Form C := alpha*A*B**T + beta*C
  336. *
  337. DO 170 J = 1,N
  338. IF (BETA.EQ.ZERO) THEN
  339. DO 130 I = 1,M
  340. C(I,J) = ZERO
  341. 130 CONTINUE
  342. ELSE IF (BETA.NE.ONE) THEN
  343. DO 140 I = 1,M
  344. C(I,J) = BETA*C(I,J)
  345. 140 CONTINUE
  346. END IF
  347. DO 160 L = 1,K
  348. TEMP = ALPHA*B(J,L)
  349. DO 150 I = 1,M
  350. C(I,J) = C(I,J) + TEMP*A(I,L)
  351. 150 CONTINUE
  352. 160 CONTINUE
  353. 170 CONTINUE
  354. ELSE
  355. *
  356. * Form C := alpha*A**T*B**T + beta*C
  357. *
  358. DO 200 J = 1,N
  359. DO 190 I = 1,M
  360. TEMP = ZERO
  361. DO 180 L = 1,K
  362. TEMP = TEMP + A(L,I)*B(J,L)
  363. 180 CONTINUE
  364. IF (BETA.EQ.ZERO) THEN
  365. C(I,J) = ALPHA*TEMP
  366. ELSE
  367. C(I,J) = ALPHA*TEMP + BETA*C(I,J)
  368. END IF
  369. 190 CONTINUE
  370. 200 CONTINUE
  371. END IF
  372. END IF
  373. *
  374. RETURN
  375. *
  376. * End of GEMM .
  377. *
  378. END
  379.  
  380.  
  381.  
  382.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales