In [ ]:

    

include("../src/HPFEM.jl")


    

In [ ]:

    

nel = 10
nnodes = nel + 1
idir = [1]#,nnodes]
M = 9
Q = M+2
bas = HPFEM.Basis1d(M,Q)
lmap = HPFEM.locmap(bas)
dof = HPFEM.DofMap1d(lmap, nnodes, idir);


    

In [ ]:

    

fun(x) = cos(5π*x)


    

In [ ]:

    

a = -0.1
b = 2.6
nodes = collect(linspace(a, b, nnodes))


    

In [ ]:

    

elems = [HPFEM.Element1d(e, nodes[e], nodes[e+1], bas) for e = 1:nel];


    

In [ ]:

    

solver = HPFEM.CholeskySC(dof, HPFEM.BBSymTri);


    

In [ ]:

    

for e = 1:nel
    Ae = HPFEM.mass_matrix(bas, elems[e])
    HPFEM.add_local_matrix(solver, e, Ae)
end


    

In [ ]:

    

Fe = zeros(HPFEM.nmodes(lmap), nel)

for e = 1:nel
    fe = fun(elems[e].x)
    HPFEM.add_rhs!(bas, elems[e], fe, view(Fe, :, e))
end

# Apply Dirichilet BCs:
Fe[1,1] = fun(a)
#Fe[2,2] = cos(b)
Fe


    

In [ ]:

    

HPFEM.solve!(solver, Fe)


    

In [ ]:

    

nξ = 101
ξ = collect(linspace(-1,1,nξ));
ϕ = zeros(nξ, M)
for i = 1:M
    ϕ[:,i] = HPFEM.basis1d(bas,ξ, i)
end

Ue = ϕ * Fe;


    

In [ ]:

    

using PyPlot
x = [(1-ξ)*el.a/2 + (1+ξ)*el.b/2 for el in elems]
maxerr = 0.0
for e = 1:nel
    uu = fun(x[e])
    err = maxabs(uu-Ue[:,e])
    if err > maxerr maxerr = err end
        
    plot(x[e], Ue[:,e], "r", x[e], uu, "b")
end
maxerr


    

In [ ]: