Use the macro kcopy! instead of broadcast (#899)

Author: amontoison

src/bicgstab.jl

@@ -151,14 +151,14 @@ kwargs_bicgstab = (:c, :M, :N, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose,
       mul!(r₀, A, Δx)
       @kaxpby!(n, one(FC), b, -one(FC), r₀)
     else
-      r₀ .= b
+      @kcopy!(n, r₀, b)  # r₀ ← b
     end
 
     @kfill!(x, zero(FC))                # x₀
     @kfill!(s, zero(FC))                # s₀
     @kfill!(v, zero(FC))                # v₀
     MisI || mulorldiv!(r, M, r₀, ldiv)  # r₀
-    p .= r                              # p₁
+    @kcopy!(n, p, r)                    # p₁
 
     α = one(FC) # α₀
     ω = one(FC) # ω₀

src/car.jl

@@ -127,29 +127,29 @@ kwargs_car = (:M, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose, :history, :ca
       mul!(r, A, Δx)
       @kaxpby!(n, one(FC), b, -one(FC), r)
     else
-      r .= b
+      @kcopy!(n, r, b)  # r ← b
     end
 
     # p₀ = r₀ = M(b - Ax₀)
     if MisI
-      p .= r
+      @kcopy!(n, p, r)  # p ← r
     else
       mulorldiv!(p, M, r, ldiv)
-      r .= p
+      @kcopy!(n, r, p)  # r ← p
     end
 
     mul!(s, A, r)  # s₀ = Ar₀
 
     # q₀ = MAp₀ and s₀ = MAr₀
     if MisI
-      q .= s
+      @kcopy!(n, q, s)  # q ← s
     else
       mulorldiv!(q, M, s, ldiv)
-      s .= q
+      @kcopy!(n, s, q)  # s ← q
     end
 
     mul!(t, A, s)        # t₀ = As₀
-    u .= t               # u₀ = Aq₀
+    @kcopy!(n, u, t)     # u₀ = Aq₀
     ρ = @kdotr(n, t, s)  # ρ₀ = ⟨t₀ , s₀⟩
 
     # Compute ‖r₀‖

src/cg.jl

@@ -138,10 +138,10 @@ kwargs_cg = (:M, :ldiv, :radius, :linesearch, :atol, :rtol, :itmax, :timemax, :v
       mul!(r, A, Δx)
       @kaxpby!(n, one(FC), b, -one(FC), r)
     else
-      r .= b
+      @kcopy!(n, r, b)  # r ← b
     end
     MisI || mulorldiv!(z, M, r, ldiv)
-    p .= z
+    @kcopy!(n, p, z)  # p ← z
     γ = @kdotr(n, r, z)
     rNorm = sqrt(γ)
     history && push!(rNorms, rNorm)
@@ -182,7 +182,7 @@ kwargs_cg = (:M, :ldiv, :radius, :linesearch, :atol, :rtol, :itmax, :timemax, :v
           inconsistent = !linesearch
         end
         if linesearch
-          iter == 0 && (x .= b)
+          iter == 0 && @kcopy!(n, x, b)  # x ← b
           solved = true
         end
       end

src/cg_lanczos.jl

@@ -135,7 +135,7 @@ kwargs_cg_lanczos = (:M, :ldiv, :check_curvature, :atol, :rtol, :itmax, :timemax
       mul!(Mv, A, Δx)
       @kaxpby!(n, one(FC), b, -one(FC), Mv)
     else
-      Mv .= b
+      @kcopy!(n, Mv, b)  # Mv ← b
     end
     MisI || mulorldiv!(v, M, Mv, ldiv)  # v₁ = M⁻¹r₀
     β = sqrt(@kdotr(n, v, Mv))          # β₁ = v₁ᴴ M v₁
@@ -152,13 +152,13 @@ kwargs_cg_lanczos = (:M, :ldiv, :check_curvature, :atol, :rtol, :itmax, :timemax
       solver.warm_start = false
       return solver
     end
-    p .= v
+    @kcopy!(n, p, v)  # p ← v
 
     # Initialize Lanczos process.
     # β₁Mv₁ = b
     @kscal!(n, one(FC) / β, v)           # v₁  ←  v₁ / β₁
     MisI || @kscal!(n, one(FC) / β, Mv)  # Mv₁ ← Mv₁ / β₁
-    Mv_prev .= Mv
+    @kcopy!(n, Mv_prev, Mv)              # Mv_prev ← Mv
 
     iter = 0
     itmax == 0 && (itmax = 2 * n)

src/cg_lanczos_shift.jl

@@ -131,10 +131,10 @@ kwargs_cg_lanczos_shift = (:M, :ldiv, :check_curvature, :atol, :rtol, :itmax, :t
     for i = 1 : nshifts
       @kfill!(x[i], zero(FC))  # x₀
     end
-    Mv .= b                             # Mv₁ ← b
+    @kcopy!(n, Mv, b)                   # Mv₁ ← b
     MisI || mulorldiv!(v, M, Mv, ldiv)  # v₁ = M⁻¹ * Mv₁
     β = sqrt(@kdotr(n, v, Mv))          # β₁ = v₁ᴴ M v₁
-    rNorms .= β
+    @kfill!(rNorms, β)
     if history
       for i = 1 : nshifts
         push!(rNorms_history[i], rNorms[i])
@@ -155,14 +155,14 @@ kwargs_cg_lanczos_shift = (:M, :ldiv, :check_curvature, :atol, :rtol, :itmax, :t
 
     # Initialize each p to v.
     for i = 1 : nshifts
-      p[i] .= v
+      @kcopy!(n, p[i], v)  # pᵢ ← v
     end
 
     # Initialize Lanczos process.
     # β₁Mv₁ = b
     @kscal!(n, one(FC) / β, v)           # v₁  ←  v₁ / β₁
     MisI || @kscal!(n, one(FC) / β, Mv)  # Mv₁ ← Mv₁ / β₁
-    Mv_prev .= Mv
+    @kcopy!(n, Mv_prev, Mv)              # Mv_prev ← Mv
 
     # Initialize some constants used in recursions below.
     ρ = one(T)

src/cgls.jl

@@ -147,7 +147,7 @@ kwargs_cgls = (:M, :ldiv, :radius, :λ, :atol, :rtol, :itmax, :timemax, :verbose
     Mq = MisI ? q : solver.Mr
 
     @kfill!(x, zero(FC))
-    r .= b
+    @kcopy!(m, r, b)  # r ← b
     bNorm = @knrm2(m, r)   # Marginally faster than norm(b)
     if bNorm == 0
       stats.niter = 0
@@ -160,7 +160,7 @@ kwargs_cgls = (:M, :ldiv, :radius, :λ, :atol, :rtol, :itmax, :timemax, :verbose
     end
     MisI || mulorldiv!(Mr, M, r, ldiv)
     mul!(s, Aᴴ, Mr)
-    p .= s
+    @kcopy!(n, p, s)  # p ← s
     γ = @kdotr(n, s, s)  # γ = sᴴs
     iter = 0
     itmax == 0 && (itmax = m + n)

src/cgls_lanczos_shift.jl

@@ -137,11 +137,11 @@ kwargs_cgls_lanczos_shift = (:M, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose
       @kfill!(x[i], zero(FC))  # x₀
     end
 
-    u .= b
+    @kcopy!(m, u, b)             # u ← b
     @kfill!(u_prev, zero(FC))
-    mul!(v, Aᴴ, u)                      # v₁ ← Aᴴ * b
-    β = sqrt(@kdotr(n, v, v))           # β₁ = v₁ᵀ M v₁
-    rNorms .= β
+    mul!(v, Aᴴ, u)               # v₁ ← Aᴴ * b
+    β = sqrt(@kdotr(n, v, v))    # β₁ = v₁ᵀ M v₁
+    @kfill!(rNorms, β)
     if history
       for i = 1 : nshifts
         push!(rNorms_history[i], rNorms[i])
@@ -158,17 +158,17 @@ kwargs_cgls_lanczos_shift = (:M, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose
 
     # Initialize each p to v.
     for i = 1 : nshifts
-      p[i] .= v
+      @kcopy!(n, p[i], v)  # pᵢ ← v
     end
 
     # Initialize Lanczos process.
     # β₁v₁ = b
-    @kscal!(n, one(FC) / β, v)          # v₁  ←  v₁ / β₁
+    @kscal!(n, one(FC) / β, v)  # v₁ ← v₁ / β₁
     @kscal!(m, one(FC) / β, u)
 
     # Initialize some constants used in recursions below.
     ρ = one(T)
-    σ .= β
+    @kfill!(σ, β)
     @kfill!(δhat, zero(T))
     @kfill!(ω, zero(T))
     @kfill!(γ, one(T))
@@ -206,8 +206,8 @@ kwargs_cgls_lanczos_shift = (:M, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose
       β = sqrt(@kdotr(n, v, v))            # βₖ₊₁ = vₖ₊₁ᵀ M vₖ₊₁
       @kscal!(n, one(FC) / β, v)           # vₖ₊₁  ←  vₖ₊₁ / βₖ₊₁
       @kscal!(m, one(FC) / β, utilde)      # uₖ₊₁ = uₖ₊₁ / βₖ₊₁
-      u_prev .= u
-      u .= utilde
+      @kcopy!(m, u_prev, u)  # u_prev ← u
+      @kcopy!(m, u, utilde)  # u ← utilde
 
       MisI || (ρ = @kdotr(n, v, v))
       for i = 1 : nshifts

src/cgne.jl

@@ -152,7 +152,7 @@ kwargs_cgne = (:N, :ldiv, :λ, :atol, :rtol, :itmax, :timemax, :verbose, :histor
     z = NisI ? r : solver.z
 
     @kfill!(x, zero(FC))
-    r .= b
+    @kcopy!(m, r, b)  # r ← b
     NisI || mulorldiv!(z, N, r, ldiv)
     rNorm = @knrm2(m, r)   # Marginally faster than norm(r)
     history && push!(rNorms, rNorm)
@@ -163,7 +163,7 @@ kwargs_cgne = (:N, :ldiv, :λ, :atol, :rtol, :itmax, :timemax, :verbose, :histor
       stats.status = "x = 0 is a zero-residual solution"
       return solver
     end
-    λ > 0 && (s .= r)
+    λ > 0 && @kcopy!(m, s, r)  # s ← r
     mul!(p, Aᴴ, z)
 
     # Use ‖p‖ to detect inconsistent system.

src/cgs.jl

@@ -153,7 +153,7 @@ kwargs_cgs = (:c, :M, :N, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose, :hist
       mul!(r₀, A, Δx)
       @kaxpby!(n, one(FC), b, -one(FC), r₀)
     else
-      r₀ .= b
+      @kcopy!(n, r₀, b)  # r₀ ← b
     end
 
     @kfill!(x, zero(FC))                # x₀
@@ -189,9 +189,9 @@ kwargs_cgs = (:c, :M, :N, :ldiv, :atol, :rtol, :itmax, :timemax, :verbose, :hist
     (verbose > 0) && @printf(iostream, "%5s  %7s  %5s\n", "k", "‖rₖ‖", "timer")
     kdisplay(iter, verbose) && @printf(iostream, "%5d  %7.1e  %.2fs\n", iter, rNorm, ktimer(start_time))
 
-    u .= r               # u₀
-    p .= r               # p₀
-    @kfill!(q, zero(FC)) # q₋₁
+    @kcopy!(n, u, r)      # u₀
+    @kcopy!(n, p, r)      # p₀
+    @kfill!(q, zero(FC))  # q₋₁
 
     # Stopping criterion.
     solved = rNorm ≤ ε

src/cr.jl

@@ -146,9 +146,9 @@ kwargs_cr = (:M, :ldiv, :radius, :linesearch, :γ, :atol, :rtol, :itmax, :timema
       mul!(p, A, Δx)
       @kaxpby!(n, one(FC), b, -one(FC), p)
     else
-      p .= b
+      @kcopy!(n, p, b)  # p ← b
     end
-    MisI && (r .= p)
+    MisI && @kcopy!(n, r, p)  # r ← p
     MisI || mulorldiv!(r, M, p, ldiv)
     mul!(Ar, A, r)
     ρ = @kdotr(n, r, Ar)
@@ -165,8 +165,8 @@ kwargs_cr = (:M, :ldiv, :radius, :linesearch, :γ, :atol, :rtol, :itmax, :timema
       solver.warm_start = false
       return solver
     end
-    p .= r
-    q .= Ar
+    @kcopy!(n, p, r)   # p ← r
+    @kcopy!(n, q, Ar)  # q ← Ar
     (verbose > 0) && (m = zero(T)) # quadratic model
 
     iter = 0