level2_real_impl.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2009-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #include "common.h"
11 
12 // y = alpha*A*x + beta*y
13 int EIGEN_BLAS_FUNC(symv) (const char *uplo, const int *n, const RealScalar *palpha, const RealScalar *pa, const int *lda,
14  const RealScalar *px, const int *incx, const RealScalar *pbeta, RealScalar *py, const int *incy)
15 {
16  typedef void (*functype)(int, const Scalar*, int, const Scalar*, Scalar*, Scalar);
17  static const functype func[2] = {
18  // array index: UP
19  (internal::selfadjoint_matrix_vector_product<Scalar,int,ColMajor,Upper,false,false>::run),
20  // array index: LO
21  (internal::selfadjoint_matrix_vector_product<Scalar,int,ColMajor,Lower,false,false>::run),
22  };
23 
24  const Scalar* a = reinterpret_cast<const Scalar*>(pa);
25  const Scalar* x = reinterpret_cast<const Scalar*>(px);
26  Scalar* y = reinterpret_cast<Scalar*>(py);
27  Scalar alpha = *reinterpret_cast<const Scalar*>(palpha);
28  Scalar beta = *reinterpret_cast<const Scalar*>(pbeta);
29 
30  // check arguments
31  int info = 0;
32  if(UPLO(*uplo)==INVALID) info = 1;
33  else if(*n<0) info = 2;
34  else if(*lda<std::max(1,*n)) info = 5;
35  else if(*incx==0) info = 7;
36  else if(*incy==0) info = 10;
37  if(info)
38  return xerbla_(SCALAR_SUFFIX_UP"SYMV ",&info,6);
39 
40  if(*n==0)
41  return 0;
42 
43  const Scalar* actual_x = get_compact_vector(x,*n,*incx);
44  Scalar* actual_y = get_compact_vector(y,*n,*incy);
45 
46  if(beta!=Scalar(1))
47  {
48  if(beta==Scalar(0)) make_vector(actual_y, *n).setZero();
49  else make_vector(actual_y, *n) *= beta;
50  }
51 
52  int code = UPLO(*uplo);
53  if(code>=2 || func[code]==0)
54  return 0;
55 
56  func[code](*n, a, *lda, actual_x, actual_y, alpha);
57 
58  if(actual_x!=x) delete[] actual_x;
59  if(actual_y!=y) delete[] copy_back(actual_y,y,*n,*incy);
60 
61  return 1;
62 }
63 
64 // C := alpha*x*x' + C
65 int EIGEN_BLAS_FUNC(syr)(const char *uplo, const int *n, const RealScalar *palpha, const RealScalar *px, const int *incx, RealScalar *pc, const int *ldc)
66 {
67 
68  typedef void (*functype)(int, Scalar*, int, const Scalar*, const Scalar*, const Scalar&);
69  static const functype func[2] = {
70  // array index: UP
71  (selfadjoint_rank1_update<Scalar,int,ColMajor,Upper,false,Conj>::run),
72  // array index: LO
73  (selfadjoint_rank1_update<Scalar,int,ColMajor,Lower,false,Conj>::run),
74  };
75 
76  const Scalar* x = reinterpret_cast<const Scalar*>(px);
77  Scalar* c = reinterpret_cast<Scalar*>(pc);
78  Scalar alpha = *reinterpret_cast<const Scalar*>(palpha);
79 
80  int info = 0;
81  if(UPLO(*uplo)==INVALID) info = 1;
82  else if(*n<0) info = 2;
83  else if(*incx==0) info = 5;
84  else if(*ldc<std::max(1,*n)) info = 7;
85  if(info)
86  return xerbla_(SCALAR_SUFFIX_UP"SYR ",&info,6);
87 
88  if(*n==0 || alpha==Scalar(0)) return 1;
89 
90  // if the increment is not 1, let's copy it to a temporary vector to enable vectorization
91  const Scalar* x_cpy = get_compact_vector(x,*n,*incx);
92 
93  int code = UPLO(*uplo);
94  if(code>=2 || func[code]==0)
95  return 0;
96 
97  func[code](*n, c, *ldc, x_cpy, x_cpy, alpha);
98 
99  if(x_cpy!=x) delete[] x_cpy;
100 
101  return 1;
102 }
103 
104 // C := alpha*x*y' + alpha*y*x' + C
105 int EIGEN_BLAS_FUNC(syr2)(const char *uplo, const int *n, const RealScalar *palpha, const RealScalar *px, const int *incx, const RealScalar *py, const int *incy, RealScalar *pc, const int *ldc)
106 {
107  typedef void (*functype)(int, Scalar*, int, const Scalar*, const Scalar*, Scalar);
108  static const functype func[2] = {
109  // array index: UP
110  (internal::rank2_update_selector<Scalar,int,Upper>::run),
111  // array index: LO
112  (internal::rank2_update_selector<Scalar,int,Lower>::run),
113  };
114 
115  const Scalar* x = reinterpret_cast<const Scalar*>(px);
116  const Scalar* y = reinterpret_cast<const Scalar*>(py);
117  Scalar* c = reinterpret_cast<Scalar*>(pc);
118  Scalar alpha = *reinterpret_cast<const Scalar*>(palpha);
119 
120  int info = 0;
121  if(UPLO(*uplo)==INVALID) info = 1;
122  else if(*n<0) info = 2;
123  else if(*incx==0) info = 5;
124  else if(*incy==0) info = 7;
125  else if(*ldc<std::max(1,*n)) info = 9;
126  if(info)
127  return xerbla_(SCALAR_SUFFIX_UP"SYR2 ",&info,6);
128 
129  if(alpha==Scalar(0))
130  return 1;
131 
132  const Scalar* x_cpy = get_compact_vector(x,*n,*incx);
133  const Scalar* y_cpy = get_compact_vector(y,*n,*incy);
134 
135  int code = UPLO(*uplo);
136  if(code>=2 || func[code]==0)
137  return 0;
138 
139  func[code](*n, c, *ldc, x_cpy, y_cpy, alpha);
140 
141  if(x_cpy!=x) delete[] x_cpy;
142  if(y_cpy!=y) delete[] y_cpy;
143 
144 // int code = UPLO(*uplo);
145 // if(code>=2 || func[code]==0)
146 // return 0;
147 
148 // func[code](*n, a, *inca, b, *incb, c, *ldc, alpha);
149  return 1;
150 }
151 
159 // int EIGEN_BLAS_FUNC(sbmv)( char *uplo, int *n, int *k, RealScalar *alpha, RealScalar *a, int *lda,
160 // RealScalar *x, int *incx, RealScalar *beta, RealScalar *y, int *incy)
161 // {
162 // return 1;
163 // }
164 
165 
174 // int EIGEN_BLAS_FUNC(spmv)(char *uplo, int *n, RealScalar *alpha, RealScalar *ap, RealScalar *x, int *incx, RealScalar *beta, RealScalar *y, int *incy)
175 // {
176 // return 1;
177 // }
178 
186 int EIGEN_BLAS_FUNC(spr)(char *uplo, int *n, Scalar *palpha, Scalar *px, int *incx, Scalar *pap)
187 {
188  typedef void (*functype)(int, Scalar*, const Scalar*, Scalar);
189  static const functype func[2] = {
190  // array index: UP
191  (internal::selfadjoint_packed_rank1_update<Scalar,int,ColMajor,Upper,false,false>::run),
192  // array index: LO
193  (internal::selfadjoint_packed_rank1_update<Scalar,int,ColMajor,Lower,false,false>::run),
194  };
195 
196  Scalar* x = reinterpret_cast<Scalar*>(px);
197  Scalar* ap = reinterpret_cast<Scalar*>(pap);
198  Scalar alpha = *reinterpret_cast<Scalar*>(palpha);
199 
200  int info = 0;
201  if(UPLO(*uplo)==INVALID) info = 1;
202  else if(*n<0) info = 2;
203  else if(*incx==0) info = 5;
204  if(info)
205  return xerbla_(SCALAR_SUFFIX_UP"SPR ",&info,6);
206 
207  if(alpha==Scalar(0))
208  return 1;
209 
210  Scalar* x_cpy = get_compact_vector(x, *n, *incx);
211 
212  int code = UPLO(*uplo);
213  if(code>=2 || func[code]==0)
214  return 0;
215 
216  func[code](*n, ap, x_cpy, alpha);
217 
218  if(x_cpy!=x) delete[] x_cpy;
219 
220  return 1;
221 }
222 
230 int EIGEN_BLAS_FUNC(spr2)(char *uplo, int *n, RealScalar *palpha, RealScalar *px, int *incx, RealScalar *py, int *incy, RealScalar *pap)
231 {
232  typedef void (*functype)(int, Scalar*, const Scalar*, const Scalar*, Scalar);
233  static const functype func[2] = {
234  // array index: UP
235  (internal::packed_rank2_update_selector<Scalar,int,Upper>::run),
236  // array index: LO
237  (internal::packed_rank2_update_selector<Scalar,int,Lower>::run),
238  };
239 
240  Scalar* x = reinterpret_cast<Scalar*>(px);
241  Scalar* y = reinterpret_cast<Scalar*>(py);
242  Scalar* ap = reinterpret_cast<Scalar*>(pap);
243  Scalar alpha = *reinterpret_cast<Scalar*>(palpha);
244 
245  int info = 0;
246  if(UPLO(*uplo)==INVALID) info = 1;
247  else if(*n<0) info = 2;
248  else if(*incx==0) info = 5;
249  else if(*incy==0) info = 7;
250  if(info)
251  return xerbla_(SCALAR_SUFFIX_UP"SPR2 ",&info,6);
252 
253  if(alpha==Scalar(0))
254  return 1;
255 
256  Scalar* x_cpy = get_compact_vector(x, *n, *incx);
257  Scalar* y_cpy = get_compact_vector(y, *n, *incy);
258 
259  int code = UPLO(*uplo);
260  if(code>=2 || func[code]==0)
261  return 0;
262 
263  func[code](*n, ap, x_cpy, y_cpy, alpha);
264 
265  if(x_cpy!=x) delete[] x_cpy;
266  if(y_cpy!=y) delete[] y_cpy;
267 
268  return 1;
269 }
270 
278 int EIGEN_BLAS_FUNC(ger)(int *m, int *n, Scalar *palpha, Scalar *px, int *incx, Scalar *py, int *incy, Scalar *pa, int *lda)
279 {
280  Scalar* x = reinterpret_cast<Scalar*>(px);
281  Scalar* y = reinterpret_cast<Scalar*>(py);
282  Scalar* a = reinterpret_cast<Scalar*>(pa);
283  Scalar alpha = *reinterpret_cast<Scalar*>(palpha);
284 
285  int info = 0;
286  if(*m<0) info = 1;
287  else if(*n<0) info = 2;
288  else if(*incx==0) info = 5;
289  else if(*incy==0) info = 7;
290  else if(*lda<std::max(1,*m)) info = 9;
291  if(info)
292  return xerbla_(SCALAR_SUFFIX_UP"GER ",&info,6);
293 
294  if(alpha==Scalar(0))
295  return 1;
296 
297  Scalar* x_cpy = get_compact_vector(x,*m,*incx);
298  Scalar* y_cpy = get_compact_vector(y,*n,*incy);
299 
300  internal::general_rank1_update<Scalar,int,ColMajor,false,false>::run(*m, *n, a, *lda, x_cpy, y_cpy, alpha);
301 
302  if(x_cpy!=x) delete[] x_cpy;
303  if(y_cpy!=y) delete[] y_cpy;
304 
305  return 1;
306 }
gtsam.examples.DogLegOptimizerExample.int
int
Definition: DogLegOptimizerExample.py:111
internal::selfadjoint_packed_rank1_update
Definition: PackedSelfadjointProduct.h:19
UPLO
#define UPLO(X)
Definition: gtsam/3rdparty/Eigen/blas/common.h:56
alpha
RealScalar alpha
Definition: level1_cplx_impl.h:147
c
Scalar Scalar * c
Definition: benchVecAdd.cpp:17
x
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
Definition: gnuplot_common_settings.hh:12
Eigen::selfadjoint_rank1_update
Definition: GeneralMatrixMatrixTriangular.h:16
make_vector
Map< Matrix< T, Dynamic, 1 >, 0, InnerStride< Dynamic > > make_vector(T *data, int size, int incr)
Definition: gtsam/3rdparty/Eigen/blas/common.h:123
get_compact_vector
T * get_compact_vector(T *x, int n, int incx)
Definition: gtsam/3rdparty/Eigen/blas/common.h:147
ger
int EIGEN_BLAS_FUNC() ger(int *m, int *n, Scalar *palpha, Scalar *px, int *incx, Scalar *py, int *incy, Scalar *pa, int *lda)
Definition: level2_real_impl.h:278
internal::rank2_update_selector
Definition: Rank2Update.h:19
beta
double beta(double a, double b)
Definition: beta.c:61
symv
int EIGEN_BLAS_FUNC() symv(const char *uplo, const int *n, const RealScalar *palpha, const RealScalar *pa, const int *lda, const RealScalar *px, const int *incx, const RealScalar *pbeta, RealScalar *py, const int *incy)
Definition: level2_real_impl.h:13
n
int n
Definition: BiCGSTAB_simple.cpp:1
common.h
incy
int RealScalar int RealScalar int * incy
Definition: level1_cplx_impl.h:119
syr
int EIGEN_BLAS_FUNC() syr(const char *uplo, const int *n, const RealScalar *palpha, const RealScalar *px, const int *incx, RealScalar *pc, const int *ldc)
Definition: level2_real_impl.h:65
EIGEN_BLAS_FUNC
#define EIGEN_BLAS_FUNC(X)
Definition: gtsam/3rdparty/Eigen/blas/common.h:173
incx
RealScalar RealScalar int * incx
Definition: level1_cplx_impl.h:29
info
else if n * info
Definition: 3rdparty/Eigen/lapack/cholesky.cpp:18
SCALAR_SUFFIX_UP
#define SCALAR_SUFFIX_UP
Definition: blas/complex_double.cpp:12
gtsam.examples.DogLegOptimizerExample.run
def run(args)
Definition: DogLegOptimizerExample.py:21
m
Matrix3f m
Definition: AngleAxis_mimic_euler.cpp:1
spr2
int EIGEN_BLAS_FUNC() spr2(char *uplo, int *n, RealScalar *palpha, RealScalar *px, int *incx, RealScalar *py, int *incy, RealScalar *pap)
Definition: level2_real_impl.h:230
pc
int RealScalar int RealScalar int RealScalar * pc
Definition: level1_cplx_impl.h:119
lda
* lda
Definition: eigenvalues.cpp:59
y
Scalar * y
Definition: level1_cplx_impl.h:124
RealScalar
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:47
a
ArrayXXi a
Definition: Array_initializer_list_23_cxx11.cpp:1
xerbla_
EIGEN_WEAK_LINKING int xerbla_(const char *msg, int *info, int)
Definition: xerbla.cpp:15
syr2
int EIGEN_BLAS_FUNC() syr2(const char *uplo, const int *n, const RealScalar *palpha, const RealScalar *px, const int *incx, const RealScalar *py, const int *incy, RealScalar *pc, const int *ldc)
Definition: level2_real_impl.h:105
internal::packed_rank2_update_selector
Definition: Rank2Update.h:37
copy_back
T * copy_back(T *x_cpy, T *x, int n, int incx)
Definition: gtsam/3rdparty/Eigen/blas/common.h:159
py
int RealScalar int RealScalar * py
Definition: level1_cplx_impl.h:119
spr
int EIGEN_BLAS_FUNC() spr(char *uplo, int *n, Scalar *palpha, Scalar *px, int *incx, Scalar *pap)
Definition: level2_real_impl.h:186
palpha
int RealScalar * palpha
Definition: level1_cplx_impl.h:142
func
Definition: benchGeometry.cpp:23
max
#define max(a, b)
Definition: datatypes.h:20
conf.code
code
Definition: wrap/pybind11/docs/conf.py:73
INVALID
#define INVALID
Definition: gtsam/3rdparty/Eigen/blas/common.h:45
Scalar
SCALAR Scalar
Definition: bench_gemm.cpp:46
px
RealScalar RealScalar * px
Definition: level1_cplx_impl.h:28

gtsam
Author(s):
autogenerated on Sun Dec 22 2024 04:11:54