add pairwise! for Covariance types

src/covariance.jl

@@ -18,15 +18,31 @@ Evaluate the covariance at objects `x₁` and `x₁`.
 
 """
     pairwise(cov, domain)
-    
+
 Evaluate covariance `cov` between all elements in the `domain`.
-    
+
     pairwise(cov, domain₁, domain₂)
 
 Evaluate covariance `cov` between all elements of `domain₁` and `domain₂`.
 """
 pairwise(cov::Covariance, args...) = sill(cov.γ) .- pairwise(cov.γ, args...)
 
+"""
+    pairwise!(Γ, cov, domain)
+
+Evaluates covariance `cov` between all elements in the `domain` in-place, filling the matrix `Γ`.
+
+    pairwise!(Γ, cov, domain₁, domain₂)
+
+Evaluates covariance `cov` between all elements of `domain₁` and `domain₂` in-place, filling the matrix `Γ`."""
+function pairwise!(Γ, cov::Covariance, args...)
+  pairwise!(Γ, cov.γ, args...)
+  for i ∈ eachindex(Γ)
+    Γ[i] = sill(cov.γ) - Γ[i]
+  end
+  Γ
+end
+
 # -----------
 # IO METHODS
 # -----------

src/variogram.jl

@@ -152,6 +152,11 @@ function pairwise(γ::Variogram, domain)
   pairwise!(Γ, γ, domain)
 end
 
+"""
+    pairwise!(Γ, γ, domain)
+
+Evaluates covariance `γ` between all elements in the `domain` in-place, filling the matrix `Γ`.
+"""
 function pairwise!(Γ, γ::Variogram, domain)
   n = length(domain)
   @inbounds for j in 1:n
@@ -185,6 +190,11 @@ function pairwise(γ::Variogram, domain₁, domain₂)
   pairwise!(Γ, γ, domain₁, domain₂)
 end
 
+"""
+    pairwise!(Γ, γ, domain)
+
+Evaluates covariance `γ` between all elements of `domain₁` and `domain₂` in-place, filling the matrix `Γ`.
+"""
 function pairwise!(Γ, γ::Variogram, domain₁, domain₂)
   m = length(domain₁)
   n = length(domain₂)

test/covariance.jl

@@ -40,6 +40,11 @@
   @test eltype(Γ_f) == Float32
   @test issymmetric(Γ_f)
 
+  𝒟 = PointSet(Matrix(1.0I, 3, 3))
+  Γ = Matrix{Float64}(undef, 3, 3)
+  GeoStatsFunctions.pairwise!(Γ, GaussianCovariance(range=1.0, sill=1.0, nugget=1.0), 𝒟)
+  @test issymmetric(Γ)
+
   # shows
   cov = CircularCovariance()
   @test sprint(show, cov) == "CircularCovariance(sill: 1.0, nugget: 0.0, range: 1.0, distance: Euclidean)"