some comments and coefficient type change

Author: ajozefiak

src/derivative_operators/derivative_operator.jl

@@ -117,7 +117,7 @@ function UpwindDifference{N}(derivative_order::Int,
     low_boundary_coefs      = convert(SVector{boundary_point_count},_low_boundary_coefs)
     high_boundary_coefs     = convert(SVector{boundary_point_count},reverse(SVector{boundary_stencil_length, T}[reverse(low_boundary_coefs[i]) for i in 1:boundary_point_count]))
 
-    coefficients            = Vector{T}(undef,len)
+    coefficients            = nothing
 
     DerivativeOperator{T,N,true,T,typeof(stencil_coefs),
         typeof(low_boundary_coefs),Vector{T},

src/derivative_operators/derivative_operator_functions.jl

@@ -1,3 +1,5 @@
+
+# Fallback mul! implementation for a single DerivativeOperator operating on an AbstractArray
 function LinearAlgebra.mul!(x_temp::AbstractArray{T}, A::DerivativeOperator{T,N}, M::AbstractArray{T}) where {T<:Real,N}
 
     # Check that x_temp has correct dimensions
@@ -25,6 +27,8 @@ function LinearAlgebra.mul!(x_temp::AbstractArray{T}, A::DerivativeOperator{T,N}
     end
 end
 
+# A more efficient mul! implementation for a single, regular-grid, centered difference
+# DerivativeOperator operating on a 2D or 3D AbstractArray
 for MT in [2,3]
     @eval begin
         function LinearAlgebra.mul!(x_temp::AbstractArray{T,$MT}, A::DerivativeOperator{T,N,false,T2,S1}, M::AbstractArray{T,$MT}) where {T<:Real,N,T2,SL,S1<:SArray{Tuple{SL},T,1,SL}}
@@ -92,9 +96,17 @@ function *(A::DerivativeOperator{T,N},M::AbstractArray{T}) where {T<:Real,N}
     return x_temp
 end
 
-# Will eventually implement proper dispatch
+# A more efficient mul! implementation for a composition of regular-grid, centered difference
+# DerivativeOperator operating on a 2D or 3D AbstractArray
 function LinearAlgebra.mul!(x_temp::AbstractArray{T,2}, A::AbstractDiffEqCompositeOperator, M::AbstractArray{T,2}) where {T}
 
+    #Check that composite operator satisfies: regular-grid, centered difference:
+    #for L in A.ops
+    #    if L ... (does not satisfy conditions)
+    #        return (call fall back for multiplication of composite operators)
+    #    end
+    #end
+
     ndimsM = ndims(M)
     Wdims = ones(Int64,ndimsM)
     pad = zeros(Int64, ndimsM)
@@ -121,7 +133,7 @@ function LinearAlgebra.mul!(x_temp::AbstractArray{T,2}, A::AbstractDiffEqComposi
         for i in offset+1:Wdims[axis]-offset
             idx[axis]=i
 
-            W[idx...] += s[i-offset] 
+            W[idx...] += s[i-offset]
 
             idx[axis] = mid_Wdims[axis]
         end