Pullback API working

src/ReferenceFEs/Pullbacks.jl

@@ -31,7 +31,7 @@ end
 function Arrays.lazy_map(
   ::Broadcasting{typeof(gradient)}, a::LazyArray{<:Fill{Broadcasting{Operation{<:Pushforward}}}}
 )
-  cell_ref_fields, args = a.args
+  cell_ref_fields, args... = a.args
   cell_ref_gradient = lazy_map(Broadcasting(∇),cell_ref_fields)
   return lazy_map(a.maps.value,cell_ref_gradient,args...)
 end
@@ -123,15 +123,16 @@ function Arrays.lazy_map(
   ::typeof(evaluate),k::LazyArray{<:Fill{<:Pullback}},ref_cell_fields::AbstractArray
 )
   pb = k.maps.value
-  phys_cell_dofs, pf_args = k.args
+  phys_cell_dofs, pf_args... = k.args
   phys_cell_fields = lazy_map(pb.pushforward,ref_cell_fields,pf_args...)
   return lazy_map(evaluate,phys_cell_dofs,phys_cell_fields)
 end
 
 function evaluate!(
-  cache, k::PullBack, σ_phys::MomentBasedDofBasis, args...
+  cache, k::Pullback, σ_phys::AbstractVector{<:Dof}, args...
 )
-  Broadcasting(Operation(k))(f_phys,args...)
+  pf(f_ref) = evaluate(k.pushforward,f_ref,args...) # TODO: Can we avoid this? 
+  return Arrays.OperationMap(σ_phys,pf)
 end
 
 # InversePullback
@@ -145,7 +146,7 @@ where
   - V̂* is a dof space on the reference cell K̂ and 
   - V* is a dof space on the physical cell K.
 Its action on reference dofs ̂σ : V̂ -> R is defined in terms of the pushforward map F* as
-  σ = F**(̂σ) := ̂σ∘(F*)^-1 : V -> R
+  σ = (F**)^-1(̂σ) := ̂σ∘(F*)^-1 : V -> R
 """
 struct InversePullback{PF} <: Map
   pushforward::PF
@@ -162,11 +163,18 @@ function Arrays.lazy_map(
   ::typeof(evaluate),k::LazyArray{<:Fill{<:InversePullback}},phys_cell_fields::AbstractArray
 )
   pb = inverse_map(k.maps.value)
-  ref_cell_dofs, pf_args = k.args
+  ref_cell_dofs, pf_args... = k.args
   ref_cell_fields = lazy_map(inverse_map(pb.pushforward),phys_cell_fields,pf_args...)
   return lazy_map(evaluate,ref_cell_dofs,ref_cell_fields)
 end
 
+function evaluate!(
+  cache, k::InversePullback, σ_ref::AbstractVector{<:Dof}, args...
+)
+  ipf(f_phys) = evaluate(inverse_map(k.pushforward),f_phys,args...) # TODO: Can we avoid this? 
+  return Arrays.OperationMap(σ_ref,ipf)
+end
+
 # ContraVariantPiolaMap
 
 struct ContraVariantPiolaMap <: Pushforward end

test/pullbacks.jl

@@ -1,7 +1,8 @@
 
 using FillArrays
 using Gridap
-using Gridap.ReferenceFEs, Gridap.FESpaces
+using Gridap.ReferenceFEs, Gridap.FESpaces, Gridap.CellData
+using Gridap.Fields
 
 using Gridap.ReferenceFEs: Pullback
 
@@ -29,15 +30,25 @@ Arr = lazy_map(evaluate,σ_ref,φ_ref)[1]
 
 pf  = ContraVariantPiolaMap()
 
+# Inverse Pushforward
 ipf = inverse_map(pf)
 f_ref = lazy_map(ipf,φ_phys,Jt,sign)
 Brr = lazy_map(evaluate,σ_ref,f_ref)[1]
 Brr == Arr
 φ_ref_x = lazy_map(evaluate,φ_ref,pts)[1]
 f_ref_x = lazy_map(evaluate,f_ref,pts)[1]
-f_ref_x == φ_ref_x
+f_ref_x ≈ φ_ref_x
 
+# Pullback
 pb  = Pullback(pf)
 θ_ref = lazy_map(pb,σ_phys,Jt,sign)
+θ_ref_x = lazy_map(evaluate,θ_ref,φ_ref)
+σ_ref_x = lazy_map(evaluate,σ_ref,φ_ref)
+θ_ref_x ≈ σ_ref_x
 
+# Inverse Pullback
 ipb = inverse_map(pb)
+θ_phys = lazy_map(ipb,σ_ref,Jt,sign)
+θ_phys_x = lazy_map(evaluate,θ_phys,φ_phys)
+σ_phys_x = lazy_map(evaluate,σ_phys,φ_phys)
+θ_phys_x ≈ σ_phys_x