From 2caea13914c721ba87ce77bc03fb9e1b0ea0a7e9 Mon Sep 17 00:00:00 2001 From: Oriol Colomes Date: Tue, 5 Oct 2021 14:43:26 +0200 Subject: [PATCH 1/6] Fixing TransientFEOperator for autoDiff --- src/TransientFETools/TransientFEOperators.jl | 13 ++++++++++--- test/TransientFEsTests/HeatEquationAutoDiffTests.jl | 8 ++++---- 2 files changed, 14 insertions(+), 7 deletions(-) diff --git a/src/TransientFETools/TransientFEOperators.jl b/src/TransientFETools/TransientFEOperators.jl index 8a1a629..02b4f7f 100644 --- a/src/TransientFETools/TransientFEOperators.jl +++ b/src/TransientFETools/TransientFEOperators.jl @@ -193,11 +193,18 @@ function TransientFEOperator(res::Function,jac::Function,jac_t::Function, end function TransientFEOperator(res::Function,trial,test;order::Integer=1) - jacs = () - for i in 1:order+1 + function jac_0(t,x,dx0,dv) + function res_0(y) + x0 = (y,x[2:end]...) + res(t,x0,dv) + end + jacobian(res_0,x[1]) + end + jacs = (jac_0,) + for i in 2:order+1 function jac_i(t,x,dxi,dv) function res_i(y) - x[i] = y + x = (x[1:i-1]...,y,x[i+1:end]...) res(t,x,dv) end jacobian(res_i,x[i]) diff --git a/test/TransientFEsTests/HeatEquationAutoDiffTests.jl b/test/TransientFEsTests/HeatEquationAutoDiffTests.jl index 3ec19eb..e548305 100644 --- a/test/TransientFEsTests/HeatEquationAutoDiffTests.jl +++ b/test/TransientFEsTests/HeatEquationAutoDiffTests.jl @@ -50,11 +50,11 @@ dv = get_fe_basis(V0) du = get_trial_fe_basis(U₀) uh = FEFunction(U₀,rand(num_free_dofs(U₀))) -cell_j = get_array(jac(0.5,uh,uh,du,dv)) -cell_j_t = get_array(jac_t(0.5,uh,uh,du,dv)) +cell_j = get_array(jac(0.5,(uh,uh),du,dv)) +cell_j_t = get_array(jac_t(0.5,(uh,uh),du,dv)) -cell_j_auto = get_array(jacobian(x->res(0.5,x,uh,dv),uh)) -cell_j_t_auto = get_array(jacobian(x->res(0.5,uh,x,dv),uh)) +cell_j_auto = get_array(jacobian(x->res(0.5,(x,uh),dv),uh)) +cell_j_t_auto = get_array(jacobian(x->res(0.5,(uh,x),dv),uh)) test_array(cell_j_auto,cell_j,≈) test_array(cell_j_t_auto,cell_j_t,≈) From c866ee8048a65b87b36de2ffae3bc852c202800c Mon Sep 17 00:00:00 2001 From: Oriol Colomes Date: Tue, 5 Oct 2021 17:31:12 +0200 Subject: [PATCH 2/6] Updated StokesAutodiff --- .../TransientFEsTests/StokesEquationAutoDiffTests.jl | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/test/TransientFEsTests/StokesEquationAutoDiffTests.jl b/test/TransientFEsTests/StokesEquationAutoDiffTests.jl index b4d9fa0..a1d0a25 100644 --- a/test/TransientFEsTests/StokesEquationAutoDiffTests.jl +++ b/test/TransientFEsTests/StokesEquationAutoDiffTests.jl @@ -1,4 +1,4 @@ -module StokesEquationTests +module StokesEquationAutoDiffTests using Gridap using ForwardDiff @@ -6,7 +6,7 @@ using LinearAlgebra using Test using GridapODEs.ODETools using GridapODEs.TransientFETools -using Gridap.FESpaces: get_algebraic_operator +using Gridap.FESpaces using Gridap.Arrays: test_array # using GridapODEs.ODETools: ThetaMethodLinear @@ -76,11 +76,11 @@ dy = get_fe_basis(Y) dx = get_trial_fe_basis(X₀) xh = FEFunction(X₀,rand(num_free_dofs(X₀))) -cell_j = get_array(jac(0.5,xh,xh,dx,dy)) -cell_j_t = get_array(jac_t(0.5,xh,xh,dx,dy)) +cell_j = get_array(jac(0.5,(xh,xh),dx,dy)) +cell_j_t = get_array(jac_t(0.5,(xh,xh),dx,dy)) -cell_j_auto = get_array(jacobian(x->res(0.5,x,xh,dy),xh)) -cell_j_t_auto = get_array(jacobian(x->res(0.5,xh,x,dy),xh)) +cell_j_auto = get_array(jacobian(x->res(0.5,(x,xh),dy),xh)) +cell_j_t_auto = get_array(jacobian(x->res(0.5,(xh,x),dy),xh)) test_array(cell_j_auto,cell_j,≈) test_array(cell_j_t_auto,cell_j_t,≈) From 0ceabc11c7967486db45f3f689710be1390ba7d8 Mon Sep 17 00:00:00 2001 From: Oriol Colomes Date: Tue, 5 Oct 2021 17:32:35 +0200 Subject: [PATCH 3/6] Added HeatEquationAutoDiffTests.jl to runtests.jl --- test/TransientFEsTests/runtests.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/test/TransientFEsTests/runtests.jl b/test/TransientFEsTests/runtests.jl index 95a0ca2..05797ad 100644 --- a/test/TransientFEsTests/runtests.jl +++ b/test/TransientFEsTests/runtests.jl @@ -24,4 +24,6 @@ using Test @testset "DGHeatEquationTests" begin include("DGHeatEquationTests.jl") end +@testset "HeatEquationAutoDiffTests" begin include("HeatEquationAutoDiffTests.jl") end + end # module From 218203d1b7fba57bd4f04ea181ae55b366651ca1 Mon Sep 17 00:00:00 2001 From: Oriol Colomes Date: Mon, 25 Oct 2021 17:12:34 +0200 Subject: [PATCH 4/6] debugging Forward Euler --- .../ForwardEulerHeatEquationTests.jl | 22 +++++++++++-------- 1 file changed, 13 insertions(+), 9 deletions(-) diff --git a/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl b/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl index 63427bc..af78fe1 100644 --- a/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl +++ b/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl @@ -1,4 +1,4 @@ -module HeatEquationTests +module ForwardEulerHeatEquationTests using Gridap using ForwardDiff @@ -30,25 +30,29 @@ model = CartesianDiscreteModel(domain,partition) order = 2 +reffe = ReferenceFE(lagrangian,Float64,order) V0 = FESpace( - reffe=lagrangian, order=order, valuetype=Float64, - conformity=:H1, model=model, dirichlet_tags="boundary") + model, + reffe, + conformity=:H1, + dirichlet_tags="boundary" +) U = TransientTrialFESpace(V0,u) -trian = Triangulation(model) +Ω = Triangulation(model) degree = 2*order -quad = CellQuadrature(trian,degree) +dΩ = Measure(Ω,degree) # a(u,v) = ∇(v)⋅∇(u) b(v,t) = v*f(t) -res(t,(u,ut),v) = a(u,v) + ut*v - b(v,t) +res(t,(u,ut),v) = ∫( a(u,v) + ut*v - b(v,t) )dΩ jac(t,(u,ut),du,v) = a(du,v) -jac_t(t,(u,ut),dut,v) = dut*v +jac_t(t,(u,ut),dut,v) = ∫( dut*v )dΩ -t_Ω = FETerm(res,jac,jac_t,trian,quad) -op = TransientFEOperator(U,V0,t_Ω) +#t_Ω = FETerm(res,jac,jac_t,trian,quad) +op = TransientFEOperator(res,jac,jac_t,U,V0) t0 = 0.0 tF = 1.0 From 9be6d01e2c59c0ad59c28513c0c88f3ba60a4c71 Mon Sep 17 00:00:00 2001 From: Oriol Colomes Date: Thu, 9 Dec 2021 20:24:20 +0100 Subject: [PATCH 5/6] Fixing Forward Euler test --- results.vtu | Bin 1207 -> 0 bytes src/ODETools/ForwardEuler.jl | 3 +-- src/TransientFETools/TransientFEOperators.jl | 4 +++- .../ForwardEulerHeatEquationTests.jl | 7 ++----- test/TransientFEsTests/runtests.jl | 2 ++ 5 files changed, 8 insertions(+), 8 deletions(-) delete mode 100644 results.vtu diff --git a/results.vtu b/results.vtu deleted file mode 100644 index 6a70a188c9905ba97455588f2eeb6635ca4eba93..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 1207 zcmbVM%}&EG3?4rR*ojAoa)++#Ks%(Hroj*b!4AN6WUUL4+N4U;LEoB3AnDIms?kkE zN>aOy{e6xdGX9D}v;vth!Ch=$^e_VM37_$I7pE#XAL8+-OJ1MvXDkGyvKULpy8F+vm(`X6uDVPh!RpR0S zUaO`MhFW9iUWy~RQE@(3TvJ6SQqrshItv9=SH=g8H15JFMuCtKRW7D+4BQ8eUj#t{ zTBGS#2&N*61&17kgRLqGu?jD1IfqqAnnx4E9am1jLf{F`!BcF-RMv*W?)`zp-XPf( ztKgKhhUxz+%!}qUtSFh1*b&6omiJoEYK;rAex;EZ$t1EpubZ3ow0$d4$JrA@Ju1ns zx~|!h^G1!OAFVhBmwVu7gw-Wydc)*`SlhaXOv*wwP;31w>8Ea%#an0B_#U@gH2+s} qwRbvarGo>#AL)6h=W#iW?({@2=EDq>ukA}tw$4kxLRR? 0.0) + _matdata = (_matdata...,matdata_jacobian(op,t,xh,i,γ[i])) + end end matdata = vcat_matdata(_matdata) assemble_matrix_add!(A,op.assem_t, matdata) diff --git a/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl b/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl index 6cf75a3..759f014 100644 --- a/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl +++ b/test/TransientFEsTests/ForwardEulerHeatEquationTests.jl @@ -48,10 +48,9 @@ a(u,v) = ∇(v)⋅∇(u) b(v,t) = v*f(t) res(t,(u,ut),v) = ∫( a(u,v) + ut*v - b(v,t) )dΩ -jac(t,(u,ut),du,v) = a(du,v) +jac(t,(u,ut),du,v) = ∫( a(du,v) )dΩ jac_t(t,(u,ut),dut,v) = ∫( dut*v )dΩ -#t_Ω = FETerm(res,jac,jac_t,trian,quad) op = TransientFEOperator(res,jac,jac_t,U,V0) t0 = 0.0 @@ -62,8 +61,6 @@ U0 = U(0.0) uh0 = interpolate_everywhere(u(0.0),U0) ls = LUSolver() -using Gridap.Algebra: NewtonRaphsonSolver -nls = NLSolver(ls;show_trace=true,method=:newton) #linesearch=BackTracking()) ode_solver = ThetaMethod(ls,dt,θ) sol_t = solve(ode_solver,op,uh0,t0,tF) @@ -77,7 +74,7 @@ for (uh_tn, tn) in sol_t global _t_n _t_n += dt e = u(tn) - uh_tn - el2 = sqrt(sum( integrate(l2(e),trian,quad) )) + el2 = sqrt(sum( ∫(l2(e))dΩ )) @test el2 < tol end diff --git a/test/TransientFEsTests/runtests.jl b/test/TransientFEsTests/runtests.jl index 05797ad..aec4392 100644 --- a/test/TransientFEsTests/runtests.jl +++ b/test/TransientFEsTests/runtests.jl @@ -26,4 +26,6 @@ using Test @testset "HeatEquationAutoDiffTests" begin include("HeatEquationAutoDiffTests.jl") end +@testset "ForwardEulerHeatEquationTests" begin include("ForwardEulerHeatEquationTests.jl") end + end # module From 7cb1ecc0db7fe3514cd2de5caab4e7a7a2d65312 Mon Sep 17 00:00:00 2001 From: Oriol Colomes Date: Thu, 9 Dec 2021 21:19:32 +0100 Subject: [PATCH 6/6] Adding Stokes autodiff test --- test/TransientFEsTests/StokesEquationAutoDiffTests.jl | 8 ++++++-- test/TransientFEsTests/runtests.jl | 2 ++ 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/test/TransientFEsTests/StokesEquationAutoDiffTests.jl b/test/TransientFEsTests/StokesEquationAutoDiffTests.jl index 65a7f34..1876df5 100644 --- a/test/TransientFEsTests/StokesEquationAutoDiffTests.jl +++ b/test/TransientFEsTests/StokesEquationAutoDiffTests.jl @@ -82,8 +82,12 @@ cell_j_t = get_array(jac_t(0.5,(xh,xh),dx,dy)) cell_j_auto = get_array(jacobian(x->res(0.5,(x,xh),dy),xh)) cell_j_t_auto = get_array(jacobian(x->res(0.5,(xh,x),dy),xh)) -test_array(cell_j_auto,cell_j,≈) -test_array(cell_j_t_auto,cell_j_t,≈) +for i in 1:length(cell_j) + test_array(cell_j[i].array[1,1],cell_j_auto[i].array[1,1],≈) + test_array(cell_j[i].array[1,2],cell_j_auto[i].array[1,2],≈) + test_array(cell_j[i].array[2,1],cell_j_auto[i].array[2,1],≈) + test_array(cell_j_t[i].array[1,1],cell_j_t_auto[i].array[1,1],≈) +end op = TransientFEOperator(res,X,Y) diff --git a/test/TransientFEsTests/runtests.jl b/test/TransientFEsTests/runtests.jl index aec4392..80dc0b1 100644 --- a/test/TransientFEsTests/runtests.jl +++ b/test/TransientFEsTests/runtests.jl @@ -26,6 +26,8 @@ using Test @testset "HeatEquationAutoDiffTests" begin include("HeatEquationAutoDiffTests.jl") end +@testset "StokesEquationAutoDiffTests" begin include("StokesEquationAutoDiffTests.jl") end + @testset "ForwardEulerHeatEquationTests" begin include("ForwardEulerHeatEquationTests.jl") end end # module