From 981a832a1d40de0997b6e313e9305925ec6053af Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Wed, 12 Oct 2022 14:43:05 +1100
Subject: [PATCH 01/17] added some missing Indexing tests with `end` usage

---
 test/TensorValuesTests/IndexingTests.jl | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/test/TensorValuesTests/IndexingTests.jl b/test/TensorValuesTests/IndexingTests.jl
index 4f7c20c1e..5f6951076 100644
--- a/test/TensorValuesTests/IndexingTests.jl
+++ b/test/TensorValuesTests/IndexingTests.jl
@@ -14,6 +14,7 @@ v = VectorValue{4}(a)
 @test size(v) == (4,)
 @test length(v) == 4
 @test lastindex(v) == length(v)
+@test v[end] == a[end]
 
 for (k,i) in enumerate(eachindex(v))
   @test v[i] == a[k]
@@ -24,6 +25,7 @@ t = TensorValue{2}(a)
 @test size(t) == (2,2)
 @test length(t) == 4
 @test lastindex(t) == length(t)
+@test t[end] == a[end]
 
 for (k,i) in enumerate(eachindex(t))
   @test t[i] == a[k]
@@ -43,6 +45,7 @@ t = TensorValue(convert(SMatrix{2,2,Int},s))
 @test size(s) == (2,2)
 @test length(s) == 4
 @test lastindex(s) == length(s)
+@test s[end] == 22 
 
 for (k,i) in enumerate(eachindex(t))
     @test s[i] == t[k]

From ac342e1a19613279a3d44af8e96617f6422f95c0 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Wed, 12 Oct 2022 16:23:10 +1100
Subject: [PATCH 02/17] first attempt with CartesianDiscreteModel, works only
 in 1D

---
 src/CellData/DiracDeltas.jl | 30 ++++++++++++++++++++++++++++++
 1 file changed, 30 insertions(+)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 4cfea2f3f..f4d6a9da9 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -76,3 +76,33 @@ end
 function get_triangulation(d::DiracDelta)
   d.Γ
 end
+
+# For handling DiracDelta at a generic Point in the domain #
+
+#=
+The idea would be to use the existing structure of DiracDelta.
+In place of triangulation, we need something similar for Measure
+it should be in PhysicalSpace to not repeat integration
+reuse the cache from KDTree , if possible.
+
+we also need to add the contribution to the respective cell, which apparently
+was not done in the existing implementation of DiracDelta. One way to this is
+to add the contribution to the respective cell where the point belongs to!
+So that there is no problem with AD happening cell wise
+=#
+function DiracDelta(x::Point)
+  DiracDelta{0}(x)
+end
+
+function DiracDelta{0}(x::Point{D,T}) where {D,T}
+  desc = CartesianDescriptor(x, ntuple(i->zero(T),Val(D)), ntuple(i->1,Val(D)))
+  point_model = CartesianDiscreteModel(desc)
+  point_Ω = Triangulation(point_model)
+  point_Γ = BoundaryTriangulation(point_model)
+  quad = GenericQuadrature([x],ones(T,1),"DiracDelta point node")
+  cell_quad = CellQuadrature(SA[quad],SA[quad.coordinates],SA[quad.weights],point_Ω, PhysicalDomain(),PhysicalDomain())
+  dx = Measure(cell_quad)
+  DiracDelta{0}(point_Γ,dx)
+end
+
+# goal is to build a DiscreteModel of just a single Point 

From 4cda93dd1d7c8e3b468f7f690e33fe5931977536 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Wed, 12 Oct 2022 18:41:31 +1100
Subject: [PATCH 03/17] general point DiracDelta working, need to fix
 application with CellField

---
 src/CellData/DiracDeltas.jl | 37 ++++++++++++++++---------------------
 1 file changed, 16 insertions(+), 21 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index f4d6a9da9..f92846558 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -79,30 +79,25 @@ end
 
 # For handling DiracDelta at a generic Point in the domain #
 
-#=
-The idea would be to use the existing structure of DiracDelta.
-In place of triangulation, we need something similar for Measure
-it should be in PhysicalSpace to not repeat integration
-reuse the cache from KDTree , if possible.
-
-we also need to add the contribution to the respective cell, which apparently
+#= TO DO
+We also need to add the contribution to the respective cell, which apparently
 was not done in the existing implementation of DiracDelta. One way to this is
 to add the contribution to the respective cell where the point belongs to!
-So that there is no problem with AD happening cell wise
+So that there is no problem with AD happening cell wise.
 =#
-function DiracDelta(x::Point)
-  DiracDelta{0}(x)
+function DiracDelta(x::Point{D}, model::DiscreteModel{D}) where D
+  DiracDelta{0}(x, model)
 end
 
