Skip to content

Commit

Permalink
Now collect used for closure dofs
Browse files Browse the repository at this point in the history
  • Loading branch information
santiagobadia committed Sep 5, 2019
1 parent 529c764 commit 6b1c62f
Show file tree
Hide file tree
Showing 2 changed files with 54 additions and 49 deletions.
2 changes: 1 addition & 1 deletion src/RefFEs/RefFEs.jl
Original file line number Diff line number Diff line change
Expand Up @@ -110,7 +110,7 @@ nfacedofs(this::LagrangianRefFE{D,T} where {D,T}) = this.nfacedofs
function closurenfacedofs(this::LagrangianRefFE{D,T} where {D,T})
cv1 = CellValueFromArray(this.polytope.nf_nfs) # cell to index
cv2 = CellValueFromArray(this.nfacedofs)
return CellVectorByComposition(cv1,cv2)
return collect(CellVectorByComposition(cv1,cv2))
end

# Generate the linear nodes by computing the polytope vertices. Only the
Expand Down
101 changes: 53 additions & 48 deletions test/temp_devs/HighOrderTets.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,71 +3,76 @@ module HighOrderTets
##
using Gridap, Test
using Gridap.CellValuesGallery
# 1) Clean constructors without D or T

# Create dofbasis using node array for Lagrangian FEs
ufun(x) = x[1] + x[2]
ufun_grad(x) = VectorValue(1.0,1.0)
(::typeof(ufun)) = ufun_grad
bfun(x) = 0.0

# Create BasisWithChangeOfBasis
# i.e., CanonicalBasis given DOFs
# Construct the discrete model
model = CartesianDiscreteModel(domain=(0.0,1.0,0.0,1.0), partition=(2,2))
model = simplexify(model)

# nfacetoowndofs
order = 2
diritag = "boundary"
fespace = ConformingFESpace(Float64,model,order,diritag)

# D = 1
#
# Define test and trial
V = TestFESpace(fespace)
U = TrialFESpace(fespace,ufun)

# Define integration mesh and quadrature
trian = Triangulation(model)
quad = CellQuadrature(trian,order=6)

# Closure n-face nodes
# Method that given the set of nodes and the nfacedofs, returns the
# nface dofs on the closure of the n-face (only sense for Lagrangian FEs)
##
cell_to_x_l = [2,3,1,3,4,4,3,2,5,4,3,4]
cell_to_x_p = [1,4,4,7,13]
cell_to_x = CellVectorFromDataAndPtrs(cell_to_x_l,cell_to_x_p)
x_to_vals_l = [5,4,1,2,3,6,7,8,9,10]
x_to_vals_p = [1,3,5,6,10,11]
x_to_vals = CellVectorFromDataAndPtrs(x_to_vals_l,x_to_vals_p)
cell_to_vals = CellVectorByComposition(cell_to_x, x_to_vals)
a = [
[1, 2, 3, 5, 4],
Int[],
[3, 6, 7, 8, 9, 6, 7, 8, 9],
[3, 1, 2, 10, 6, 7, 8, 9, 3, 6, 7, 8, 9]]
##
cell_to_vals
D = 2
orders=[3,3]
p = Gridap.Polytope(1,1)
reffe = LagrangianRefFE(Float64,p,orders)
# Define forms
a(v,u) = inner((v), (u))
b(v) = inner(v,bfun)

reffe.nfacedofs
p.nf_nfs
reffe.nfacedofs[p.nf_nfs[i_nf]]
# Define Assembler
assem = SparseMatrixAssembler(V,U)

reffe.polytope
cv1 = CellValueFromArray(reffe.polytope.nf_nfs) # cell to index
cv2 = CellValueFromArray(reffe.nfacedofs)
a = CellVectorByComposition(cv1,cv2)
# Define the FEOperator
op = LinearFEOperator(a,b,V,U,assem,trian,quad)

# Define the FESolver
ls = LUSolver()
solver = LinearFESolver(ls)

# Solve!
uh = solve(solver,op)

# Define exact solution and error
u = CellField(trian,ufun)
e = u - uh

i_nf = 9
nf_nfs = reffe.nfacedofs[p.nf_nfs[i_nf]]
# Define norms to measure the error
l2(u) = inner(u,u)
h1(u) = a(u,u) + l2(u)

# Compute errors
el2 = sqrt(sum( integrate(l2(e),trian,quad) ))
eh1 = sqrt(sum( integrate(h1(e),trian,quad) ))

length(orders) == 2
#
p = Polytope(1,2)
@test el2 < 1.e-8
@test eh1 < 1.e-8

orders=[2,3]
# 1) Clean constructors without D or T

nodes, nfacenodes = Gridap.RefFEs._high_order_lagrangian_nodes_polytope(p,orders)
# Create dofbasis using node array for Lagrangian FEs

orders = [1,1,1]
# Create BasisWithChangeOfBasis
# i.e., CanonicalBasis given DOFs

if !( all(orders.==1))
b=10
end
b
# nfacetoowndofs

# D = 1
#


# Closure n-face nodes
# Method that given the set of nodes and the nfacedofs, returns the
# nface dofs on the closure of the n-face (only sense for Lagrangian FEs)
##

end # module

0 comments on commit 6b1c62f

Please sign in to comment.