xref: /libCEED/julia/LibCEED.jl/test/buildmats.jl (revision 2efebffe38fa9227caaeab504b43e3a698cb86d7)
1*b2e3f8ecSSebastian Grimbergfunction build_mats_hdiv(qref, qweight, ::Type{T}=Float64) where {T}
2*b2e3f8ecSSebastian Grimberg    P, Q, dim = 4, 4, 2
3*b2e3f8ecSSebastian Grimberg    interp = Array{T}(undef, dim, Q, P)
4*b2e3f8ecSSebastian Grimberg    div = Array{T}(undef, Q, P)
5*b2e3f8ecSSebastian Grimberg
6*b2e3f8ecSSebastian Grimberg    qref[1, 1] = -1.0/sqrt(3.0)
7*b2e3f8ecSSebastian Grimberg    qref[1, 2] = qref[1, 1]
8*b2e3f8ecSSebastian Grimberg    qref[1, 3] = qref[1, 1]
9*b2e3f8ecSSebastian Grimberg    qref[1, 4] = -qref[1, 1]
10*b2e3f8ecSSebastian Grimberg    qref[2, 1] = -qref[1, 1]
11*b2e3f8ecSSebastian Grimberg    qref[2, 2] = -qref[1, 1]
12*b2e3f8ecSSebastian Grimberg    qref[2, 3] = qref[1, 1]
13*b2e3f8ecSSebastian Grimberg    qref[2, 4] = qref[1, 1]
14*b2e3f8ecSSebastian Grimberg    qweight[1] = 1.0
15*b2e3f8ecSSebastian Grimberg    qweight[2] = 1.0
16*b2e3f8ecSSebastian Grimberg    qweight[3] = 1.0
17*b2e3f8ecSSebastian Grimberg    qweight[4] = 1.0
18*b2e3f8ecSSebastian Grimberg
19*b2e3f8ecSSebastian Grimberg    # Loop over quadrature points
20*b2e3f8ecSSebastian Grimberg    for i = 1:Q
21*b2e3f8ecSSebastian Grimberg        x1 = qref[1, i]
22*b2e3f8ecSSebastian Grimberg        x2 = qref[2, i]
23*b2e3f8ecSSebastian Grimberg        # Interp
24*b2e3f8ecSSebastian Grimberg        interp[1, i, 1] = 0.0
25*b2e3f8ecSSebastian Grimberg        interp[2, i, 1] = 1.0 - x2
26*b2e3f8ecSSebastian Grimberg        interp[1, i, 2] = x1 - 1.0
27*b2e3f8ecSSebastian Grimberg        interp[2, i, 2] = 0.0
28*b2e3f8ecSSebastian Grimberg        interp[1, i, 3] = -x1
29*b2e3f8ecSSebastian Grimberg        interp[2, i, 3] = 0.0
30*b2e3f8ecSSebastian Grimberg        interp[1, i, 4] = 0.0
31*b2e3f8ecSSebastian Grimberg        interp[2, i, 4] = x2
32*b2e3f8ecSSebastian Grimberg        # Div
33*b2e3f8ecSSebastian Grimberg        div[i, 1] = -1.0
34*b2e3f8ecSSebastian Grimberg        div[i, 2] = 1.0
35*b2e3f8ecSSebastian Grimberg        div[i, 3] = -1.0
36*b2e3f8ecSSebastian Grimberg        div[i, 4] = 1.0
37*b2e3f8ecSSebastian Grimberg    end
38*b2e3f8ecSSebastian Grimberg
39*b2e3f8ecSSebastian Grimberg    return interp, div
40*b2e3f8ecSSebastian Grimbergend
41*b2e3f8ecSSebastian Grimberg
42*b2e3f8ecSSebastian Grimbergfunction build_mats_hcurl(qref, qweight, ::Type{T}=Float64) where {T}
43*b2e3f8ecSSebastian Grimberg    P, Q, dim = 3, 4, 2
44*b2e3f8ecSSebastian Grimberg    interp = Array{T}(undef, dim, Q, P)
45*b2e3f8ecSSebastian Grimberg    curl = Array{T}(undef, 1, Q, P)
46*b2e3f8ecSSebastian Grimberg
47*b2e3f8ecSSebastian Grimberg    qref[1, 1] = 0.2
48*b2e3f8ecSSebastian Grimberg    qref[1, 2] = 0.6
49*b2e3f8ecSSebastian Grimberg    qref[1, 3] = 1.0/3.0
50*b2e3f8ecSSebastian Grimberg    qref[1, 4] = 0.2
51*b2e3f8ecSSebastian Grimberg    qref[2, 1] = 0.2
52*b2e3f8ecSSebastian Grimberg    qref[2, 2] = 0.2
53*b2e3f8ecSSebastian Grimberg    qref[2, 3] = 1.0/3.0
54*b2e3f8ecSSebastian Grimberg    qref[2, 4] = 0.6
55*b2e3f8ecSSebastian Grimberg    qweight[1] = 25.0/96.0
56*b2e3f8ecSSebastian Grimberg    qweight[2] = 25.0/96.0
57*b2e3f8ecSSebastian Grimberg    qweight[3] = -27.0/96.0
58*b2e3f8ecSSebastian Grimberg    qweight[4] = 25.0/96.0
59*b2e3f8ecSSebastian Grimberg
60*b2e3f8ecSSebastian Grimberg    # Loop over quadrature points
61*b2e3f8ecSSebastian Grimberg    for i = 1:Q
62*b2e3f8ecSSebastian Grimberg        x1 = qref[1, i]
63*b2e3f8ecSSebastian Grimberg        x2 = qref[2, i]
64*b2e3f8ecSSebastian Grimberg        # Interp
65*b2e3f8ecSSebastian Grimberg        interp[1, i, 1] = -x2
66*b2e3f8ecSSebastian Grimberg        interp[2, i, 1] = x1
67*b2e3f8ecSSebastian Grimberg        interp[1, i, 2] = x2
68*b2e3f8ecSSebastian Grimberg        interp[2, i, 2] = 1.0 - x1
69*b2e3f8ecSSebastian Grimberg        interp[1, i, 3] = 1.0 - x2
70*b2e3f8ecSSebastian Grimberg        interp[2, i, 3] = x1
71*b2e3f8ecSSebastian Grimberg        # Curl
72*b2e3f8ecSSebastian Grimberg        curl[1, i, 1] = 2.0
73*b2e3f8ecSSebastian Grimberg        curl[1, i, 2] = -2.0
74*b2e3f8ecSSebastian Grimberg        curl[1, i, 3] = -2.0
75*b2e3f8ecSSebastian Grimberg    end
76*b2e3f8ecSSebastian Grimberg
77*b2e3f8ecSSebastian Grimberg    return interp, curl
78*b2e3f8ecSSebastian Grimbergend
79