ergo
template_blas_ger.h
Go to the documentation of this file.
1/* Ergo, version 3.8.2, a program for linear scaling electronic structure
2 * calculations.
3 * Copyright (C) 2023 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek,
4 * and Anastasia Kruchinina.
5 *
6 * This program is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * (at your option) any later version.
10 *
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 *
19 * Primary academic reference:
20 * Ergo: An open-source program for linear-scaling electronic structure
21 * calculations,
22 * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia
23 * Kruchinina,
24 * SoftwareX 7, 107 (2018),
25 * <http://dx.doi.org/10.1016/j.softx.2018.03.005>
26 *
27 * For further information about Ergo, see <http://www.ergoscf.org>.
28 */
29
30 /* This file belongs to the template_lapack part of the Ergo source
31 * code. The source files in the template_lapack directory are modified
32 * versions of files originally distributed as CLAPACK, see the
33 * Copyright/license notice in the file template_lapack/COPYING.
34 */
35
36
37#ifndef TEMPLATE_BLAS_GER_HEADER
38#define TEMPLATE_BLAS_GER_HEADER
39
40
41template<class Treal>
42int template_blas_ger(const integer *m, const integer *n, const Treal *alpha,
43 const Treal *x, const integer *incx, const Treal *y, const integer *incy,
44 Treal *a, const integer *lda)
45{
46 /* System generated locals */
47 integer a_dim1, a_offset, i__1, i__2;
48 /* Local variables */
49 integer info;
50 Treal temp;
51 integer i__, j, ix, jy, kx;
52#define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1]
53/* Purpose
54 =======
55 DGER performs the rank 1 operation
56 A := alpha*x*y' + A,
57 where alpha is a scalar, x is an m element vector, y is an n element
58 vector and A is an m by n matrix.
59 Parameters
60 ==========
61 M - INTEGER.
62 On entry, M specifies the number of rows of the matrix A.
63 M must be at least zero.
64 Unchanged on exit.
65 N - INTEGER.
66 On entry, N specifies the number of columns of the matrix A.
67 N must be at least zero.
68 Unchanged on exit.
69 ALPHA - DOUBLE PRECISION.
70 On entry, ALPHA specifies the scalar alpha.
71 Unchanged on exit.
72 X - DOUBLE PRECISION array of dimension at least
73 ( 1 + ( m - 1 )*abs( INCX ) ).
74 Before entry, the incremented array X must contain the m
75 element vector x.
76 Unchanged on exit.
77 INCX - INTEGER.
78 On entry, INCX specifies the increment for the elements of
79 X. INCX must not be zero.
80 Unchanged on exit.
81 Y - DOUBLE PRECISION array of dimension at least
82 ( 1 + ( n - 1 )*abs( INCY ) ).
83 Before entry, the incremented array Y must contain the n
84 element vector y.
85 Unchanged on exit.
86 INCY - INTEGER.
87 On entry, INCY specifies the increment for the elements of
88 Y. INCY must not be zero.
89 Unchanged on exit.
90 A - DOUBLE PRECISION array of DIMENSION ( LDA, n ).
91 Before entry, the leading m by n part of the array A must
92 contain the matrix of coefficients. On exit, A is
93 overwritten by the updated matrix.
94 LDA - INTEGER.
95 On entry, LDA specifies the first dimension of A as declared
96 in the calling (sub) program. LDA must be at least
97 max( 1, m ).
98 Unchanged on exit.
99 Level 2 Blas routine.
100 -- Written on 22-October-1986.
101 Jack Dongarra, Argonne National Lab.
102 Jeremy Du Croz, Nag Central Office.
103 Sven Hammarling, Nag Central Office.
104 Richard Hanson, Sandia National Labs.
105 Test the input parameters.
106 Parameter adjustments */
107 --x;
108 --y;
109 a_dim1 = *lda;
110 a_offset = 1 + a_dim1 * 1;
111 a -= a_offset;
112 /* Function Body */
113 info = 0;
114 if (*m < 0) {
115 info = 1;
116 } else if (*n < 0) {
117 info = 2;
118 } else if (*incx == 0) {
119 info = 5;
120 } else if (*incy == 0) {
121 info = 7;
122 } else if (*lda < maxMACRO(1,*m)) {
123 info = 9;
124 }
125 if (info != 0) {
126 template_blas_erbla("GER ", &info);
127 return 0;
128 }
129/* Quick return if possible. */
130 if (*m == 0 || *n == 0 || *alpha == 0.) {
131 return 0;
132 }
133/* Start the operations. In this version the elements of A are
134 accessed sequentially with one pass through A. */
135 if (*incy > 0) {
136 jy = 1;
137 } else {
138 jy = 1 - (*n - 1) * *incy;
139 }
140 if (*incx == 1) {
141 i__1 = *n;
142 for (j = 1; j <= i__1; ++j) {
143 if (y[jy] != 0.) {
144 temp = *alpha * y[jy];
145 i__2 = *m;
146 for (i__ = 1; i__ <= i__2; ++i__) {
147 a_ref(i__, j) = a_ref(i__, j) + x[i__] * temp;
148/* L10: */
149 }
150 }
151 jy += *incy;
152/* L20: */
153 }
154 } else {
155 if (*incx > 0) {
156 kx = 1;
157 } else {
158 kx = 1 - (*m - 1) * *incx;
159 }
160 i__1 = *n;
161 for (j = 1; j <= i__1; ++j) {
162 if (y[jy] != 0.) {
163 temp = *alpha * y[jy];
164 ix = kx;
165 i__2 = *m;
166 for (i__ = 1; i__ <= i__2; ++i__) {
167 a_ref(i__, j) = a_ref(i__, j) + x[ix] * temp;
168 ix += *incx;
169/* L30: */
170 }
171 }
172 jy += *incy;
173/* L40: */
174 }
175 }
176 return 0;
177/* End of DGER . */
178} /* dger_ */
179#undef a_ref
180
181#endif
int template_blas_erbla(const char *srname, integer *info)
Definition template_blas_common.cc:146
int integer
Definition template_blas_common.h:40
#define maxMACRO(a, b)
Definition template_blas_common.h:45
int template_blas_ger(const integer *m, const integer *n, const Treal *alpha, const Treal *x, const integer *incx, const Treal *y, const integer *incy, Treal *a, const integer *lda)
Definition template_blas_ger.h:42
#define a_ref(a_1, a_2)