diff --git a/src/FESpaces/FESpaceFactories.jl b/src/FESpaces/FESpaceFactories.jl index 28c810aa9..ae6a1b276 100644 --- a/src/FESpaces/FESpaceFactories.jl +++ b/src/FESpaces/FESpaceFactories.jl @@ -4,13 +4,15 @@ function FESpace(;kwargs...) constraint = _get_kwarg(:constraint,kwargs,nothing) - # constraint = nothing reffe = _get_kwarg(:reffe,kwargs) - @notimplementedif !isa(reffe,Symbol) "For the moment, reffe can only be a symbol" fespace = nothing - if reffe in [:Lagrangian,:QLagrangian,:PLagrangian,:SLagrangian] + if isa(reffe,ReferenceFE) + + fespace = _setup_generic_space(kwargs) + + elseif reffe in [:Lagrangian,:QLagrangian,:PLagrangian,:SLagrangian] fespace = _setup_lagrange_spaces(kwargs) @@ -96,6 +98,40 @@ function TestFESpace(;kwargs...) FESpace(;kwargs...) end +function _setup_generic_space(kwargs) + + reffe = _get_kwarg(:reffe,kwargs) + model = _get_kwarg(:model,kwargs) + labels = _get_kwarg(:labels,kwargs,get_face_labeling(model)) + conformity = _get_kwarg(:conformity,kwargs,true) + diritags = _get_kwarg(:dirichlet_tags,kwargs,Int[]) + order = _get_kwarg(:order,kwargs,nothing) + dofspace = _get_kwarg(:dof_space,kwargs,:reference) + ( dofspace == :reference ? true : false ) + + is_ref = (dofspace==:reference) + + Tf = _get_kwarg(:valuetype,kwargs,VectorValue{1,Float64}) + T = eltype(Tf) + + + polytopes = get_polytopes(model) + if length(polytopes) > 2 || polytopes[1] != get_polytope(reffe) + @unreachable "RefFE and geometrical model are not compatible in generic FESpace constructor" + end + + reffes = [ reffe ] + + if conformity in [true, :default] + V = ConformingFESpace(reffes,model,labels,diritags,nothing,is_ref) + else + s = "Conformity is determined by $reffe reference FE and cannot be imposed in generic FESpace constructor, leave it `true` or simply do not specify it" + @unreachable s + end + + V +end + function _setup_hdiv_space(kwargs) reffe = _get_kwarg(:reffe,kwargs) @@ -213,7 +249,7 @@ function _setup_lagrange_spaces(kwargs) elseif reffe == :QLagrangian _reffes = [QDiscRefFE(T,p,order) for p in polytopes] else - @unreachable "Not possible to use a $reffe reffe on polytopoes $(polytopes...)" + @unreachable "Not possible to use a $reffe reffe on polytopes $(polytopes...)" end end diff --git a/test/FESpacesTests/CurlConformingFESpacesTests.jl b/test/FESpacesTests/CurlConformingFESpacesTests.jl index 9a976f248..2b81540de 100644 --- a/test/FESpacesTests/CurlConformingFESpacesTests.jl +++ b/test/FESpacesTests/CurlConformingFESpacesTests.jl @@ -4,6 +4,7 @@ using Test using Gridap using LinearAlgebra using Gridap.FESpaces +using Gridap.ReferenceFEs domain =(0,1,0,1,0,1) partition = (3,3,3) @@ -28,8 +29,6 @@ test_single_field_fe_space(V) U = TrialFESpace(V,u) uh = interpolate(U,u) -length(uh.free_values) -uh.free_values e = u - uh @@ -39,6 +38,18 @@ quad = CellQuadrature(trian,order) el2 = sqrt(sum(integrate(inner(e,e),trian,quad))) @test el2 < 1.0e-10 +T = Float64 +reffe = NedelecRefFE(T,HEX,order) +V = FESpace(model=model,reffe=reffe,dirichlet_tags = [21,22]) +test_single_field_fe_space(V) + +U = TrialFESpace(V,u) + +uh = interpolate(U,u) + +el2 = sqrt(sum(integrate(inner(e,e),trian,quad))) +@test el2 < 1.0e-10 + # using Gridap.Visualization # writevtk(trian,"trian",nsubcells=10,cellfields=["uh"=>uh]) diff --git a/test/FESpacesTests/DivConformingFESpacesTests.jl b/test/FESpacesTests/DivConformingFESpacesTests.jl index 4c8e5e437..1c4134f15 100644 --- a/test/FESpacesTests/DivConformingFESpacesTests.jl +++ b/test/FESpacesTests/DivConformingFESpacesTests.jl @@ -35,6 +35,18 @@ quad = CellQuadrature(trian,order) el2 = sqrt(sum(integrate(inner(e,e),trian,quad))) @test el2 < 1.0e-10 +T = Float64 +reffe = RaviartThomasRefFE(T,QUAD,order) +V = FESpace(model=model,reffe=reffe,dirichlet_tags = [1,6]) +test_single_field_fe_space(V) + +U = TrialFESpace(V,u) + +uh = interpolate(U,u) + +el2 = sqrt(sum(integrate(inner(e,e),trian,quad))) +@test el2 < 1.0e-10 + #using Gridap.Visualization # #writevtk(trian,"trian",nsubcells=10,cellfields=["uh"=>uh]) diff --git a/test/FESpacesTests/FESpaceFactoriesTests.jl b/test/FESpacesTests/FESpaceFactoriesTests.jl index 4167e406f..967b52074 100644 --- a/test/FESpacesTests/FESpaceFactoriesTests.jl +++ b/test/FESpacesTests/FESpaceFactoriesTests.jl @@ -13,6 +13,11 @@ partition = (3,3) model = CartesianDiscreteModel(domain,partition) order = 3 +T = VectorValue{2,Float64} +reffe = QDiscRefFE(T,QUAD,2) +V = FESpace(model=model,reffe=reffe) +@test isa(V,UnconstrainedFESpace) + V = FESpace( model=model, valuetype=Float64,