-function DiracDelta{0}(x::Point{D,T}) where {D,T}
-  desc = CartesianDescriptor(x, ntuple(i->zero(T),Val(D)), ntuple(i->1,Val(D)))
-  point_model = CartesianDiscreteModel(desc)
-  point_Ω = Triangulation(point_model)
-  point_Γ = BoundaryTriangulation(point_model)
-  quad = GenericQuadrature([x],ones(T,1),"DiracDelta point node")
-  cell_quad = CellQuadrature(SA[quad],SA[quad.coordinates],SA[quad.weights],point_Ω, PhysicalDomain(),PhysicalDomain())
-  dx = Measure(cell_quad)
-  DiracDelta{0}(point_Γ,dx)
+function DiracDelta{0}(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
+  # check if the point is inside an active cell and save the caches and cell
+  # as wouldn't be caught for user-defined functions (i.e. not CellFields)
+  trian = Triangulation(model)
+  cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
+  cell = _point_to_cell!(cache1, x) # throws error if Point not in domain
+  point_grid = UnstructuredGrid([x])
+  point_model = UnstructuredDiscreteModel(point_grid)
+  point_trian = Triangulation(point_model)
+  dx = Measure(point_trian,1)
+  DiracDelta{0}(point_trian,dx)
 end
-
-# goal is to build a DiscreteModel of just a single Point 

From 5bca9855a974edadc2908fb51c248c6fa53c62e8 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 02:07:06 +1100
Subject: [PATCH 04/17] Added support for DiracDelta at generic Points

---
 src/CellData/DiracDeltas.jl            | 47 +++++++++++++++++++-------
 test/CellDataTests/DiracDeltasTests.jl | 16 +++++++++
 2 files changed, 50 insertions(+), 13 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index f92846558..258cd20e8 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -1,7 +1,15 @@
+abstract type DiracDeltaSupportType <: GridapType end
+struct IsGridEntity <: DiracDeltaSupportType end
+struct NotGridEntity <: DiracDeltaSupportType end
 
-struct DiracDelta{D} <: GridapType
+const _is_grid_entity = IsGridEntity()
+const _not_grid_entity = NotGridEntity()
+
+
+struct DiracDelta{D,S<:DiracDeltaSupportType} <: GridapType
   Γ::Triangulation{D}
   dΓ::Measure
+  supp::S
 end
 
 function DiracDelta{D}(
@@ -24,7 +32,7 @@ function DiracDelta{D}(
   trian = BodyFittedTriangulation(model,face_grid,face_to_bgface)
   Γ = BoundaryTriangulation(trian,glue)
   dΓ = Measure(Γ,degree)
-  DiracDelta(Γ,dΓ)
+  DiracDelta(Γ,dΓ,_is_grid_entity)
 end
 
 function DiracDelta{D}(
@@ -79,19 +87,13 @@ end
 
 # For handling DiracDelta at a generic Point in the domain #
 
-#= TO DO
-We also need to add the contribution to the respective cell, which apparently
-was not done in the existing implementation of DiracDelta. One way to this is
-to add the contribution to the respective cell where the point belongs to!
-So that there is no problem with AD happening cell wise.
-=#
 function DiracDelta(x::Point{D}, model::DiscreteModel{D}) where D
-  DiracDelta{0}(x, model)
+  DiracDelta(x, model)
 end
 
-function DiracDelta{0}(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
-  # check if the point is inside an active cell and save the caches and cell
-  # as wouldn't be caught for user-defined functions (i.e. not CellFields)
+function DiracDelta(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
+  # check if the point is inside an active cell, as it wouldn't be caught for
+  # user-defined functions (i.e. which are not CellFields)
   trian = Triangulation(model)
   cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
   cell = _point_to_cell!(cache1, x) # throws error if Point not in domain
@@ -99,5 +101,24 @@ function DiracDelta{0}(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
   point_model = UnstructuredDiscreteModel(point_grid)
   point_trian = Triangulation(point_model)
   dx = Measure(point_trian,1)
-  DiracDelta{0}(point_trian,dx)
+  DiracDelta(point_trian,dx,_not_grid_entity)
+end
+
+function DiracDelta(v::Vector{Point{D,T}},model::DiscreteModel{D}) where {D,T}
+  # check if the points are inside an active cell
+  trian = Triangulation(model)
+  cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
+  cell = map(x -> _point_to_cell!(cache1, x), v)
+  # throws error if any Point not in domain
+  point_grid = UnstructuredGrid(v)
+  point_model = UnstructuredDiscreteModel(point_grid)
+  point_trian = Triangulation(point_model)
+  dx = Measure(point_trian,1)
+  DiracDelta(point_trian,dx,_not_grid_entity)
+end
+
+function evaluate!(cache,d::DiracDelta{0,NotGridEntity},f::CellField)
+  eval = lazy_map(f,d.Γ.grid.node_coordinates)
+  dc = DomainContribution()
+  add_contribution!(dc, d.Γ, eval)
 end
diff --git a/test/CellDataTests/DiracDeltasTests.jl b/test/CellDataTests/DiracDeltasTests.jl
index 051779abe..f2732cc25 100644
--- a/test/CellDataTests/DiracDeltasTests.jl
+++ b/test/CellDataTests/DiracDeltasTests.jl
@@ -37,6 +37,22 @@ degree = 2
 δ = DiracDelta{0}(model,tags=[2,4])
 @test sum(δ(v)) ≈ 5
 
+# Tests for DiracDelta at a generic Point in the domain #
+
+p = VectorValue(0.2,0.3)
+δ = DiracDelta(p,model)
+@test sum(δ(v)) == v(p)
+
+pvec = [p,π*p]
+δ = DiracDelta(pvec,model)
+@test sum(δ(v)) == sum(v(pvec))
+
+# below passes but cannot import FESpace and FEFunction here
+# reffe = ReferenceFE(lagrangian,Float64,3)
+# V = FESpace(model,reffe,conformity=:L2)
+# uh = FEFunction(V,rand(num_free_dofs(V)))
+# @test sum(d(uh)) == uh(p)
+
 #using Gridap
 #
 #order = 2

From b13aa3261a51ea3f00add45c89d2e7384b765e76 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 02:47:06 +1100
Subject: [PATCH 05/17] fixed recursive DiracDelta calls

---
 src/CellData/DiracDeltas.jl | 4 ----
 1 file changed, 4 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 258cd20e8..e0b76e576 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -87,10 +87,6 @@ end
 
 # For handling DiracDelta at a generic Point in the domain #
 
-function DiracDelta(x::Point{D}, model::DiscreteModel{D}) where D
-  DiracDelta(x, model)
-end
-
 function DiracDelta(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
   # check if the point is inside an active cell, as it wouldn't be caught for
   # user-defined functions (i.e. which are not CellFields)

From db790bfcdc31e5690223f95caa513c2a69e4f833 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 10:21:06 +1100
Subject: [PATCH 06/17] Added DiracDelta for generic points PR to NEWS.md [skip
 ci]

---
 NEWS.md | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/NEWS.md b/NEWS.md
index e97c1dab8..d443ddbdd 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -10,7 +10,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - `lastindex` for `MultiValue`s for consistent usage of `[end]` as per `length`. Since PR [#834](https://github.com/gridap/Gridap.jl/pull/834)
 - BDM (Brezzi-Douglas-Marini) ReferenceFEs in PR [#823](https://github.com/gridap/Gridap.jl/pull/823)
 - Implemented `ConstantFESpace`. This space allows, e.g., to impose the mean value of the pressure via an additional unknown/equation (a Lagrange multiplier). Since PR [#836](https://github.com/gridap/Gridap.jl/pull/836)
-
+- Added support to construct `DiracDelta` at generic `Point`(s) in the domain. Since PR [#837](https://github.com/gridap/Gridap.jl/pull/837)
 ## [0.17.14] - 2022-07-29
 
 ### Added

From 546ca0f09f48f57379015ccddae54811b47e6e4d Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 10:31:47 +1100
Subject: [PATCH 07/17] removing some supporting comments

---
 src/CellData/DiracDeltas.jl            | 6 +-----
 test/CellDataTests/DiracDeltasTests.jl | 2 +-
 2 files changed, 2 insertions(+), 6 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index e0b76e576..04e43c1d5 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -88,11 +88,9 @@ end
 # For handling DiracDelta at a generic Point in the domain #
 
 function DiracDelta(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
-  # check if the point is inside an active cell, as it wouldn't be caught for
-  # user-defined functions (i.e. which are not CellFields)
   trian = Triangulation(model)
   cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
-  cell = _point_to_cell!(cache1, x) # throws error if Point not in domain
+  cell = _point_to_cell!(cache1, x)
   point_grid = UnstructuredGrid([x])
   point_model = UnstructuredDiscreteModel(point_grid)
   point_trian = Triangulation(point_model)
@@ -101,11 +99,9 @@ function DiracDelta(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
 end
 
 function DiracDelta(v::Vector{Point{D,T}},model::DiscreteModel{D}) where {D,T}
-  # check if the points are inside an active cell
   trian = Triangulation(model)
   cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
   cell = map(x -> _point_to_cell!(cache1, x), v)
-  # throws error if any Point not in domain
   point_grid = UnstructuredGrid(v)
   point_model = UnstructuredDiscreteModel(point_grid)
   point_trian = Triangulation(point_model)
diff --git a/test/CellDataTests/DiracDeltasTests.jl b/test/CellDataTests/DiracDeltasTests.jl
index f2732cc25..3875904ff 100644
--- a/test/CellDataTests/DiracDeltasTests.jl
+++ b/test/CellDataTests/DiracDeltasTests.jl
@@ -47,7 +47,7 @@ pvec = [p,π*p]
 δ = DiracDelta(pvec,model)
 @test sum(δ(v)) == sum(v(pvec))
 
-# below passes but cannot import FESpace and FEFunction here
+# below passes but cannot use FESpace and FEFunction here
 # reffe = ReferenceFE(lagrangian,Float64,3)
 # V = FESpace(model,reffe,conformity=:L2)
 # uh = FEFunction(V,rand(num_free_dofs(V)))

From 2a4921b473fe113b9d4b6168b8728a7a708f7bc6 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 11:03:35 +1100
Subject: [PATCH 08/17] removed Grid support attribute from DiracDelta, using
 just template parameter

---
 src/CellData/DiracDeltas.jl | 15 +++++----------
 1 file changed, 5 insertions(+), 10 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 04e43c1d5..29b9b2235 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -2,14 +2,9 @@ abstract type DiracDeltaSupportType <: GridapType end
 struct IsGridEntity <: DiracDeltaSupportType end
 struct NotGridEntity <: DiracDeltaSupportType end
 
-const _is_grid_entity = IsGridEntity()
-const _not_grid_entity = NotGridEntity()
-
-
 struct DiracDelta{D,S<:DiracDeltaSupportType} <: GridapType
   Γ::Triangulation{D}
   dΓ::Measure
-  supp::S
 end
 
 function DiracDelta{D}(
@@ -32,7 +27,7 @@ function DiracDelta{D}(
   trian = BodyFittedTriangulation(model,face_grid,face_to_bgface)
   Γ = BoundaryTriangulation(trian,glue)
   dΓ = Measure(Γ,degree)
-  DiracDelta(Γ,dΓ,_is_grid_entity)
+  DiracDelta{D,IsGridEntity}(Γ,dΓ)
 end
 
 function DiracDelta{D}(
@@ -95,7 +90,7 @@ function DiracDelta(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
   point_model = UnstructuredDiscreteModel(point_grid)
   point_trian = Triangulation(point_model)
   dx = Measure(point_trian,1)
-  DiracDelta(point_trian,dx,_not_grid_entity)
+  DiracDelta{0,NotGridEntity}(point_trian,dx)
 end
 
 function DiracDelta(v::Vector{Point{D,T}},model::DiscreteModel{D}) where {D,T}
@@ -106,11 +101,11 @@ function DiracDelta(v::Vector{Point{D,T}},model::DiscreteModel{D}) where {D,T}
   point_model = UnstructuredDiscreteModel(point_grid)
   point_trian = Triangulation(point_model)
   dx = Measure(point_trian,1)
-  DiracDelta(point_trian,dx,_not_grid_entity)
+  DiracDelta{0,NotGridEntity}(point_trian,dx)
 end
 
 function evaluate!(cache,d::DiracDelta{0,NotGridEntity},f::CellField)
-  eval = lazy_map(f,d.Γ.grid.node_coordinates)
+  d_f = lazy_map(f,d.Γ.grid.node_coordinates)
   dc = DomainContribution()
-  add_contribution!(dc, d.Γ, eval)
+  add_contribution!(dc, d.Γ, d_f)
 end

From b99d484ca479dcf9fd16d2b03a640fcaa7f00a5c Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 18:36:15 +1100
Subject: [PATCH 09/17] added Generic point DiracDelta tests for FEFunction and
 FEBasis

---
 test/CellDataTests/DiracDeltasTests.jl | 38 ++++++++++++++++++++------
 1 file changed, 30 insertions(+), 8 deletions(-)

diff --git a/test/CellDataTests/DiracDeltasTests.jl b/test/CellDataTests/DiracDeltasTests.jl
index 3875904ff..51a02fb96 100644
--- a/test/CellDataTests/DiracDeltasTests.jl
+++ b/test/CellDataTests/DiracDeltasTests.jl
@@ -8,6 +8,10 @@ using Gridap.ReferenceFEs
 using Gridap.Geometry
 using Gridap.CellData
 
+using Gridap.Fields
+using Gridap.FESpaces
+using Gridap.ReferenceFEs
+
 domain = (0,1,0,1)
 cells = (30,30)
 model = CartesianDiscreteModel(domain,cells)
@@ -39,19 +43,37 @@ degree = 2
 
 # Tests for DiracDelta at a generic Point in the domain #
 
-p = VectorValue(0.2,0.3)
+p = Point(0.2,0.3)
+δ = DiracDelta(p,model)
+@test sum(δ(v)) ≈ v(p)
+
+pvec = [p,π*p]
+δ = DiracDelta(pvec,model)
+@test sum(δ(v)) ≈ sum(v(pvec))
+
+p = Point(0.2,0.3)
 δ = DiracDelta(p,model)
-@test sum(δ(v)) == v(p)
+reffe = ReferenceFE(lagrangian,Float64,3)
+V = FESpace(model,reffe,conformity=:L2)
+V0 = TestFESpace(model,reffe,conformity=:H1,dirichlet_tags="boundary")
+U0 = TrialFESpace(V0)
+uh = FEFunction(U0,rand(num_free_dofs(U0)))
+vh = FEFunction(V,rand(num_free_dofs(V)))
+fe_basis = get_fe_basis(U0)
+dg_basis = get_fe_basis(V)
+
+@test sum(δ(uh)) ≈ uh(p)
+@test sum(δ(vh)) ≈ vh(p)
+@test norm(sum(δ(fe_basis)) .- fe_basis(p)) ≈ 0
+@test norm(sum(δ(dg_basis)) .- dg_basis(p)) ≈ 0
 
 pvec = [p,π*p]
 δ = DiracDelta(pvec,model)
-@test sum(δ(v)) == sum(v(pvec))
 
-# below passes but cannot use FESpace and FEFunction here
-# reffe = ReferenceFE(lagrangian,Float64,3)
-# V = FESpace(model,reffe,conformity=:L2)
-# uh = FEFunction(V,rand(num_free_dofs(V)))
-# @test sum(d(uh)) == uh(p)
+@test sum(δ(uh)) ≈ sum(map(uh,pvec))
+@test sum(δ(vh)) ≈ sum(map(vh,pvec))
+@test norm(sum(δ(fe_basis)) .- sum(map(fe_basis,pvec))) ≈ 0
+@test norm(sum(δ(dg_basis)) .- sum(map(dg_basis,pvec))) ≈ 0
 
 #using Gridap
 #

From cf74c49b65faade66d3c2960c0f7eb904fa70147 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Thu, 13 Oct 2022 18:38:08 +1100
Subject: [PATCH 10/17] Mentioned info on `lastindex` missing tests added later

---
 NEWS.md | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/NEWS.md b/NEWS.md
index d443ddbdd..a675316ac 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -10,7 +10,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - `lastindex` for `MultiValue`s for consistent usage of `[end]` as per `length`. Since PR [#834](https://github.com/gridap/Gridap.jl/pull/834)
 - BDM (Brezzi-Douglas-Marini) ReferenceFEs in PR [#823](https://github.com/gridap/Gridap.jl/pull/823)
 - Implemented `ConstantFESpace`. This space allows, e.g., to impose the mean value of the pressure via an additional unknown/equation (a Lagrange multiplier). Since PR [#836](https://github.com/gridap/Gridap.jl/pull/836)
-- Added support to construct `DiracDelta` at generic `Point`(s) in the domain. Since PR [#837](https://github.com/gridap/Gridap.jl/pull/837)
+- Methods to construct `DiracDelta` at generic `Point`(s) in the domain. Since PR [#837](https://github.com/gridap/Gridap.jl/pull/837)
+- some key missing `lastindex` and `end` tests belonging to PR [#834]. Since PR [#837](https://github.com/gridap/Gridap.jl/pull/837)
+
 ## [0.17.14] - 2022-07-29
 
 ### Added

From 69a727ed25e75e6e35afa79bc1f32bfdac12d013 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Sun, 16 Oct 2022 16:33:13 +1100
Subject: [PATCH 11/17] fully working implementation of DiracDelta at random
 points

---
 src/CellData/DiracDeltas.jl | 80 +++++++++++++++++++++++--------------
 src/Geometry/Geometry.jl    |  1 +
 2 files changed, 52 insertions(+), 29 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 29b9b2235..1fd847268 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -1,12 +1,14 @@
 abstract type DiracDeltaSupportType <: GridapType end
-struct IsGridEntity <: DiracDeltaSupportType end
-struct NotGridEntity <: DiracDeltaSupportType end
+abstract type IsGridEntity <: DiracDeltaSupportType end
+abstract type NotGridEntity <: DiracDeltaSupportType end
 
-struct DiracDelta{D,S<:DiracDeltaSupportType} <: GridapType
-  Γ::Triangulation{D}
+struct GenericDiracDelta{D,Dt,S<:DiracDeltaSupportType} <: GridapType
+  Γ::Triangulation{Dt}
   dΓ::Measure
 end
 
+const DiracDelta{D} = GenericDiracDelta{D,D,IsGridEntity}
+
 function DiracDelta{D}(
   model::DiscreteModel,
   face_to_bgface::AbstractVector{<:Integer},
@@ -27,7 +29,7 @@ function DiracDelta{D}(
   trian = BodyFittedTriangulation(model,face_grid,face_to_bgface)
   Γ = BoundaryTriangulation(trian,glue)
   dΓ = Measure(Γ,degree)
-  DiracDelta{D,IsGridEntity}(Γ,dΓ)
+  DiracDelta{D}(Γ,dΓ)
 end
 
 function DiracDelta{D}(
@@ -68,44 +70,64 @@ function DiracDelta{0}(model::DiscreteModel;tags)
   DiracDelta{0}(model,degree;tags=tags)
 end
 
-function (d::DiracDelta)(f)
+function (d::GenericDiracDelta)(f)
   evaluate(d,f)
 end
 
-function evaluate!(cache,d::DiracDelta,f)
+function evaluate!(cache,d::GenericDiracDelta,f)
   ∫(f)*d.dΓ
 end
 
-function get_triangulation(d::DiracDelta)
+function get_triangulation(d::GenericDiracDelta)
   d.Γ
 end
 
 # For handling DiracDelta at a generic Point in the domain #
 
-function DiracDelta(x::Point{D,T}, model::DiscreteModel{D}) where {D,T}
-  trian = Triangulation(model)
-  cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
-  cell = _point_to_cell!(cache1, x)
-  point_grid = UnstructuredGrid([x])
-  point_model = UnstructuredDiscreteModel(point_grid)
-  point_trian = Triangulation(point_model)
-  dx = Measure(point_trian,1)
-  DiracDelta{0,NotGridEntity}(point_trian,dx)
+function _cell_to_pindex(pvec::Vector{<:Point},trian::Triangulation)
+  cache = _point_to_cell_cache(KDTreeSearch(),trian)
+  cell_to_pindex = Dict{Int, Vector{Int32}}()
+  for i in 1:length(pvec)
+    cell = _point_to_cell!(cache, pvec[i])
+    push!(get!(() -> valtype(cell_to_pindex)[], cell_to_pindex, cell), i)
+  end
+  cell_to_pindex
 end
 
-function DiracDelta(v::Vector{Point{D,T}},model::DiscreteModel{D}) where {D,T}
+function DiracDelta(p::Point{D,T}, model::DiscreteModel{D}) where {D,T}
   trian = Triangulation(model)
-  cache1 = _point_to_cell_cache(KDTreeSearch(),trian)
-  cell = map(x -> _point_to_cell!(cache1, x), v)
-  point_grid = UnstructuredGrid(v)
-  point_model = UnstructuredDiscreteModel(point_grid)
-  point_trian = Triangulation(point_model)
-  dx = Measure(point_trian,1)
-  DiracDelta{0,NotGridEntity}(point_trian,dx)
+  cache = _point_to_cell_cache(KDTreeSearch(),trian)
+  cell = _point_to_cell!(cache, p)
+  trianv = TriangulationView(trian,[cell])
+  pquad = GenericQuadrature([p],[one(T)])
+  pmeas = Measure(CellQuadrature([pquad],[[p]],[[one(T)]],trianv,PhysicalDomain(),PhysicalDomain()))
+  GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
 
-function evaluate!(cache,d::DiracDelta{0,NotGridEntity},f::CellField)
-  d_f = lazy_map(f,d.Γ.grid.node_coordinates)
-  dc = DomainContribution()
-  add_contribution!(dc, d.Γ, d_f)
+function DiracDelta(pvec::Vector{Point{D,T}}, model::DiscreteModel{D}) where {D,T}
+  trian = Triangulation(model)
+  cell_to_pindices = _cell_to_pindex(pvec,trian)
+  cell_ids = collect(keys(cell_to_pindices))
+  cell_points = collect(values(cell_to_pindices))
+  points = map(i->pvec[cell_points[i]], 1:length(cell_ids))
+  weights_x_cell = collect.(Fill.(one(T),length.(cell_points)))
+  pquad = map(i -> GenericQuadrature(pvec[cell_points[i]],weights_x_cell[i]), 1:length(cell_ids))
+  trianv = Triangulation(trian,cell_ids)
+  pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
+  GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
+
+# TO DO: point_to_cell search cache needs to reused somehow
+# for Finite Element Function build using a fixed FE basis
+# since the basis remain fixed for each cell unless p or h
+# adapted, we can pre-compute them for given DiracDelta and
+# store them in cache. But for this we have to diverge from
+# the current Struct and add more attributes - like cache
+
+# GenericQuadrature doesn't taken FillArrays only Vectors
+# would be good to add this, as it would be useful here
+# The question is if it creates a problem down somewhere?
+# ofcourse excluding the tests which taken into account the types
+# like test_quadrature functions
+# if type independent problems exist for this generalization then
+# this is a more important problem
diff --git a/src/Geometry/Geometry.jl b/src/Geometry/Geometry.jl
index bf6f5f16a..bfca9bd9c 100644
--- a/src/Geometry/Geometry.jl
+++ b/src/Geometry/Geometry.jl
@@ -102,6 +102,7 @@ export compress_ids
 export UnstructuredGridTopology
 
 export Triangulation
+export TriangulationView
 export get_reffes
 export get_cell_coordinates
 export get_cell_ref_coordinates

From d967ea6e5597915d8f465952d32b6055c42a3509 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Mon, 17 Oct 2022 11:04:58 +1100
Subject: [PATCH 12/17] added DiracDelta tests for generic points

---
 src/CellData/DiracDeltas.jl            | 19 +-----
 test/CellDataTests/DiracDeltasTests.jl | 80 ++++++++++++++++++++++++--
 2 files changed, 78 insertions(+), 21 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 1fd847268..421102689 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -94,7 +94,7 @@ function _cell_to_pindex(pvec::Vector{<:Point},trian::Triangulation)
   cell_to_pindex
 end
 
-function DiracDelta(p::Point{D,T}, model::DiscreteModel{D}) where {D,T}
+function DiracDelta(model::DiscreteModel{D}, p::Point{D,T}) where {D,T}
   trian = Triangulation(model)
   cache = _point_to_cell_cache(KDTreeSearch(),trian)
   cell = _point_to_cell!(cache, p)
@@ -104,7 +104,7 @@ function DiracDelta(p::Point{D,T}, model::DiscreteModel{D}) where {D,T}
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
 
-function DiracDelta(pvec::Vector{Point{D,T}}, model::DiscreteModel{D}) where {D,T}
+function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,T}
   trian = Triangulation(model)
   cell_to_pindices = _cell_to_pindex(pvec,trian)
   cell_ids = collect(keys(cell_to_pindices))
@@ -116,18 +116,3 @@ function DiracDelta(pvec::Vector{Point{D,T}}, model::DiscreteModel{D}) where {D,
   pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
-
-# TO DO: point_to_cell search cache needs to reused somehow
-# for Finite Element Function build using a fixed FE basis
-# since the basis remain fixed for each cell unless p or h
-# adapted, we can pre-compute them for given DiracDelta and
-# store them in cache. But for this we have to diverge from
-# the current Struct and add more attributes - like cache
-
-# GenericQuadrature doesn't taken FillArrays only Vectors
-# would be good to add this, as it would be useful here
-# The question is if it creates a problem down somewhere?
-# ofcourse excluding the tests which taken into account the types
-# like test_quadrature functions
-# if type independent problems exist for this generalization then
-# this is a more important problem
diff --git a/test/CellDataTests/DiracDeltasTests.jl b/test/CellDataTests/DiracDeltasTests.jl
index 51a02fb96..3a5057744 100644
--- a/test/CellDataTests/DiracDeltasTests.jl
+++ b/test/CellDataTests/DiracDeltasTests.jl
@@ -8,6 +8,7 @@ using Gridap.ReferenceFEs
 using Gridap.Geometry
 using Gridap.CellData
 
+using Gridap.Algebra
 using Gridap.Fields
 using Gridap.FESpaces
 using Gridap.ReferenceFEs
@@ -44,15 +45,15 @@ degree = 2
 # Tests for DiracDelta at a generic Point in the domain #
 
 p = Point(0.2,0.3)
-δ = DiracDelta(p,model)
+δ = DiracDelta(model,p)
 @test sum(δ(v)) ≈ v(p)
 
 pvec = [p,π*p]
-δ = DiracDelta(pvec,model)
+δ = DiracDelta(model,pvec)
 @test sum(δ(v)) ≈ sum(v(pvec))
 
 p = Point(0.2,0.3)
-δ = DiracDelta(p,model)
+δ = DiracDelta(model,p)
 reffe = ReferenceFE(lagrangian,Float64,3)
 V = FESpace(model,reffe,conformity=:L2)
 V0 = TestFESpace(model,reffe,conformity=:H1,dirichlet_tags="boundary")
@@ -68,13 +69,84 @@ dg_basis = get_fe_basis(V)
 @test norm(sum(δ(dg_basis)) .- dg_basis(p)) ≈ 0
 
 pvec = [p,π*p]
-δ = DiracDelta(pvec,model)
+δ = DiracDelta(model,pvec)
 
 @test sum(δ(uh)) ≈ sum(map(uh,pvec))
 @test sum(δ(vh)) ≈ sum(map(vh,pvec))
 @test norm(sum(δ(fe_basis)) .- sum(map(fe_basis,pvec))) ≈ 0
 @test norm(sum(δ(dg_basis)) .- sum(map(dg_basis,pvec))) ≈ 0
 
+# comparing exact solution with that of GenericDiracDelta
+
+# 1D Test
+model = CartesianDiscreteModel((-1.,1.),(2,))
+Ω = Triangulation(model)
+dΩ = Measure(Ω,2)
+
+# exact solution
+function u(x)
+  if x[1] < 0.0
+    return 0.0
+  else
+    return -x[1]
+  end
+end
+
+ucf = CellField(u,Ω)
+
+order = 1
+reffe = ReferenceFE(lagrangian,Float64,order)
+V = FESpace(model,reffe,conformity=:L2)
+V0 = TestFESpace(model,reffe,conformity=:H1,dirichlet_tags="boundary")
+U0 = TrialFESpace(V0,u)
+
+p = Point(0.0)
+
+δ_p = DiracDelta(model,p)
+δ_tag = DiracDelta{0}(model,tags=3)
+
+a(u,v) = ∫(∇(u)⋅∇(v))*dΩ
+l(v) = ∫(0.0*v)*dΩ + δ_p(1.0*v)
+
+op = AffineFEOperator(a,l,U0,V0)
+uh_p = solve(op)
+
+l(v) = ∫(0.0*v)*dΩ + δ_tag(1.0*v)
+op = AffineFEOperator(a,l,U0,V0)
+uh_tag = solve(op)
+
+e_p = u - uh_p
+e = uh_p - uh_tag
+err_p = ∫(e_p*e_p)*dΩ
+err = ∫(e*e)*dΩ
+@test sum(err_p) < eps()
+@test sum(err) < eps()
+
+# 2D Test - comparing with tag version and point version
+
+model = CartesianDiscreteModel((-1.,1.,-1.,1.),(3,3))
+Ω = Triangulation(model)
+dΩ = Measure(Ω,2)
+
+pvec = [Point(-1.0/3,-1.0/3), Point(-1.0/3,1.0/3), Point(1.0/3,1.0/3), Point(1.0/3,-1.0/3)]
+δ_p = DiracDelta(model,pvec)
+
+a(u,v) = ∫(∇(u)⋅∇(v))*dΩ
+l(v) = δ_p(v)
+V0 = TestFESpace(model,reffe,conformity=:H1,dirichlet_tags="boundary")
+U0 = TrialFESpace(V0,u)
+op = AffineFEOperator(a,l,U0,V0)
+uh_p = solve(op)
+
+δ_tag = DiracDelta{0}(model,tags=9)
+l(v) = δ_tag(v)
+op = AffineFEOperator(a,l,U0,V0)
+uh_tag = solve(op)
+
+err = uh_p - uh_tag
+@test sum(∫(err*err)*dΩ) < eps()
+
+
 #using Gridap
 #
 #order = 2

From 633326e02bd8aec4d64bd0a9e757c72dfd2c4f7e Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Mon, 17 Oct 2022 12:26:35 +1100
Subject: [PATCH 13/17] added a specialized evaluate for DiracDelta with
 Function

---
 src/CellData/DiracDeltas.jl | 18 ++++++++++++++++--
 1 file changed, 16 insertions(+), 2 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 421102689..ce2f6f9fc 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -84,7 +84,7 @@ end
 
 # For handling DiracDelta at a generic Point in the domain #
 
-function _cell_to_pindex(pvec::Vector{<:Point},trian::Triangulation)
+function _cell_to_pindices(pvec::Vector{<:Point},trian::Triangulation)
   cache = _point_to_cell_cache(KDTreeSearch(),trian)
   cell_to_pindex = Dict{Int, Vector{Int32}}()
   for i in 1:length(pvec)
@@ -106,7 +106,7 @@ end
 
 function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,T}
   trian = Triangulation(model)
-  cell_to_pindices = _cell_to_pindex(pvec,trian)
+  cell_to_pindices = _cell_to_pindices(pvec,trian)
   cell_ids = collect(keys(cell_to_pindices))
   cell_points = collect(values(cell_to_pindices))
   points = map(i->pvec[cell_points[i]], 1:length(cell_ids))
@@ -116,3 +116,17 @@ function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,
   pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
+
+function evaluate!(cache,d::GenericDiracDelta{0,Dt,NotGridEntity},f::Function) where Dt
+  dc = DomainContribution()
+  quad_points = d.dΓ.quad.cell_point
+  weights = d.dΓ.quad.cell_weight
+  dc_Γ = zeros(eltype(eltype(weights)), length(weights))
+  for i in eachindex(weights)
+    for j in eachindex(weights[i])
+      @inbounds dc_Γ[i] += f(quad_points[i][j])*weights[i][j]
+    end
+  end
+  add_contribution!(dc,d.Γ,dc_Γ)
+  dc
+end

From d3f26fb45288b7020c569516627ba7dd57faf1fb Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Mon, 17 Oct 2022 15:20:17 +1100
Subject: [PATCH 14/17] made changes for memory reuse and added tests

---
 src/CellData/DiracDeltas.jl            | 8 +++++---
 test/CellDataTests/DiracDeltasTests.jl | 6 ++++++
 2 files changed, 11 insertions(+), 3 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index ce2f6f9fc..825dd5ec7 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -99,8 +99,10 @@ function DiracDelta(model::DiscreteModel{D}, p::Point{D,T}) where {D,T}
   cache = _point_to_cell_cache(KDTreeSearch(),trian)
   cell = _point_to_cell!(cache, p)
   trianv = TriangulationView(trian,[cell])
-  pquad = GenericQuadrature([p],[one(T)])
-  pmeas = Measure(CellQuadrature([pquad],[[p]],[[one(T)]],trianv,PhysicalDomain(),PhysicalDomain()))
+  point = [p]
+  weight = [one(T)]
+  pquad = GenericQuadrature(point,weight)
+  pmeas = Measure(CellQuadrature([pquad],[pquad.coordinates],[pquad.weights],trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
 
@@ -111,7 +113,7 @@ function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,
   cell_points = collect(values(cell_to_pindices))
   points = map(i->pvec[cell_points[i]], 1:length(cell_ids))
   weights_x_cell = collect.(Fill.(one(T),length.(cell_points)))
-  pquad = map(i -> GenericQuadrature(pvec[cell_points[i]],weights_x_cell[i]), 1:length(cell_ids))
+  pquad = map(i -> GenericQuadrature(points[i],weights_x_cell[i]), 1:length(cell_ids))
   trianv = Triangulation(trian,cell_ids)
   pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
diff --git a/test/CellDataTests/DiracDeltasTests.jl b/test/CellDataTests/DiracDeltasTests.jl
index 3a5057744..16501a1ed 100644
--- a/test/CellDataTests/DiracDeltasTests.jl
+++ b/test/CellDataTests/DiracDeltasTests.jl
@@ -146,6 +146,12 @@ uh_tag = solve(op)
 err = uh_p - uh_tag
 @test sum(∫(err*err)*dΩ) < eps()
 
+# testing the faster and direct functional evalulation method
+
+f(x) = sin(norm(x))
+fcf = CellField(f,Ω)
+@test sum(δ_p(f)) ≈ sum(δ_p(fcf))
+
 
 #using Gridap
 #

From c5d38b15f2480be977c244067cb453f3a4c208db Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Mon, 17 Oct 2022 15:29:12 +1100
Subject: [PATCH 15/17] mistake of having Triangulation instead of View (was
 still working)

---
 src/CellData/DiracDeltas.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 825dd5ec7..a3ed844d2 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -114,7 +114,7 @@ function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,
   points = map(i->pvec[cell_points[i]], 1:length(cell_ids))
   weights_x_cell = collect.(Fill.(one(T),length.(cell_points)))
   pquad = map(i -> GenericQuadrature(points[i],weights_x_cell[i]), 1:length(cell_ids))
-  trianv = Triangulation(trian,cell_ids)
+  trianv = TriangulationView(trian,cell_ids)
   pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end

From b3b867c41b14827ad21c1afb10f24690a4fc5e9a Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Mon, 17 Oct 2022 17:14:11 +1100
Subject: [PATCH 16/17] removed specialized evaluate! method for Function

Need to analyze why this is much faster than the usual way and
if we can generalize it to Gridap to make it faster, for such
situations involving Functions
---
 src/CellData/DiracDeltas.jl | 14 --------------
 1 file changed, 14 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index a3ed844d2..83825f663 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -118,17 +118,3 @@ function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,
   pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
-
-function evaluate!(cache,d::GenericDiracDelta{0,Dt,NotGridEntity},f::Function) where Dt
-  dc = DomainContribution()
-  quad_points = d.dΓ.quad.cell_point
-  weights = d.dΓ.quad.cell_weight
-  dc_Γ = zeros(eltype(eltype(weights)), length(weights))
-  for i in eachindex(weights)
-    for j in eachindex(weights[i])
-      @inbounds dc_Γ[i] += f(quad_points[i][j])*weights[i][j]
-    end
-  end
-  add_contribution!(dc,d.Γ,dc_Γ)
-  dc
-end

From 4b841f13c6e7cabee29013a7b0c3cb7a0d9152b9 Mon Sep 17 00:00:00 2001
From: Kishore Nori <sai.nori@monash.edu>
Date: Mon, 17 Oct 2022 17:22:36 +1100
Subject: [PATCH 17/17] Replaced `TriangulationView` by `view`

---
 src/CellData/DiracDeltas.jl | 4 ++--
 src/Geometry/Geometry.jl    | 1 -
 2 files changed, 2 insertions(+), 3 deletions(-)

diff --git a/src/CellData/DiracDeltas.jl b/src/CellData/DiracDeltas.jl
index 83825f663..db9cf9946 100644
--- a/src/CellData/DiracDeltas.jl
+++ b/src/CellData/DiracDeltas.jl
@@ -98,7 +98,7 @@ function DiracDelta(model::DiscreteModel{D}, p::Point{D,T}) where {D,T}
   trian = Triangulation(model)
   cache = _point_to_cell_cache(KDTreeSearch(),trian)
   cell = _point_to_cell!(cache, p)
-  trianv = TriangulationView(trian,[cell])
+  trianv = view(trian,[cell])
   point = [p]
   weight = [one(T)]
   pquad = GenericQuadrature(point,weight)
@@ -114,7 +114,7 @@ function DiracDelta(model::DiscreteModel{D}, pvec::Vector{Point{D,T}}) where {D,
   points = map(i->pvec[cell_points[i]], 1:length(cell_ids))
   weights_x_cell = collect.(Fill.(one(T),length.(cell_points)))
   pquad = map(i -> GenericQuadrature(points[i],weights_x_cell[i]), 1:length(cell_ids))
-  trianv = TriangulationView(trian,cell_ids)
+  trianv = view(trian,cell_ids)
   pmeas = Measure(CellQuadrature(pquad,points,weights_x_cell,trianv,PhysicalDomain(),PhysicalDomain()))
   GenericDiracDelta{0,D,NotGridEntity}(trianv,pmeas)
 end
diff --git a/src/Geometry/Geometry.jl b/src/Geometry/Geometry.jl
index bfca9bd9c..bf6f5f16a 100644
--- a/src/Geometry/Geometry.jl
+++ b/src/Geometry/Geometry.jl
@@ -102,7 +102,6 @@ export compress_ids
 export UnstructuredGridTopology
 
 export Triangulation
-export TriangulationView
 export get_reffes
 export get_cell_coordinates
 export get_cell_ref_coordinates