diff --git a/.gitignore b/.gitignore index 079af296..838b6d86 100644 --- a/.gitignore +++ b/.gitignore @@ -9,3 +9,4 @@ /tmp/ *.vtu *.pvtu +*.pvd diff --git a/Project.toml b/Project.toml index f537399a..67bd66be 100644 --- a/Project.toml +++ b/Project.toml @@ -14,10 +14,10 @@ SparseMatricesCSR = "a0a7dd2c-ebf4-11e9-1f05-cf50bc540ca1" WriteVTK = "64499a7a-5c06-52f2-abe2-ccb03c286192" [compat] -Gridap = "0.17.8" -PartitionedArrays = "0.2.10" FillArrays = "0.8.4, 0.9, 0.10, 0.11, 0.12" +Gridap = "0.17.9" MPI = "0.16, 0.17, 0.18, 0.19" +PartitionedArrays = "0.2.10" SparseMatricesCSR = "0.6.6" WriteVTK = "1.12.0" julia = "1.3" diff --git a/src/GridapDistributed.jl b/src/GridapDistributed.jl index 9baba964..4f666597 100644 --- a/src/GridapDistributed.jl +++ b/src/GridapDistributed.jl @@ -11,6 +11,8 @@ using Gridap.CellData using Gridap.Visualization using Gridap.FESpaces using Gridap.MultiField +using Gridap.ODEs.TransientFETools +using Gridap.ODEs.ODETools using PartitionedArrays const PArrays = PartitionedArrays @@ -23,6 +25,7 @@ import Gridap.TensorValues: inner, outer, double_contraction, symmetric_part import LinearAlgebra: det, tr, cross, dot, ⋅ import Base: inv, abs, abs2, *, +, -, /, adjoint, transpose, real, imag, conj import Gridap.Fields: grad2curl +import Gridap.ODEs.ODETools: ∂t, ∂tt export FullyAssembledRows export SubAssembledRows @@ -41,4 +44,10 @@ include("DivConformingFESpaces.jl") include("MultiField.jl") +include("TransientDistributedCellField.jl") + +include("TransientMultiFieldDistributedCellField.jl") + +include("TransientFESpaces.jl") + end # module diff --git a/src/TransientDistributedCellField.jl b/src/TransientDistributedCellField.jl new file mode 100644 index 00000000..b82e6787 --- /dev/null +++ b/src/TransientDistributedCellField.jl @@ -0,0 +1,99 @@ +# Transient Distributed CellField +abstract type TransientDistributedCellField <: DistributedCellDatum end + +# Transient SingleField +struct TransientSingleFieldDistributedCellField{A} <: TransientDistributedCellField + cellfield::A + derivatives::Tuple +end + +# Constructors +function TransientFETools.TransientCellField(single_field::DistributedSingleFieldFEFunction,derivatives::Tuple) + TransientSingleFieldDistributedCellField(single_field,derivatives) +end + +function TransientFETools.TransientCellField(single_field::DistributedCellField,derivatives::Tuple) +TransientSingleFieldDistributedCellField(single_field,derivatives) +end + +# Time derivative +function ∂t(f::TransientDistributedCellField) + cellfield, derivatives = first_and_tail(f.derivatives) + TransientCellField(cellfield,derivatives) +end + +∂tt(f::TransientDistributedCellField) = ∂t(∂t(f)) + +# Integration related +function Fields.integrate(f::TransientDistributedCellField,b::DistributedMeasure) + integrate(f.cellfield,b) +end + +# Differential Operations +Fields.gradient(f::TransientDistributedCellField) = gradient(f.cellfield) +Fields.∇∇(f::TransientDistributedCellField) = ∇∇(f.cellfield) + +# Unary ops +for op in (:symmetric_part,:inv,:det,:abs,:abs2,:+,:-,:tr,:transpose,:adjoint,:grad2curl,:real,:imag,:conj) + @eval begin + ($op)(a::TransientDistributedCellField) = ($op)(a.cellfield) + end +end + +# Binary ops +for op in (:inner,:outer,:double_contraction,:+,:-,:*,:cross,:dot,:/) + @eval begin + ($op)(a::TransientDistributedCellField,b::TransientDistributedCellField) = ($op)(a.cellfield,b.cellfield) + ($op)(a::TransientDistributedCellField,b::DistributedCellField) = ($op)(a.cellfield,b) + ($op)(a::DistributedCellField,b::TransientDistributedCellField) = ($op)(a,b.cellfield) + ($op)(a::TransientDistributedCellField,b::Number) = ($op)(a.cellfield,b) + ($op)(a::Number,b::TransientDistributedCellField) = ($op)(a,b.cellfield) + ($op)(a::TransientDistributedCellField,b::Function) = ($op)(a.cellfield,b) + ($op)(a::Function,b::TransientDistributedCellField) = ($op)(a,b.cellfield) + end +end + +Base.broadcasted(f,a::TransientDistributedCellField,b::TransientDistributedCellField) = broadcasted(f,a.cellfield,b.cellfield) +Base.broadcasted(f,a::TransientDistributedCellField,b::DistributedCellField) = broadcasted(f,a.cellfield,b) +Base.broadcasted(f,a::DistributedCellField,b::TransientDistributedCellField) = broadcasted(f,a,b.cellfield) +Base.broadcasted(f,a::Number,b::TransientDistributedCellField) = broadcasted(f,a,b.cellfield) +Base.broadcasted(f,a::TransientDistributedCellField,b::Number) = broadcasted(f,a.cellfield,b) +Base.broadcasted(f,a::Function,b::TransientDistributedCellField) = broadcasted(f,a,b.cellfield) +Base.broadcasted(f,a::TransientDistributedCellField,b::Function) = broadcasted(f,a.cellfield,b) +Base.broadcasted(a::typeof(*),b::typeof(∇),f::TransientDistributedCellField) = broadcasted(a,b,f.cellfield) +Base.broadcasted(a::typeof(*),s::Fields.ShiftedNabla,f::TransientDistributedCellField) = broadcasted(a,s,f.cellfield) + +dot(::typeof(∇),f::TransientDistributedCellField) = dot(∇,f.cellfield) +outer(::typeof(∇),f::TransientDistributedCellField) = outer(∇,f.cellfield) +outer(f::TransientDistributedCellField,::typeof(∇)) = outer(f.cellfield,∇) +cross(::typeof(∇),f::TransientDistributedCellField) = cross(∇,f.cellfield) + +# Skeleton related +function Base.getproperty(f::TransientDistributedCellField, sym::Symbol) + if sym in (:⁺,:plus,:⁻, :minus) + derivatives = () + cellfield = DistributedCellField(f.cellfield,sym) + for iderivative in f.derivatives + derivatives = (derivatives...,DistributedCellField(iderivative)) + end + return TransientSingleFieldCellField(cellfield,derivatives) + else + return getfield(f, sym) + end +end + +Base.propertynames(x::TransientDistributedCellField, private::Bool=false) = propertynames(x.cellfield, private) + +for op in (:outer,:*,:dot) + @eval begin + ($op)(a::TransientDistributedCellField,b::SkeletonPair{<:DistributedCellField}) = ($op)(a.cellfield,b) + ($op)(a::SkeletonPair{<:DistributedCellField},b::TransientDistributedCellField) = ($op)(a,b.cellfield) + end +end + +Arrays.evaluate!(cache,k::Operation,a::TransientDistributedCellField,b::SkeletonPair{<:DistributedCellField}) = evaluate!(cache,k,a.cellfield,b) + +Arrays.evaluate!(cache,k::Operation,a::SkeletonPair{<:DistributedCellField},b::TransientDistributedCellField) = evaluate!(cache,k,a,b.cellfield) + +CellData.jump(a::TransientDistributedCellField) = jump(a.cellfield) +CellData.mean(a::TransientDistributedCellField) = mean(a.cellfield) diff --git a/src/TransientFESpaces.jl b/src/TransientFESpaces.jl new file mode 100644 index 00000000..6a2f7c0a --- /dev/null +++ b/src/TransientFESpaces.jl @@ -0,0 +1,66 @@ +# Functions for transient FE spaces + +Fields.evaluate(U::DistributedSingleFieldFESpace,t::Nothing) = U + +(U::DistributedSingleFieldFESpace)(t) = U + +∂t(U::DistributedSingleFieldFESpace) = HomogeneousTrialFESpace(U) +∂t(U::DistributedMultiFieldFESpace) = MultiFieldFESpace(∂t.(U.field_fe_space)) +∂tt(U::DistributedSingleFieldFESpace) = HomogeneousTrialFESpace(U) +∂tt(U::DistributedMultiFieldFESpace) = MultiFieldFESpace(∂tt.(U.field_fe_spaces)) + +function TransientMultiFieldFESpace(spaces::Vector{<:DistributedSingleFieldFESpace}) + MultiFieldFESpace(spaces) +end + +# Functions for transient FE Functions + +function ODETools.allocate_jacobian( + op::TransientFETools.TransientFEOperatorFromWeakForm, + duh::DistributedCellField, + cache) + _matdata_jacobians = TransientFETools.fill_initial_jacobians(op,duh) + matdata = _vcat_distributed_matdata(_matdata_jacobians) + allocate_matrix(op.assem_t,matdata) +end + +function ODETools.allocate_jacobian( + op::TransientFETools.TransientFEOperatorFromWeakForm, + duh::DistributedMultiFieldFEFunction, + cache) + _matdata_jacobians = TransientFETools.fill_initial_jacobians(op,duh) + matdata = _vcat_distributed_matdata(_matdata_jacobians) + allocate_matrix(op.assem_t,matdata) +end + +function ODETools.jacobians!( + A::AbstractMatrix, + op::TransientFETools.TransientFEOperatorFromWeakForm, + t::Real, + xh::TransientDistributedCellField, + γ::Tuple{Vararg{Real}}, + cache) + _matdata_jacobians = TransientFETools.fill_jacobians(op,t,xh,γ) + matdata = _vcat_distributed_matdata(_matdata_jacobians) + assemble_matrix_add!(A,op.assem_t, matdata) + A +end + +function _vcat_distributed_matdata(_matdata) + term_to_cellmat = map_parts(a->a[1],local_views(_matdata[1])) + term_to_cellidsrows = map_parts(a->a[2],local_views(_matdata[1])) + term_to_cellidscols = map_parts(a->a[3],local_views(_matdata[1])) + for j in 2:length(_matdata) + term_to_cellmat_j = map_parts(a->a[1],local_views(_matdata[j])) + term_to_cellidsrows_j = map_parts(a->a[2],local_views(_matdata[j])) + term_to_cellidscols_j = map_parts(a->a[3],local_views(_matdata[j])) + term_to_cellmat = map_parts((a,b)->vcat(a,b),local_views(term_to_cellmat),local_views(term_to_cellmat_j)) + term_to_cellidsrows = map_parts((a,b)->vcat(a,b),local_views(term_to_cellidsrows),local_views(term_to_cellidsrows_j)) + term_to_cellidscols = map_parts((a,b)->vcat(a,b),local_views(term_to_cellidscols),local_views(term_to_cellidscols_j)) + end + map_parts( (a,b,c) -> (a,b,c), + local_views(term_to_cellmat), + local_views(term_to_cellidsrows), + local_views(term_to_cellidscols) + ) +end diff --git a/src/TransientMultiFieldDistributedCellField.jl b/src/TransientMultiFieldDistributedCellField.jl new file mode 100644 index 00000000..41822882 --- /dev/null +++ b/src/TransientMultiFieldDistributedCellField.jl @@ -0,0 +1,61 @@ +# Transient MultiField +struct TransientMultiFieldDistributedCellField{A} <: TransientDistributedCellField + cellfield::A + derivatives::Tuple + transient_single_fields::Vector{<:TransientDistributedCellField} # used to iterate +end + +# Constructors +function TransientFETools.TransientCellField(multi_field::DistributedMultiFieldFEFunction,derivatives::Tuple) + transient_single_fields = _to_transient_single_distributed_fields(multi_field,derivatives) + TransientMultiFieldDistributedCellField(multi_field,derivatives,transient_single_fields) +end + +# Get single index +function Base.getindex(f::TransientMultiFieldDistributedCellField,ifield::Integer) + single_field = f.cellfield[ifield] + single_derivatives = () + for ifield_derivatives in f.derivatives + single_derivatives = (single_derivatives...,getindex(ifield_derivatives,ifield)) + end + TransientSingleFieldDistributedCellField(single_field,single_derivatives) +end + +# Get multiple indices +function Base.getindex(f::TransientMultiFieldDistributedCellField,indices::Vector{<:Int}) + cellfield = DistributedMultiFieldCellField(f.cellfield[indices],DomainStyle(f.cellfield)) + derivatives = () + for derivative in f.derivatives + derivatives = (derivatives...,DistributedMultiFieldCellField(derivative[indices],DomainStyle(derivative))) + end + transient_single_fields = _to_transient_single_distributed_fields(cellfield,derivatives) + TransientMultiFieldDistributedCellField(cellfield,derivatives,transient_single_fields) +end + +function _to_transient_single_distributed_fields(multi_field,derivatives) + transient_single_fields = TransientDistributedCellField[] + for ifield in 1:num_fields(multi_field) + single_field = multi_field[ifield] + single_derivatives = () + for ifield_derivatives in derivatives + single_derivatives = (single_derivatives...,getindex(ifield_derivatives,ifield)) + end + transient_single_field = TransientSingleFieldDistributedCellField(single_field,single_derivatives) + push!(transient_single_fields,transient_single_field) + end + transient_single_fields +end + +# Iterate functions +Base.iterate(f::TransientMultiFieldDistributedCellField) = iterate(f.transient_single_fields) +Base.iterate(f::TransientMultiFieldDistributedCellField,state) = iterate(f.transient_single_fields,state) + +# Time derivative +function ∂t(f::TransientMultiFieldDistributedCellField) + cellfield, derivatives = first_and_tail(f.derivatives) + transient_single_field_derivatives = TransientDistributedCellField[] + for transient_single_field in f.transient_single_fields + push!(transient_single_field_derivatives,∂t(transient_single_field)) + end + TransientMultiFieldDistributedCellField(cellfield,derivatives,transient_single_field_derivatives) +end diff --git a/test/HeatEquationTests.jl b/test/HeatEquationTests.jl new file mode 100644 index 00000000..940ac7e8 --- /dev/null +++ b/test/HeatEquationTests.jl @@ -0,0 +1,77 @@ +module HeatEquationTests + +using Gridap +using Gridap.ODEs +using GridapDistributed +using PartitionedArrays +using Test + +function main(parts) + θ = 0.2 + + u(x,t) = (1.0-x[1])*x[1]*(1.0-x[2])*x[2]*t + u(t::Real) = x -> u(x,t) + f(t) = x -> ∂t(u)(x,t)-Δ(u(t))(x) + + domain = (0,1,0,1) + partition = (4,4) + model = CartesianDiscreteModel(parts,domain,partition) + + order = 2 + + reffe = ReferenceFE(lagrangian,Float64,order) + V0 = FESpace( + model, + reffe, + conformity=:H1, + dirichlet_tags="boundary" + ) + U = TransientTrialFESpace(V0,u) + + Ω = Triangulation(model) + degree = 2*order + dΩ = Measure(Ω,degree) + + # + m(u,v) = ∫(u*v)dΩ + a(u,v) = ∫(∇(v)⋅∇(u))dΩ + b(t,v) = ∫(v*f(t))dΩ + + res(t,u,v) = a(u,v) + m(∂t(u),v) - b(t,v) + jac(t,u,du,v) = a(du,v) + jac_t(t,u,dut,v) = m(dut,v) + + op = TransientFEOperator(res,jac,jac_t,U,V0) + op_constant = TransientConstantMatrixFEOperator(m,a,b,U,V0) + + t0 = 0.0 + tF = 1.0 + dt = 0.1 + + U0 = U(0.0) + uh0 = interpolate_everywhere(u(0.0),U0) + + ls = LUSolver() + ode_solver = ThetaMethod(ls,dt,θ) + + sol_t = solve(ode_solver,op,uh0,t0,tF) + sol_t_const = solve(ode_solver,op_constant,uh0,t0,tF) + + l2(w) = w*w + + tol = 1.0e-6 + + for (uh_tn, tn) in sol_t + e = u(tn) - uh_tn + el2 = sqrt(sum( ∫(l2(e))dΩ )) + @test el2 < tol + end + + for (uh_tn, tn) in sol_t_const + e = u(tn) - uh_tn + el2 = sqrt(sum( ∫(l2(e))dΩ )) + @test el2 < tol + end +end + +end #module diff --git a/test/StokesOpenBoundaryTests.jl b/test/StokesOpenBoundaryTests.jl new file mode 100644 index 00000000..90af1b35 --- /dev/null +++ b/test/StokesOpenBoundaryTests.jl @@ -0,0 +1,115 @@ +module StokesOpenBoundaryTests + +using Gridap +using LinearAlgebra +using Test +using Gridap.ODEs +using GridapDistributed +using PartitionedArrays + +function main(parts) + + θ = 0.5 + + u(x,t) = VectorValue(x[1],x[2])*t + u(t::Real) = x -> u(x,t) + + p(x,t) = (x[1]-x[2])*t + p(t::Real) = x -> p(x,t) + q(x) = t -> p(x,t) + + f(t) = x -> ∂t(u)(t)(x)-Δ(u(t))(x)+ ∇(p(t))(x) + g(t) = x -> (∇⋅u(t))(x) + h(t) = x -> ∇(u(t))(x)⋅VectorValue(0.0,1.0) - p(t)(x)*VectorValue(0.0,1.0) + + parts = get_part_ids(sequential,(2,2)) + domain = (0,1,0,1) + partition = (4,4) + model = CartesianDiscreteModel(parts,domain,partition) + labels = get_face_labeling(model) + add_tag_from_tags!(labels,"neumann",6) + add_tag_from_tags!(labels,"dirichlet",[1,2,3,4,5,7,8]) + + order = 2 + + reffeᵤ = ReferenceFE(lagrangian,VectorValue{2,Float64},order) + V0 = FESpace( + model, + reffeᵤ, + conformity=:H1, + dirichlet_tags="dirichlet" + ) + + reffeₚ = ReferenceFE(lagrangian,Float64,order-1) + Q = TestFESpace( + model, + reffeₚ, + conformity=:H1 + ) + + U = TransientTrialFESpace(V0,u) + + P = TrialFESpace(Q) + + Ω = Triangulation(model) + degree = 2*order + dΩ = Measure(Ω,degree) + + Γ = Boundary(model,tags=["neumann"]) + dΓ = Measure(Γ,degree) + + # + a(u,v) = ∫(∇(u)⊙∇(v))dΩ + b((v,q),t) = ∫(v⋅f(t))dΩ + ∫(q*g(t))dΩ + ∫(v⋅h(t))dΓ + m(ut,v) = ∫(ut⋅v)dΩ + + X = TransientMultiFieldFESpace([U,P]) + Y = MultiFieldFESpace([V0,Q]) + + res(t,(u,p),(v,q)) = a(u,v) + m(∂t(u),v) - ∫((∇⋅v)*p)dΩ + ∫(q*(∇⋅u))dΩ - b((v,q),t) + jac(t,(u,p),(du,dp),(v,q)) = a(du,v) - ∫((∇⋅v)*dp)dΩ + ∫(q*(∇⋅du))dΩ + jac_t(t,(u,p),(dut,dpt),(v,q)) = m(dut,v) + + b((v,q)) = b((v,q),0.0) + + mat((du1,du2),(v1,v2)) = a(du1,v1)+a(du2,v2) + + U0 = U(0.0) + P0 = P(0.0) + X0 = X(0.0) + uh0 = interpolate_everywhere(u(0.0),U0) + ph0 = interpolate_everywhere(p(0.0),P0) + xh0 = interpolate_everywhere([uh0,ph0],X0) + + op = TransientFEOperator(res,jac,jac_t,X,Y) + + t0 = 0.0 + tF = 1.0 + dt = 0.1 + + ls = LUSolver() + ode_solver = ThetaMethod(ls,dt,θ) + + sol_t = solve(ode_solver,op,xh0,t0,tF) + + l2(w) = w⋅w + + + tol = 1.0e-6 + _t_n = t0 + + result = Base.iterate(sol_t) + + for (xh_tn, tn) in sol_t + uh_tn = xh_tn[1] + ph_tn = xh_tn[2] + writevtk(Ω,"output/tmp_stokes_OB_sol_$tn",cellfields=["u"=>uh_tn,"p"=>ph_tn]) + e = u(tn) - uh_tn + el2 = sqrt(sum( ∫(l2(e))dΩ )) + e = p(tn) - ph_tn + el2 = sqrt(sum( ∫(l2(e))dΩ )) + @test el2 < tol + end +end + +end #module diff --git a/test/TestApp/src/TestApp.jl b/test/TestApp/src/TestApp.jl index d240749b..3c5ed61e 100644 --- a/test/TestApp/src/TestApp.jl +++ b/test/TestApp/src/TestApp.jl @@ -6,4 +6,8 @@ module TestApp include("../../PLaplacianTests.jl") include("../../PoissonTests.jl") include("../../PeriodicBCsTests.jl") -end + include("../../TransientDistributedCellFieldTests.jl") + include("../../TransientMultiFieldDistributedCellFieldTests.jl") + include("../../HeatEquationTests.jl") + include("../../StokesOpenBoundaryTests.jl") +end diff --git a/test/TransientDistributedCellFieldTests.jl b/test/TransientDistributedCellFieldTests.jl new file mode 100644 index 00000000..314bee04 --- /dev/null +++ b/test/TransientDistributedCellFieldTests.jl @@ -0,0 +1,49 @@ +module TransientDistributedCellFieldTests + +using Gridap +using GridapDistributed +using Gridap.ODEs.ODETools: ∂t, ∂tt +using Gridap.ODEs.TransientFETools: TransientCellField +using PartitionedArrays +using Test + +function main(parts) + domain = (0,1,0,1) + cells = (4,4) + 𝒯 = CartesianDiscreteModel(parts,domain,cells) + Ω = Interior(𝒯) + dΩ = Measure(Ω,2) + + f(t) = t^2 + df(t) = 2t + ddf(t) = 2 + + a(t) = CellField(f(t),Ω) + da(t) = CellField(df(t),Ω) + dda(t) = CellField(ddf(t),Ω) + @test isa(a(0),GridapDistributed.DistributedCellField) + @test isa(da(0),GridapDistributed.DistributedCellField) + @test isa(dda(0),GridapDistributed.DistributedCellField) + + b(t) = TransientCellField(a(t),(da(t),dda(t))) + @test isa(b(0),GridapDistributed.TransientDistributedCellField) + @test isa(b(0),GridapDistributed.TransientSingleFieldDistributedCellField) + + db(t) = ∂t(b(t)) + @test isa(db(0),GridapDistributed.TransientDistributedCellField) + @test isa(db(0),GridapDistributed.TransientSingleFieldDistributedCellField) + + ddb(t) = ∂t(db(t)) + @test isa(ddb(0),GridapDistributed.TransientDistributedCellField) + @test isa(ddb(0),GridapDistributed.TransientSingleFieldDistributedCellField) + + @test (∑(∫(a(0.5))dΩ)) ≈ 0.25 + @test (∑(∫(da(0.5))dΩ)) ≈ 1.0 + @test (∑(∫(dda(0.5))dΩ)) ≈ 2.0 + @test (∑(∫(b(0.5))dΩ)) ≈ 0.25 + @test (∑(∫(db(0.5))dΩ)) ≈ 1.0 + @test (∑(∫(ddb(0.5))dΩ)) ≈ 2.0 + @test (∑(∫(∂tt(b(0.5)))dΩ)) ≈ 2.0 +end + +end diff --git a/test/TransientMultiFieldDistributedCellFieldTests.jl b/test/TransientMultiFieldDistributedCellFieldTests.jl new file mode 100644 index 00000000..a1a73c1f --- /dev/null +++ b/test/TransientMultiFieldDistributedCellFieldTests.jl @@ -0,0 +1,66 @@ +module TransientMultiFieldDistributedCellFieldTests + +using Gridap +using GridapDistributed +using Gridap.ODEs.ODETools: ∂t, ∂tt +using Gridap.ODEs.TransientFETools: TransientCellField +using Gridap.ODEs.TransientFETools: TransientTrialFESpace, TransientMultiFieldFESpace +using PartitionedArrays +using Test + +function main(parts) + domain = (0,1,0,1) + cells = (4,4) + 𝒯 = CartesianDiscreteModel(parts,domain,cells) + Ω = Interior(𝒯) + dΩ = Measure(Ω,2) + + reffe = ReferenceFE(lagrangian,Float64,1) + V = FESpace(𝒯, reffe) + U = TrialFESpace(V) + Ut = TransientTrialFESpace(V) + Y = MultiFieldFESpace([V,V]) + X = MultiFieldFESpace([U,U]) + Xt = TransientMultiFieldFESpace([Ut,Ut]) + + f(t) = t^2 + df(t) = 2t + ddf(t) = 2 + + a(t) = interpolate([f(t),f(t)],X) + da(t) = interpolate([df(t),df(t)],X) + dda(t) = interpolate([ddf(t),ddf(t)],X) + @test isa(a(0),GridapDistributed.DistributedMultiFieldFEFunction) + @test isa(da(0),GridapDistributed.DistributedMultiFieldFEFunction) + @test isa(dda(0),GridapDistributed.DistributedMultiFieldFEFunction) + + b(t) = TransientCellField(a(t),(da(t),dda(t))) + @test isa(b(0),GridapDistributed.TransientDistributedCellField) + @test isa(b(0),GridapDistributed.TransientMultiFieldDistributedCellField) + + db(t) = ∂t(b(t)) + @test isa(db(0),GridapDistributed.TransientDistributedCellField) + @test isa(db(0),GridapDistributed.TransientMultiFieldDistributedCellField) + + ddb(t) = ∂t(db(t)) + @test isa(ddb(0),GridapDistributed.TransientDistributedCellField) + @test isa(ddb(0),GridapDistributed.TransientMultiFieldDistributedCellField) + + b1(t) = b(t)[1] + @test isa(b1(0),GridapDistributed.TransientDistributedCellField) + @test isa(b1(0),GridapDistributed.TransientSingleFieldDistributedCellField) + + db1(t) = ∂t(b1(t)) + @test isa(db1(0),GridapDistributed.TransientDistributedCellField) + @test isa(db1(0),GridapDistributed.TransientSingleFieldDistributedCellField) + + ddb1(t) = ∂t(db1(t)) + @test isa(ddb1(0),GridapDistributed.TransientDistributedCellField) + @test isa(ddb1(0),GridapDistributed.TransientSingleFieldDistributedCellField) + + @test (∑(∫(b(0.5)[1])dΩ)) == (∑(∫(b1(0.5))dΩ)) + @test (∑(∫(db(0.5)[1])dΩ)) == (∑(∫(db1(0.5))dΩ)) + @test (∑(∫(ddb(0.5)[1])dΩ)) == (∑(∫(ddb1(0.5))dΩ)) +end + +end diff --git a/test/mpi/runtests_np4_body.jl b/test/mpi/runtests_np4_body.jl index 34e9c027..cdb99e78 100644 --- a/test/mpi/runtests_np4_body.jl +++ b/test/mpi/runtests_np4_body.jl @@ -1,29 +1,41 @@ function all_tests(parts) - display(parts) + display(parts) - t = PArrays.PTimer(parts,verbose=true) - PArrays.tic!(t) + t = PArrays.PTimer(parts,verbose=true) + PArrays.tic!(t) - TestApp.GeometryTests.main(parts) - PArrays.toc!(t,"Geometry") + TestApp.GeometryTests.main(parts) + PArrays.toc!(t,"Geometry") - TestApp.CellDataTests.main(parts) - PArrays.toc!(t,"CellData") + TestApp.CellDataTests.main(parts) + PArrays.toc!(t,"CellData") - TestApp.FESpacesTests.main(parts) - PArrays.toc!(t,"FESpaces") + TestApp.FESpacesTests.main(parts) + PArrays.toc!(t,"FESpaces") - TestApp.MultiFieldTests.main(parts) - PArrays.toc!(t,"MultiField") + TestApp.MultiFieldTests.main(parts) + PArrays.toc!(t,"MultiField") - TestApp.PoissonTests.main(parts) - PArrays.toc!(t,"Poisson") + TestApp.PoissonTests.main(parts) + PArrays.toc!(t,"Poisson") TestApp.PLaplacianTests.main(parts) PArrays.toc!(t,"PLaplacian") + TestApp.TransientDistributedCellFieldTests.main(parts) + PArrays.toc!(t,"TransientDistributedCellField") + + TestApp.TransientMultiFieldDistributedCellFieldTests.main(parts) + PArrays.toc!(t,"TransientMultiFieldDistributedCellField") + TestApp.PeriodicBCsTests.main(parts) PArrays.toc!(t,"PeriodicBCs") + TestApp.HeatEquationTests.main(parts) + PArrays.toc!(t,"HeatEquation") + + TestApp.StokesOpenBoundaryTests.main(parts) + PArrays.toc!(t,"StokesOpenBoundary") + display(t) end diff --git a/test/sequential/HeatEquationTests.jl b/test/sequential/HeatEquationTests.jl new file mode 100644 index 00000000..3af75560 --- /dev/null +++ b/test/sequential/HeatEquationTests.jl @@ -0,0 +1,5 @@ +module HeatEquationTestsSeq +using PartitionedArrays +include("../HeatEquationTests.jl") +prun(HeatEquationTests.main,sequential,(2,2)) +end # module diff --git a/test/sequential/StokesOpenBoundaryTests.jl b/test/sequential/StokesOpenBoundaryTests.jl new file mode 100644 index 00000000..b2fe654c --- /dev/null +++ b/test/sequential/StokesOpenBoundaryTests.jl @@ -0,0 +1,5 @@ +module StokesOpenBoundaryTestsSeq +using PartitionedArrays +include("../StokesOpenBoundaryTests.jl") +prun(StokesOpenBoundaryTests.main,sequential,(2,2)) +end # module diff --git a/test/sequential/TransientDistributedCellFieldTests.jl b/test/sequential/TransientDistributedCellFieldTests.jl new file mode 100644 index 00000000..319a32d2 --- /dev/null +++ b/test/sequential/TransientDistributedCellFieldTests.jl @@ -0,0 +1,5 @@ +module TransientDistributedCellFieldTestsSeq +using PartitionedArrays +include("../TransientDistributedCellFieldTests.jl") +prun(TransientDistributedCellFieldTests.main,sequential,(2,2)) +end # module diff --git a/test/sequential/TransientMultiFieldDistributedCellFieldTests.jl b/test/sequential/TransientMultiFieldDistributedCellFieldTests.jl new file mode 100644 index 00000000..6e8817de --- /dev/null +++ b/test/sequential/TransientMultiFieldDistributedCellFieldTests.jl @@ -0,0 +1,5 @@ +module TransientMultiFieldDistributedCellFieldTestsSeq +using PartitionedArrays +include("../TransientMultiFieldDistributedCellFieldTests.jl") +prun(TransientMultiFieldDistributedCellFieldTests.main,sequential,(2,2)) +end # module diff --git a/test/sequential/runtests.jl b/test/sequential/runtests.jl index bc1c184e..504baf3e 100644 --- a/test/sequential/runtests.jl +++ b/test/sequential/runtests.jl @@ -18,4 +18,16 @@ using Test @time @testset "DivConformingTests" begin include("DivConformingTests.jl") end +@time @testset "TransientDistributedCellFieldTests" begin + include("TransientDistributedCellFieldTests.jl") +end + +@time @testset "TransientMultiFieldDistributedCellFieldTests" begin + include("TransientMultiFieldDistributedCellFieldTests.jl") +end + +@time @testset "HeatEquation" begin include("HeatEquationTests.jl") end + +@time @testset "StokesOpenBoundary" begin include("StokesOpenBoundaryTests.jl") end + end # module