From 9f36ea6762a0495796eed161c48ec0a4fab4e862 Mon Sep 17 00:00:00 2001
From: Jishnu Bhattacharya <jishnub.github@gmail.com>
Date: Sat, 22 Aug 2026 21:56:32 +0400
Subject: [PATCH] Define values for ArraySpace

`values(f)` errored for an array-valued `Fun`, which broke `iszero`
since it started consulting `values` in #635. Three separate causes:

* `itransform` was only defined for `VectorSpace`, so a `MatrixSpace`
  fell through to the generic `Space` method and hit
  `checkcanonicalspace`, which errors as an `ArraySpace` of canonical
  spaces is its own canonical space. Widen the method to `ArraySpace`
  and reshape the values to the shape of the space.

* The component coefficients may be views, so `pad!` may not resize
  them. Use `pad` instead.

* `_values` asserted a `Vector{float(T)}` return for any space over a
  numeric domain, which includes array-valued ones. Restrict the
  assertion to spaces with a scalar range type.

Also fix `points(f::Fun{<:ArraySpace})`, which passed the total number
of interlaced coefficients where a per-component count is expected, so
that `values(f) == f.(points(f))`.

Fixes #647

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 src/Fun.jl               |  3 ++-
 src/Spaces/ArraySpace.jl | 11 +++++++++--
 test/SpacesTest.jl       | 25 +++++++++++++++++++++++++
 3 files changed, 36 insertions(+), 3 deletions(-)

diff --git a/src/Fun.jl b/src/Fun.jl
index a6ac99c1..3bd32c64 100644
--- a/src/Fun.jl
+++ b/src/Fun.jl
@@ -426,7 +426,8 @@ true
 """
 values(f::Fun,dat...) = _values(f.space, f.coefficients, dat...)
 _values(sp, v, dat...) = itransform(sp, v, dat...)
-_values(sp::UnivariateSpace, v::Vector{T}, dat...) where {T<:Number} =
+# the return type may be asserted only for scalar-valued spaces
+_values(sp::Space{<:Domain{<:Number},<:Number}, v::Vector{T}, dat...) where {T<:Number} =
     itransform(sp, v, dat...)::Vector{float(T)}
 
 """
diff --git a/src/Spaces/ArraySpace.jl b/src/Spaces/ArraySpace.jl
index e3c60041..2c52b237 100644
--- a/src/Spaces/ArraySpace.jl
+++ b/src/Spaces/ArraySpace.jl
@@ -74,6 +74,11 @@ setdomain(A::ArraySpace,d::Domain) = ArraySpace(map(sp->setdomain(sp,d),A.spaces
 
 #TODO: rework for different spaces
 points(d::ArraySpace,n) = points(d.spaces[1],n)
+# points(::ArraySpace, n) forwards n to a component space, so n must be the number of
+# coefficients per component. The generic points(f::Fun) would pass ncoefficients(f),
+# which counts the interlaced coefficients of all the components together.
+# The count below matches the padding in itransform, so that values(f) == f.(points(f)).
+points(f::Fun{<:ArraySpace}) = points(space(f), maximum(ncoefficients, vec(f), init=0))
 
 
 transform(AS::ArraySpace{SS,1},vals::AbstractVector{Vector{V}}) where {SS,V} =
@@ -98,10 +103,12 @@ transform(AS::VectorSpace{SS},vals::AbstractVector{AV}) where {SS,AV<:AbstractVe
 transform(AS::VectorSpace{SS},vals::AbstractVector{SVector{V,n}}) where {SS,n,V} =
     transform(AS,map(Vector,vals))
 
-function itransform(AS::VectorSpace,cfs::AbstractVector)
+function itransform(AS::ArraySpace,cfs::AbstractVector)
     vf = vec(Fun(AS, cfs))
     n = maximum(ncoefficients, vf)
-    vcat.(values.(pad!.(vf, n))...)
+    vals = vcat.(values.(pad.(vf, n))...)
+    # the values of an array-valued Fun are arrays of the same shape
+    [reshape(v, size(AS)) for v in vals]
 end
 
 
diff --git a/test/SpacesTest.jl b/test/SpacesTest.jl
index 88fb2acf..46d9bdea 100644
--- a/test/SpacesTest.jl
+++ b/test/SpacesTest.jl
@@ -395,5 +395,30 @@ using Test
         A = ApproxFunBase.ArraySpace(empty!([PointSpace(1:3)]))
         @test length(A) == 0
         @test ApproxFunBase.dimension(A) == 0
+
+        @testset "values" begin
+            # the two spaces exercise different paths: an ArraySpace of
+            # HeavisideSpaces lies over a numeric domain and has a trivial
+            # interlacer, so that the component coefficients are views
+            @testset for S in (PointSpace(1:3), ApproxFunBase.HeavisideSpace([-1.0, -0.5, 0.0, 1.0]))
+                @testset for sz in ((2,), (2,2))
+                    A = ApproxFunBase.ArraySpace(S, sz...)
+                    # 3 coefficients per component, so that each component
+                    # may be evaluated on the full grid
+                    n = 3length(A)
+                    f = Fun(A, Float64.(1:n))
+                    v = values(f)
+                    @test all(x -> size(x) == sz, v)
+                    @test v == f.(points(f))
+                    @test !iszero(f)
+                    @test iszero(Fun(A, zeros(n)))
+                end
+            end
+
+            # static array space
+            A = ApproxFunBase.ArraySpace(PointSpace(1:3), Val((2,2)))
+            f = Fun(A, Float64.(1:3length(A)))
+            @test values(f) == f.(points(f))
+        end
     end
 end