vendor/gonum.org/v1/gonum/lapack/gonum/dlascl.go

   1 // Copyright ©2015 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.
   4
   5 package gonum
   6
   7 import (
   8         "math"
   9
  10         "gonum.org/v1/gonum/lapack"
  11 )
  12
  13 // Dlascl multiplies an m×n matrix by the scalar cto/cfrom.
  14 //
  15 // cfrom must not be zero, and cto and cfrom must not be NaN, otherwise Dlascl
  16 // will panic.
  17 //
  18 // Dlascl is an internal routine. It is exported for testing purposes.
  19 func (impl Implementation) Dlascl(kind lapack.MatrixType, kl, ku int, cfrom, cto float64, m, n int, a []float64, lda int) {
  20         checkMatrix(m, n, a, lda)
  21         if cfrom == 0 {
  22                 panic(zeroDiv)
  23         }
  24         if math.IsNaN(cfrom) || math.IsNaN(cto) {
  25                 panic(nanScale)
  26         }
  27         if n == 0 || m == 0 {
  28                 return
  29         }
  30         smlnum := dlamchS
  31         bignum := 1 / smlnum
  32         cfromc := cfrom
  33         ctoc := cto
  34         cfrom1 := cfromc * smlnum
  35         for {
  36                 var done bool
  37                 var mul, ctol float64
  38                 if cfrom1 == cfromc {
  39                         // cfromc is inf.
  40                         mul = ctoc / cfromc
  41                         done = true
  42                         ctol = ctoc
  43                 } else {
  44                         ctol = ctoc / bignum
  45                         if ctol == ctoc {
  46                                 // ctoc is either 0 or inf.
  47                                 mul = ctoc
  48                                 done = true
  49                                 cfromc = 1
  50                         } else if math.Abs(cfrom1) > math.Abs(ctoc) && ctoc != 0 {
  51                                 mul = smlnum
  52                                 done = false
  53                                 cfromc = cfrom1
  54                         } else if math.Abs(ctol) > math.Abs(cfromc) {
  55                                 mul = bignum
  56                                 done = false
  57                                 ctoc = ctol
  58                         } else {
  59                                 mul = ctoc / cfromc
  60                                 done = true
  61                         }
  62                 }
  63                 switch kind {
  64                 default:
  65                         panic("lapack: not implemented")
  66                 case lapack.General:
  67                         for i := 0; i < m; i++ {
  68                                 for j := 0; j < n; j++ {
  69                                         a[i*lda+j] = a[i*lda+j] * mul
  70                                 }
  71                         }
  72                 case lapack.UpperTri:
  73                         for i := 0; i < m; i++ {
  74                                 for j := i; j < n; j++ {
  75                                         a[i*lda+j] = a[i*lda+j] * mul
  76                                 }
  77                         }
  78                 case lapack.LowerTri:
  79                         for i := 0; i < m; i++ {
  80                                 for j := 0; j <= min(i, n-1); j++ {
  81                                         a[i*lda+j] = a[i*lda+j] * mul
  82                                 }
  83                         }
  84                 }
  85                 if done {
  86                         break
  87                 }
  88         }
  89 }