1 // Copyright ©2017 The Gonum Authors. All rights reserved.
2 // Use of this source code is governed by a BSD-style
3 // license that can be found in the LICENSE file.
10 "gonum.org/v1/gonum/internal/asm/c128"
13 // Dzasum returns the sum of the absolute values of the elements of x
14 // \sum_i |Re(x[i])| + |Im(x[i])|
15 // Dzasum returns 0 if incX is negative.
16 func (Implementation) Dzasum(n int, x []complex128, incX int) float64 {
31 for _, v := range x[:n] {
36 if (n-1)*incX >= len(x) {
39 for i := 0; i < n; i++ {
46 // Dznrm2 computes the Euclidean norm of the complex vector x,
47 // ‖x‖_2 = sqrt(\sum_i x[i] * conj(x[i])).
48 // This function returns 0 if incX is negative.
49 func (Implementation) Dznrm2(n int, x []complex128, incX int) float64 {
62 if (n-1)*incX >= len(x) {
70 for _, v := range x[:n] {
71 re, im := math.Abs(real(v)), math.Abs(imag(v))
74 ssq = 1 + ssq*(scale/re)*(scale/re)
77 ssq += (re / scale) * (re / scale)
82 ssq = 1 + ssq*(scale/im)*(scale/im)
85 ssq += (im / scale) * (im / scale)
89 if math.IsInf(scale, 1) {
92 return scale * math.Sqrt(ssq)
94 for ix := 0; ix < n*incX; ix += incX {
95 re, im := math.Abs(real(x[ix])), math.Abs(imag(x[ix]))
98 ssq = 1 + ssq*(scale/re)*(scale/re)
101 ssq += (re / scale) * (re / scale)
106 ssq = 1 + ssq*(scale/im)*(scale/im)
109 ssq += (im / scale) * (im / scale)
113 if math.IsInf(scale, 1) {
116 return scale * math.Sqrt(ssq)
119 // Izamax returns the index of the first element of x having largest |Re(·)|+|Im(·)|.
120 // Izamax returns -1 if n is 0 or incX is negative.
121 func (Implementation) Izamax(n int, x []complex128, incX int) int {
126 // Return invalid index.
131 // Return invalid index.
136 if len(x) <= (n-1)*incX {
142 for i, v := range x[1:n] {
152 for i := 1; i < n; i++ {
153 absV := dcabs1(x[ix])
163 // Zaxpy adds alpha times x to y:
164 // y[i] += alpha * x[i] for all i
165 func (Implementation) Zaxpy(n int, alpha complex128, x []complex128, incX int, y []complex128, incY int) {
178 if (incX > 0 && (n-1)*incX >= len(x)) || (incX < 0 && (1-n)*incX >= len(x)) {
181 if (incY > 0 && (n-1)*incY >= len(y)) || (incY < 0 && (1-n)*incY >= len(y)) {
187 if incX == 1 && incY == 1 {
188 c128.AxpyUnitary(alpha, x[:n], y[:n])
198 c128.AxpyInc(alpha, x, y, uintptr(n), uintptr(incX), uintptr(incY), uintptr(ix), uintptr(iy))
201 // Zcopy copies the vector x to vector y.
202 func (Implementation) Zcopy(n int, x []complex128, incX int, y []complex128, incY int) {
215 if (incX > 0 && (n-1)*incX >= len(x)) || (incX < 0 && (1-n)*incX >= len(x)) {
218 if (incY > 0 && (n-1)*incY >= len(y)) || (incY < 0 && (1-n)*incY >= len(y)) {
221 if incX == 1 && incY == 1 {
232 for i := 0; i < n; i++ {
239 // Zdotc computes the dot product
241 // of two complex vectors x and y.
242 func (Implementation) Zdotc(n int, x []complex128, incX int, y []complex128, incY int) complex128 {
255 if incX == 1 && incY == 1 {
262 return c128.DotcUnitary(x[:n], y[:n])
271 if ix >= len(x) || (n-1)*incX >= len(x) {
274 if iy >= len(y) || (n-1)*incY >= len(y) {
277 return c128.DotcInc(x, y, uintptr(n), uintptr(incX), uintptr(incY), uintptr(ix), uintptr(iy))
280 // Zdotu computes the dot product
282 // of two complex vectors x and y.
283 func (Implementation) Zdotu(n int, x []complex128, incX int, y []complex128, incY int) complex128 {
296 if incX == 1 && incY == 1 {
303 return c128.DotuUnitary(x[:n], y[:n])
312 if ix >= len(x) || (n-1)*incX >= len(x) {
315 if iy >= len(y) || (n-1)*incY >= len(y) {
318 return c128.DotuInc(x, y, uintptr(n), uintptr(incX), uintptr(incY), uintptr(ix), uintptr(iy))
321 // Zdscal scales the vector x by a real scalar alpha.
322 // Zdscal has no effect if incX < 0.
323 func (Implementation) Zdscal(n int, alpha float64, x []complex128, incX int) {
330 if (n-1)*incX >= len(x) {
347 for ix := 0; ix < n*incX; ix += incX {
354 for i, v := range x {
355 x[i] = complex(alpha*real(v), alpha*imag(v))
359 for ix := 0; ix < n*incX; ix += incX {
361 x[ix] = complex(alpha*real(v), alpha*imag(v))
365 // Zscal scales the vector x by a complex scalar alpha.
366 // Zscal has no effect if incX < 0.
367 func (Implementation) Zscal(n int, alpha complex128, x []complex128, incX int) {
374 if (n-1)*incX >= len(x) {
391 for ix := 0; ix < n*incX; ix += incX {
397 c128.ScalUnitary(alpha, x[:n])
400 c128.ScalInc(alpha, x, uintptr(n), uintptr(incX))
403 // Zswap exchanges the elements of two complex vectors x and y.
404 func (Implementation) Zswap(n int, x []complex128, incX int, y []complex128, incY int) {
417 if (incX > 0 && (n-1)*incX >= len(x)) || (incX < 0 && (1-n)*incX >= len(x)) {
420 if (incY > 0 && (n-1)*incY >= len(y)) || (incY < 0 && (1-n)*incY >= len(y)) {
423 if incX == 1 && incY == 1 {
425 for i, v := range x {
437 for i := 0; i < n; i++ {
438 x[ix], y[iy] = y[iy], x[ix]