From 14b54271646c78876397ce634bc72148096d0997 Mon Sep 17 00:00:00 2001 From: JordiManyer Date: Mon, 11 Dec 2023 12:38:12 +1100 Subject: [PATCH] Started porting code to new version of Gridap.ODEs --- src/GridapDistributed.jl | 4 +--- src/TransientDistributedCellField.jl | 10 +++++----- src/TransientFESpaces.jl | 20 +++++++++---------- ...TransientMultiFieldDistributedCellField.jl | 4 ++-- test/StokesOpenBoundaryTests.jl | 8 ++++++-- test/TransientDistributedCellFieldTests.jl | 4 ++-- ...ientMultiFieldDistributedCellFieldTests.jl | 6 +++--- 7 files changed, 29 insertions(+), 27 deletions(-) diff --git a/src/GridapDistributed.jl b/src/GridapDistributed.jl index 27f4fc02..93e1077c 100644 --- a/src/GridapDistributed.jl +++ b/src/GridapDistributed.jl @@ -12,8 +12,7 @@ using Gridap.CellData using Gridap.Visualization using Gridap.FESpaces using Gridap.MultiField -using Gridap.ODEs.TransientFETools -using Gridap.ODEs.ODETools +using Gridap.ODEs using MPI using PartitionedArrays @@ -29,7 +28,6 @@ import Gridap.TensorValues: inner, outer, double_contraction, symmetric_part import LinearAlgebra: det, tr, cross, dot, ⋅, diag import Base: inv, abs, abs2, *, +, -, /, adjoint, transpose, real, imag, conj, getproperty, propertynames import Gridap.Fields: grad2curl -import Gridap.ODEs.ODETools: ∂t, ∂tt export FullyAssembledRows export SubAssembledRows diff --git a/src/TransientDistributedCellField.jl b/src/TransientDistributedCellField.jl index 8934f40f..3186d986 100644 --- a/src/TransientDistributedCellField.jl +++ b/src/TransientDistributedCellField.jl @@ -10,21 +10,21 @@ end local_views(f::TransientSingleFieldDistributedCellField) = local_views(f.cellfield) # Constructors -function TransientFETools.TransientCellField(single_field::DistributedSingleFieldFEFunction,derivatives::Tuple) +function ODEs.TransientCellField(single_field::DistributedSingleFieldFEFunction,derivatives::Tuple) TransientSingleFieldDistributedCellField(single_field,derivatives) end -function TransientFETools.TransientCellField(single_field::DistributedCellField,derivatives::Tuple) -TransientSingleFieldDistributedCellField(single_field,derivatives) +function ODEs.TransientCellField(single_field::DistributedCellField,derivatives::Tuple) + TransientSingleFieldDistributedCellField(single_field,derivatives) end # Time derivative -function ∂t(f::TransientDistributedCellField) +function ODEs.∂t(f::TransientDistributedCellField) cellfield, derivatives = first_and_tail(f.derivatives) TransientCellField(cellfield,derivatives) end -∂tt(f::TransientDistributedCellField) = ∂t(∂t(f)) +ODEs.∂tt(f::TransientDistributedCellField) = ∂t(∂t(f)) # Integration related function Fields.integrate(f::TransientDistributedCellField,b::DistributedMeasure) diff --git a/src/TransientFESpaces.jl b/src/TransientFESpaces.jl index 938d86c8..6b57a014 100644 --- a/src/TransientFESpaces.jl +++ b/src/TransientFESpaces.jl @@ -4,10 +4,10 @@ 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)) +ODEs.∂t(U::DistributedSingleFieldFESpace) = HomogeneousTrialFESpace(U) +ODEs.∂t(U::DistributedMultiFieldFESpace) = MultiFieldFESpace(∂t.(U.field_fe_space)) +ODEs.∂tt(U::DistributedSingleFieldFESpace) = HomogeneousTrialFESpace(U) +ODEs.∂tt(U::DistributedMultiFieldFESpace) = MultiFieldFESpace(∂tt.(U.field_fe_spaces)) function TransientMultiFieldFESpace(spaces::Vector{<:DistributedSingleFieldFESpace}) MultiFieldFESpace(spaces) @@ -15,24 +15,24 @@ end # Functions for transient FE Functions -function ODETools.allocate_jacobian( - op::TransientFETools.TransientFEOperatorFromWeakForm, +function ODEs.allocate_jacobian( + op::ODEs.TransientFEOpFromWeakForm, t0::Real, duh::Union{DistributedCellField,DistributedMultiFieldFEFunction}, cache) - _matdata_jacobians = TransientFETools.fill_initial_jacobians(op,t0,duh) + _matdata_jacobians = ODEs.fill_initial_jacobians(op,t0,duh) matdata = _vcat_distributed_matdata(_matdata_jacobians) allocate_matrix(op.assem_t,matdata) end -function ODETools.jacobians!( +function ODEs.jacobians!( A::AbstractMatrix, - op::TransientFETools.TransientFEOperatorFromWeakForm, + op::ODEs.TransientFEOpFromWeakForm, t::Real, xh::TransientDistributedCellField, γ::Tuple{Vararg{Real}}, cache) - _matdata_jacobians = TransientFETools.fill_jacobians(op,t,xh,γ) + _matdata_jacobians = ODEs.fill_jacobians(op,t,xh,γ) matdata = _vcat_distributed_matdata(_matdata_jacobians) assemble_matrix_add!(A,op.assem_t, matdata) A diff --git a/src/TransientMultiFieldDistributedCellField.jl b/src/TransientMultiFieldDistributedCellField.jl index 4bd280db..3f4fd47f 100644 --- a/src/TransientMultiFieldDistributedCellField.jl +++ b/src/TransientMultiFieldDistributedCellField.jl @@ -8,7 +8,7 @@ end local_views(f::TransientMultiFieldDistributedCellField) = local_views(f.cellfield) # Constructors -function TransientFETools.TransientCellField(multi_field::DistributedMultiFieldFEFunction,derivatives::Tuple) +function ODEs.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 @@ -53,7 +53,7 @@ Base.iterate(f::TransientMultiFieldDistributedCellField) = iterate(f.transient_ Base.iterate(f::TransientMultiFieldDistributedCellField,state) = iterate(f.transient_single_fields,state) # Time derivative -function ∂t(f::TransientMultiFieldDistributedCellField) +function ODEs.∂t(f::TransientMultiFieldDistributedCellField) cellfield, derivatives = first_and_tail(f.derivatives) transient_single_field_derivatives = TransientDistributedCellField[] for transient_single_field in f.transient_single_fields diff --git a/test/StokesOpenBoundaryTests.jl b/test/StokesOpenBoundaryTests.jl index 8431e89b..93bf733a 100644 --- a/test/StokesOpenBoundaryTests.jl +++ b/test/StokesOpenBoundaryTests.jl @@ -81,7 +81,7 @@ function main(distribute,parts) ph0 = interpolate_everywhere(p(0.0),P0) xh0 = interpolate_everywhere([uh0,ph0],X0) - op = TransientFEOperator(res,jac,jac_t,X,Y) + op = TransientFEOperator(res,(jac,jac_t),X,Y) t0 = 0.0 tF = 1.0 @@ -90,7 +90,7 @@ function main(distribute,parts) ls = LUSolver() ode_solver = ThetaMethod(ls,dt,θ) - sol_t = solve(ode_solver,op,xh0,t0,tF) + sol_t = solve(ode_solver,op,t0,tF,xh0) l2(w) = w⋅w @@ -112,4 +112,8 @@ function main(distribute,parts) end end +with_debug() do distribute + main(distribute,(2,2)) +end + end #module diff --git a/test/TransientDistributedCellFieldTests.jl b/test/TransientDistributedCellFieldTests.jl index e356abe4..aed283d1 100644 --- a/test/TransientDistributedCellFieldTests.jl +++ b/test/TransientDistributedCellFieldTests.jl @@ -2,8 +2,8 @@ module TransientDistributedCellFieldTests using Gridap using GridapDistributed -using Gridap.ODEs.ODETools: ∂t, ∂tt -using Gridap.ODEs.TransientFETools: TransientCellField +using Gridap.ODEs: ∂t, ∂tt +using Gridap.ODEs: TransientCellField using PartitionedArrays using Test diff --git a/test/TransientMultiFieldDistributedCellFieldTests.jl b/test/TransientMultiFieldDistributedCellFieldTests.jl index bdb6eca6..76b66a76 100644 --- a/test/TransientMultiFieldDistributedCellFieldTests.jl +++ b/test/TransientMultiFieldDistributedCellFieldTests.jl @@ -2,9 +2,9 @@ module TransientMultiFieldDistributedCellFieldTests using Gridap using GridapDistributed -using Gridap.ODEs.ODETools: ∂t, ∂tt -using Gridap.ODEs.TransientFETools: TransientCellField -using Gridap.ODEs.TransientFETools: TransientTrialFESpace, TransientMultiFieldFESpace +using Gridap.ODEs: ∂t, ∂tt +using Gridap.ODEs: TransientCellField +using Gridap.ODEs: TransientTrialFESpace, TransientMultiFieldFESpace using PartitionedArrays using Test