-
Notifications
You must be signed in to change notification settings - Fork 99
/
Copy pathFESpacesWithLinearConstraintsTests.jl
104 lines (77 loc) · 2.48 KB
/
FESpacesWithLinearConstraintsTests.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
module FESpacesWithLinearConstraintsTests
using Gridap.Algebra
using Gridap.Arrays
using Gridap.Fields
using Gridap.Geometry
using Gridap.FESpaces
using Test
using LinearAlgebra
using Gridap.CellData
using Gridap.ReferenceFEs
domain = (0,1,0,1)
partition = (2,2)
model = CartesianDiscreteModel(domain,partition)
labels = get_face_labeling(model)
add_tag_from_tags!(labels,"dirichlet",[1,2,5])
add_tag_from_tags!(labels,"neumann",[6,7,8])
Ω = Triangulation(model)
Γ = BoundaryTriangulation(model,tags="neumann")
Λ = SkeletonTriangulation(model)
dΩ = Measure(Ω,2)
dΓ = Measure(Γ,2)
dΛ = Measure(Λ,2)
V = FESpace(
model,ReferenceFE(lagrangian,Float64,1), conformity=:H1, dirichlet_tags="dirichlet")
test_single_field_fe_space(V)
fdof_to_val = collect(Float64,1:num_free_dofs(V))
ddof_to_val = -collect(Float64,1:num_dirichlet_dofs(V))
vh = FEFunction(V,fdof_to_val,ddof_to_val)
sDOF_to_dof = [1,5,-2]
sDOF_to_dofs = Table([[-1,4],[4,6],[-1,-3]])
sDOF_to_coeffs = Table([[0.5,0.5],[0.5,0.5],[0.5,0.5]])
Vc = FESpaceWithLinearConstraints(
sDOF_to_dof,
sDOF_to_dofs,
sDOF_to_coeffs,
V)
test_single_field_fe_space(Vc)
@test has_constraints(Vc)
@test isa(get_cell_constraints(Vc,Λ)[1],ArrayBlock)
@test Vc.n_fdofs == 6
@test Vc.n_fmdofs == 4
fmdof_to_val = collect(Float64,1:num_free_dofs(Vc))
dmdof_to_val = -collect(Float64,1:num_dirichlet_dofs(Vc))
vch = FEFunction(Vc,fmdof_to_val,dmdof_to_val)
r = [[-1.0, -1.5, 1.0, 1.0], [-1.5, -2.0, 1.0, 2.0], [1.0, 1.0, 3.0, 3.5], [1.0, 2.0, 3.5, 4.0]]
@test get_cell_dof_values(vch) ≈ r
v(x) = sin(4*x[1]+0.4)*cos(5*x[2]+0.7)
vch = interpolate(v,Vc)
#using Gridap.Visualization
#writevtk(Ω,"trian",nsubcells=10,cellfields=["vh"=>vh,"vch"=>vch])
u(x) = x[1] + 2*x[2]
f(x) = -Δ(u)(x)
Uc = TrialFESpace(Vc,u)
@test has_constraints(Uc)
uch = interpolate(u,Uc)
n_Γ = get_normal_vector(Γ)
a(u,v) = ∫( ∇(v)⋅∇(u) )*dΩ + ∫( jump(u)*jump(v) )*dΛ
l(v) = ∫( v*f )*dΩ + ∫( v*(n_Γ⋅∇(u)) )*dΓ
op = AffineFEOperator(a,l,Uc,Vc)
uch = solve(op)
#using Gridap.Visualization
#writevtk(trian,"trian",nsubcells=10,cellfields=["uch"=>uch])
e = u - uch
e_l2 = sqrt(sum(∫( e*e )*dΩ))
e_h1 = sqrt(sum(∫( e*e + ∇(e)⋅∇(e) )*dΩ))
tol = 1.e-9
@test e_l2 < tol
@test e_h1 < tol
V2 = FESpace(
model,ReferenceFE(lagrangian,Float64,1), conformity=:H1, dirichlet_tags="dirichlet", vector_type=ComplexF64)
Vc2 = FESpaceWithLinearConstraints(
sDOF_to_dof,
sDOF_to_dofs,
sDOF_to_coeffs,
V2)
@test get_dof_value_type(Vc2) <: ComplexF64
end # module