Skip to content

Commit

Permalink
fixes
Browse files Browse the repository at this point in the history
  • Loading branch information
fredrikekre committed Oct 28, 2024
1 parent eb2b2a7 commit 4b21079
Show file tree
Hide file tree
Showing 4 changed files with 17 additions and 13 deletions.
12 changes: 8 additions & 4 deletions benchmark/benchmarks-dofs.jl
Original file line number Diff line number Diff line change
Expand Up @@ -39,23 +39,26 @@ for spatial_dim in [3] # 1:3
close_helper = function(grid, ip)
dh = DofHandler(grid)
push!(dh, :u, field_dim, ip)
return close!(dh)
close!(dh)
return
end
LAGRANGE_SUITE["DofHandler"]["one-field"] = @benchmarkable $close_helper($grid, $ip)

close_helper = function(grid, ip, ip2)
dh = DofHandler(grid)
push!(dh, :u, field_dim, ip)
push!(dh, :p, 1, ip2)
return close!(dh)
close!(dh)
return
end
LAGRANGE_SUITE["DofHandler"]["two-fields"] = @benchmarkable $close_helper($grid, $ip, $ip2)

close_helper = function(grid)
dh = DofHandler(grid)
sdh = SubDofHandler(dh, Set(1:Int(round(getncells(grid) / 2))))
add!(sdh, :u, ip^field_dim)
return close!(dh)
close!(dh)
return
end
LAGRANGE_SUITE["DofHandler"]["one-field-subdomain"] = @benchmarkable $close_helper($grid)

Expand All @@ -64,7 +67,8 @@ for spatial_dim in [3] # 1:3
sdh = SubDofHandler(dh, Set(1:Int(round(getncells(grid) / 2))))
add!(sdh, :u, ip^field_dim)
add!(sdh, :p, ip2)
return close!(dh)
close!(dh)
return
end
LAGRANGE_SUITE["DofHandler"]["two-fields-subdomain"] = @benchmarkable $close_helper($grid)
end
Expand Down
4 changes: 2 additions & 2 deletions docs/src/literate-tutorials/incompressible_elasticity.jl
Original file line number Diff line number Diff line change
Expand Up @@ -275,8 +275,8 @@ function solve(ν, interpolation_u, interpolation_p)

## Export the solution and the stress
filename = "cook_" *
(interpolation_u == Lagrange{RefTriangle, 1}()^2 ? "linear" : "quadratic") *
"_linear"
(interpolation_u == Lagrange{RefTriangle, 1}()^2 ? "linear" : "quadratic") *
"_linear"

VTKGridFile(filename, grid) do vtk
write_solution(vtk, dh, u)
Expand Down
4 changes: 2 additions & 2 deletions docs/src/literate-tutorials/linear_shell.jl
Original file line number Diff line number Diff line change
Expand Up @@ -154,10 +154,10 @@ function fiber_coordsys(Ps::Vector{Vec{3, Float64}})
a = abs.(P)
j = 1
if a[1] > a[3]
a[3] = a[1]; j = 2
a[3] = a[1]; j = 2
end
if a[2] > a[3]
j = 3
j = 3
end

