Merge pull request #771 from gridap/fix_issue_770

Fixes related with 2nd order derivatives on vector-valued fields
gridap · Mar 23, 2022 · e328bd0 · e328bd0
2 parents edf2f90 + 4aeb153
commit e328bd0
Show file tree

Hide file tree

Showing 5 changed files with 56 additions and 12 deletions.
diff --git a/src/Fields/AffineMaps.jl b/src/Fields/AffineMaps.jl
@@ -52,13 +52,29 @@ function push_∇∇(∇∇a::Field,ϕ::AffineMap)
   Operation(push_∇∇)(∇∇a, Jt_inv)
 end
 
-function push_∇∇(∇∇a::Number,Jt_inv::MultiValue{Tuple{D,D}} where D)
-  #Jt_inv⋅Jt_inv⋅∇∇a
-  Jt_inv⋅∇∇a⋅transpose(Jt_inv)
-end
+#function push_∇∇(∇∇a::Number,Jt_inv::MultiValue{Tuple{D,D}} where D)
+#  #Jt_inv⋅Jt_inv⋅∇∇a
+#  Jt_inv⋅∇∇a⋅transpose(Jt_inv)
+#end
 
 function push_∇∇(∇∇a::Number,Jt_inv::MultiValue{Tuple{D1,D2}} where {D1,D2})
-  Jt_inv⋅∇∇a⋅transpose(Jt_inv)
+  _permdims_for_∇∇(Jt_inv⋅_permdims_for_∇∇(∇∇a)⋅transpose(Jt_inv))
+end
+
+function _permdims_for_∇∇(a::MultiValue{Tuple{D1,D2}}) where {D1,D2}
+  a
+end
+@generated function _permdims_for_∇∇(a::MultiValue{Tuple{D1,D2,D3}}) where {D1,D2,D3}
+  ss = String[]
+  for k in 1:D2
+    for j in 1:D3
+      for i in 1:D1
+        push!(ss,"a[$i,$k,$j],")
+      end
+    end
+  end
+  str =  join(ss)
+  Meta.parse("ThirdOrderTensorValue{$D1,$D3,$D2}($str)")
 end
 
 function lazy_map(

diff --git a/src/Fields/FieldsInterfaces.jl b/src/Fields/FieldsInterfaces.jl
@@ -739,7 +739,7 @@ function _zero_size(a::VoidBasis{T,1} where T)
 end
 
 function _zero_size(a::VoidBasis{T,2} where T)
-  @check size(a,1) == 1
+  @check size(a,1) in (0,1)
   (1,0)
 end
 

diff --git a/src/TensorValues/Operations.jl b/src/TensorValues/Operations.jl
@@ -234,18 +234,18 @@ end
 end
 
 # a_ilm = b_ij*c_jlm
-@generated function dot(a::A,b::B) where {A<:MultiValue{Tuple{D,D}},B<:ThirdOrderTensorValue{D,D,L}} where {D,L}
+@generated function dot(a::A,b::B) where {A<:MultiValue{Tuple{D1,D2}},B<:ThirdOrderTensorValue{D2,D3,D4}} where {D1,D2,D3,D4}
   ss = String[]
-  for m in 1:L
-    for l in 1:D
-      for i in 1:D
-        s = join([ "a[$i,$j]*b[$j,$l,$m]+" for j in 1:D])
+  for m in 1:D4
+    for l in 1:D3
+      for i in 1:D1
+        s = join([ "a[$i,$j]*b[$j,$l,$m]+" for j in 1:D2])
         push!(ss,s[1:(end-1)]*", ")
       end
     end
   end
   str = join(ss)
-  Meta.parse("ThirdOrderTensorValue{$D,$D,$L}($str)")
+  Meta.parse("ThirdOrderTensorValue{$D1,$D3,$D4}($str)")
 end
 
 const ⋅¹ = dot

diff --git a/test/GridapTests/issue_770.jl b/test/GridapTests/issue_770.jl
@@ -0,0 +1,26 @@
+module Issue770
+
+using Gridap
+using Test
+
+domain = (0,1,0,1)
+cells = (2,2)
+model = CartesianDiscreteModel(domain,cells) |> simplexify
+Ω = Interior(model)
+dΩ = Measure(Ω,4)
+m = 2
+u(x) = VectorValue(sum(x)^m, sum(x)^m)
+reffe = ReferenceFE(lagrangian,VectorValue{2,Float64},m)
+V = FESpace(Ω,reffe)
+uh = interpolate(u,V)
+t = Δ(u) - Δ(uh)
+@test sum(∫(t⋅t)dΩ) < 1.0e-12
+
+s(x) = sum(x)^m
+reffe = ReferenceFE(lagrangian,Float64,m)
+V = FESpace(Ω,reffe)
+sh = interpolate(s,V)
+t = Δ(s) - Δ(sh)
+@test sum(∫(t⋅t)dΩ) < 1.0e-12
+
+end # module
diff --git a/test/GridapTests/runtests.jl b/test/GridapTests/runtests.jl
@@ -42,4 +42,6 @@ using Test
 
 @time @testset "Issue760" begin include("issue_760.jl") end
 
+@time @testset "Issue770" begin include("issue_770.jl") end
+
 end # module