e3 = P
Expand Down
10 changes: 5 additions & 5 deletions src/interpolations.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1283,7 +1283,7 @@ function reference_shape_value(ip::Lagrange{RefPyramid, 2}, ξ::Vec{3, T}, i::In
i == 2 && return zzero ? zero(T) : x * (4x ** (z - 1) + y * (6x + 2y - z) * (z² - 2z + 1) + (z - 1) * (2x + 3y - 1) * (z² - 2z + 1)) / ((z - 1) * (z² - 2z + 1))
i == 3 && return zzero ? zero(T) : y * (4* y * (z - 1) + x * (2x + 6y - z) * (z² - 2z + 1) + (z - 1) * (3x + 2y - 1) * (z² - 2z + 1)) / ((z - 1) * (z² - 2z + 1))
i == 4 && return zzero ? zero(T) : x * y * (4 * x * y + 2x * z - 2x + 2y * z - 2y + 2- 3z + 1) / (z² - 2z + 1)
i == 5 && return z * (2z - 1)
i == 5 && return z * (2z - 1)
i == 6 && return zzero ? zero(T) : 4x * (2x ** (1 - z) - y * (3x + 2y) * (z² - 2z + 1) + (z - 1) * (z² - 2z + 1) * (-x - 3y - z + 1)) / ((z - 1) * (z² - 2z + 1))
i == 7 && return zzero ? zero(T) : 4y * (2* y * (1 - z) - x * (2x + 3y) * (z² - 2z + 1) + (z - 1) * (z² - 2z + 1) * (-3x - y - z + 1)) / ((z - 1) * (z² - 2z + 1))
i == 8 && return zzero ? zero(T) : 4z * (-x * y + (z - 1) * (-x - y - z + 1)) / (z - 1)
Expand Down Expand Up @@ -1599,9 +1599,9 @@ function reference_shape_value(ip::RannacherTurek{RefQuadrilateral, 1}, ξ::Vec{
(x, y) = ξ

i == 1 && return -(x + 1)^2 / 4 + (y + 1)^2 / 4 + (x + 1) / 2 - (y + 1) + T(3) / 4
i == 2 && return (x + 1)^2 / 4 - (y + 1)^2 / 4 + (y + 1) / 2 - T(1) / 4
i == 2 && return (x + 1)^2 / 4 - (y + 1)^2 / 4 + (y + 1) / 2 - T(1) / 4
i == 3 && return -(x + 1)^2 / 4 + (y + 1)^2 / 4 + (x + 1) / 2 - T(1) / 4
i == 4 && return (x + 1)^2 / 4 - (y + 1)^2 / 4 - (x + 1) + (y + 1) / 2 + T(3) / 4
i == 4 && return (x + 1)^2 / 4 - (y + 1)^2 / 4 - (x + 1) + (y + 1) / 2 + T(3) / 4
throw(ArgumentError("no shape function $i for interpolation $ip"))
end

Expand Down Expand Up @@ -1631,9 +1631,9 @@ function reference_shape_value(ip::RannacherTurek{RefHexahedron, 1}, ξ::Vec{3,

i == 1 && return -2((x + 1))^2 / 12 + 1(x + 1) / 3 - 2((y + 1))^2 / 12 + 1(y + 1) / 3 + 4((z + 1))^2 / 12 - 7(z + 1) / 6 + T(2) / 3
i == 2 && return -2((x + 1))^2 / 12 + 1(x + 1) / 3 + 4((y + 1))^2 / 12 - 7(y + 1) / 6 - 2((z + 1))^2 / 12 + 1(z + 1) / 3 + T(2) / 3
i == 3 && return 4((x + 1))^2 / 12 - 1(x + 1) / 6 - 2((y + 1))^2 / 12 + 1(y + 1) / 3 - 2((z + 1))^2 / 12 + 1(z + 1) / 3 - T(1) / 3
i == 3 && return 4((x + 1))^2 / 12 - 1(x + 1) / 6 - 2((y + 1))^2 / 12 + 1(y + 1) / 3 - 2((z + 1))^2 / 12 + 1(z + 1) / 3 - T(1) / 3
i == 4 && return -2((x + 1))^2 / 12 + 1(x + 1) / 3 + 4((y + 1))^2 / 12 - 1(y + 1) / 6 - 2((z + 1))^2 / 12 + 1(z + 1) / 3 - T(1) / 3
i == 5 && return 4((x + 1))^2 / 12 - 7(x + 1) / 6 - 2((y + 1))^2 / 12 + 1(y + 1) / 3 - 2((z + 1))^2 / 12 + 1(z + 1) / 3 + T(2) / 3
i == 5 && return 4((x + 1))^2 / 12 - 7(x + 1) / 6 - 2((y + 1))^2 / 12 + 1(y + 1) / 3 - 2((z + 1))^2 / 12 + 1(z + 1) / 3 + T(2) / 3
i == 6 && return -2((x + 1))^2 / 12 + 1(x + 1) / 3 - 2((y + 1))^2 / 12 + 1(y + 1) / 3 + 4((z + 1))^2 / 12 - 1(z + 1) / 6 - T(1) / 3

throw(ArgumentError("no shape function $i for interpolation $ip"))
Expand Down

0 comments on commit 4b21079

Please sign in to comment.