From cbc75cd13dd1dedda203cd1c03bc4bcbe17e08e5 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Wed, 22 Feb 2023 12:55:30 +0100 Subject: [PATCH 01/19] Adding MommentFitting machinery by @ericneiva Co-authored-by: Eric Neiva --- src/GridapEmbedded.jl | 2 + src/MomentFitting/CutCellMoments.jl | 419 +++++++++++++++++++++ src/MomentFitting/JacobiPolynomialBases.jl | 360 ++++++++++++++++++ src/MomentFitting/MomentFitting.jl | 63 ++++ test/MomentFitting/MomentFittingTests.jl | 94 +++++ test/MomentFitting/runtests.jl | 7 + test/runtests.jl | 3 + 7 files changed, 948 insertions(+) create mode 100644 src/MomentFitting/CutCellMoments.jl create mode 100644 src/MomentFitting/JacobiPolynomialBases.jl create mode 100644 src/MomentFitting/MomentFitting.jl create mode 100644 test/MomentFitting/MomentFittingTests.jl create mode 100644 test/MomentFitting/runtests.jl diff --git a/src/GridapEmbedded.jl b/src/GridapEmbedded.jl index deaeaf6..d815f25 100644 --- a/src/GridapEmbedded.jl +++ b/src/GridapEmbedded.jl @@ -8,6 +8,8 @@ include("LevelSetCutters/LevelSetCutters.jl") include("AgFEM/AgFEM.jl") +include("MomentFitting/MomentFitting.jl") + include("Exports.jl") end # module diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFitting/CutCellMoments.jl new file mode 100644 index 0000000..27bce17 --- /dev/null +++ b/src/MomentFitting/CutCellMoments.jl @@ -0,0 +1,419 @@ + +struct CutCellMoments + data::Vector{Vector{Float64}} + bgcell_to_cut_cell::Vector{Int32} +end + +function CutCellMoments(trian::Triangulation, + facet_moments::DomainContribution) + fi = [ testitem(array) for (trian,array) in facet_moments.dict ] + li = map(length,fi) + @assert all(li .== first(li)) + bgmodel = get_background_model(trian) + Dm = num_dims(bgmodel) + cell_to_parent_cell = get_glue(trian,Val{Dm}()).tface_to_mface + data = [ zero(first(fi)) for i in 1:length(cell_to_parent_cell) ] + bgcell_to_cut_cell = zeros(Int32,num_cells(bgmodel)) + bgcell_to_cut_cell[cell_to_parent_cell] .= 1:length(cell_to_parent_cell) + CutCellMoments(data,bgcell_to_cut_cell) +end + +function MomentFittingMeasures(cut,degree::Int) + MomentFittingMeasures(cut,cut.geo,degree) +end + +function MomentFittingMeasures(cut,in_or_out,degree::Int) + MomentFittingMeasures(cut,cut.geo,in_or_out,degree) +end + +function MomentFittingMeasures(cut,geo::CSG.Geometry,degree::Int) + MomentFittingMeasures(cut,cut.geo,IN,degree) +end + +function MomentFittingMeasures(cut, + geo::CSG.Geometry, + in_or_out, + degree::Int) + + Ωᶜ = Triangulation(cut,CUT,geo) + Ωⁱ = Triangulation(cut,in_or_out,geo) + + ccell_to_point_vals, ccell_to_weight_vals = # + compute_lag_moments_from_leg(cut,cut.geo,in_or_out,degree) + ccell_to_weight_vals = collect(get_array(ccell_to_weight_vals)) + + nq = num_cells(Ωᶜ) + ptrs = collect(1:nq) + ccell_to_point = Fill(ccell_to_point_vals,nq) + ccell_to_weight = CompressedArray(ccell_to_weight_vals,ptrs) + ccell_to_quad = map(1:nq) do i + GenericQuadrature(ccell_to_point[i],ccell_to_weight[i]) + end + + dΩᶜ = CellQuadrature( # + ccell_to_quad,ccell_to_point,ccell_to_weight,Ωᶜ,ReferenceDomain()) + dΩⁱ = CellQuadrature(Ωⁱ,degree) + Measure(dΩᶜ), Measure(dΩⁱ), Measure(lazy_append(dΩᶜ,dΩⁱ)) + +end + +function MomentFittingQuad(active_mesh::Triangulation, + cut, + degree::Int) + MomentFittingQuad(active_mesh,cut,cut.geo,degree) +end + +function MomentFittingQuad(active_mesh::Triangulation, + cut, + in_or_out, + degree::Int) + MomentFittingQuad(active_mesh,cut,cut.geo,in_or_out,degree) +end + +function MomentFittingQuad(active_mesh::Triangulation, + cut, + geo::CSG.Geometry, + degree::Int) + + MomentFittingQuad(active_mesh,cut,geo,IN,degree) +end + +function MomentFittingQuad(active_mesh::Triangulation, + cut, + geo::CSG.Geometry, + in_or_out, + degree::Int) + + acell_to_point_vals, acell_to_weight_vals = # + compute_lag_moments_from_leg(cut,geo,in_or_out,degree) + acell_to_weight_vals = collect(get_array(acell_to_weight_vals)) + + D = num_dims(active_mesh) + bgcell_to_inoutcut = compute_bgcell_to_inoutcut(cut,geo) + acell_to_bgcell = get_glue(active_mesh,Val{D}()).tface_to_mface + acell_to_inoutcut = lazy_map(Reindex(bgcell_to_inoutcut),acell_to_bgcell) + acell_to_point_ptrs = lazy_map(i->(i == CUT ? 1 : 2),acell_to_inoutcut) + + quad = map(r->Quadrature(get_polytope(r),degree),get_reffes(active_mesh)) + @assert length(quad) == 1 + acell_to_point_vals = [acell_to_point_vals,get_coordinates(quad[1])] + + push!(acell_to_weight_vals,get_weights(quad[1])) + + acell_to_is_cut = findall(lazy_map(i->(i == CUT),acell_to_inoutcut)) + num_quads = length(acell_to_weight_vals) + acell_to_weight_ptrs = map(acell_to_inoutcut) do i + i == in_or_out ? num_quads : 0 + end + acell_to_weight_ptrs[acell_to_is_cut] .= 1:length(acell_to_is_cut) + + acell_to_point = CompressedArray(acell_to_point_vals,acell_to_point_ptrs) + acell_to_weight = CompressedArray(acell_to_weight_vals,acell_to_weight_ptrs) + acell_to_quad = map(1:length(acell_to_point)) do i + GenericQuadrature(acell_to_point[i],acell_to_weight[i]) + end + + CellQuadrature( # + acell_to_quad,acell_to_point,acell_to_weight,active_mesh,ReferenceDomain()) +end + +# function compute_lag_moments(cut::EmbeddedDiscretization{D,T}, +# deg::Int) where{D,T} +# t = Triangulation(cut,cut.geo,CUT_IN) +# b = JacobiPolynomialBasis{D}(T,deg) +# p = check_and_get_polytope(cut) +# orders = tfill(deg,Val{D}()) +# nodes, _ = compute_nodes(p,orders) +# dofs = LagrangianDofBasis(T,nodes) +# change = evaluate(dofs,b) +# println(cond(change)) +# rtol = sqrt(eps(real(float(one(eltype(change)))))) +# change = pinv(change,rtol=rtol) +# l = linear_combination(change,b) +# v = Fill(l,num_cells(t)) +# dt = CellQuadrature(t,deg*D) +# x_gp_ref_1d = dt.cell_point +# cell_map = get_cell_ref_map(t) +# x_gp_ref = lazy_map(evaluate,cell_map,x_gp_ref_1d) +# v_gp_ref = lazy_map(evaluate,v,x_gp_ref) +# cell_Jt = lazy_map(∇,cell_map) +# cell_Jtx = lazy_map(evaluate,cell_Jt,x_gp_ref_1d) +# I_v_in_t = lazy_map(IntegrationMap(),v_gp_ref,dt.cell_weight,cell_Jtx) +# cbgm = DiscreteModel(cut,cut.geo,CUT) +# moments = [ zero(first(I_v_in_t)) for i in 1:num_cells(cbgm) ] +# cell_to_bgcell = get_cell_to_bgcell(t) +# cell_to_parent_cell = get_cell_to_parent_cell(cbgm) +# bgcell_to_cell = zeros(Int32,num_cells(get_parent_model(cbgm))) +# bgcell_to_cell[cell_to_parent_cell] .= 1:length(cell_to_parent_cell) +# for i in 1:num_cells(t) +# moments[bgcell_to_cell[cell_to_bgcell[i]]] += I_v_in_t[i] +# end +# nodes, moments +# end + +function Pᵢ(i::Int) + P = [] + a = (-1)^i + for k in 0:i + push!(P,a*binomial(i,k)*binomial(i+k,k)*(-1)^k) + end + P +end + +function legendreToMonomial1D(n::Int) + B = zeros(n+1,n+1) + for i in 1:n+1 + B[i,1:i] = sqrt(2*i-1)*Pᵢ(i-1) + end + B +end + +function legendreToMonomial(n::Int,d::Int) + nt = ntuple(i->1:(n+1),d) + cis = CartesianIndices(nt) + B = zeros(length(cis),length(cis)) + B1D = legendreToMonomial1D(n) + for (i,ci) in enumerate(cis) + ti = [ B1D[j,:] for j in Tuple(ci) ] + B[i,:] = kron(ti[end:-1:1]...) + end + B +end + +function compute_lag_moments_from_leg(cut, + geo::CSG.Geometry, + in_or_out, + degree::Int) + cut_trian = Triangulation(cut,CUT,geo) + T = eltype(eltype(get_node_coordinates(cut_trian))) + D = num_dims(cut_trian) + b = MonomialBasis{D}(T,degree) + mon_contribs = compute_monomial_domain_contribution(cut,geo,in_or_out,b,degree) + mon_moments = compute_monomial_cut_cell_moments(cut_trian,mon_contribs,b) + mon_to_leg = Fill(legendreToMonomial(degree,D),num_cells(cut_trian)) + leg_moments = lazy_map(*,mon_to_leg,mon_moments) + p = JacobiPolynomialBasis{D}(T,degree) + lag_nodes, lag_to_leg = get_nodes_and_change_of_basis(cut_trian,cut,p,degree) + lag_moments = lazy_map(*,lag_to_leg,leg_moments) + lag_nodes, lag_moments +end + +# function compute_cell_moments(cut::EmbeddedDiscretization{D,T}, +# degree::Int) where{D,T} +# bgtrian = Triangulation(cut.bgmodel) +# b = MonomialBasis{D}(T,degree) +# cut_bgmodel = DiscreteModel(cut,cut.geo,CUT) +# mon_contribs = compute_monomial_domain_contribution(cut,b,degree) +# mon_moments = compute_monomial_cut_cell_moments(cut_bgmodel,mon_contribs,b) +# lag_nodes, lag_to_mon = get_nodes_and_change_of_basis(cut_bgmodel,cut,b,degree) +# lag_moments = lazy_map(*,lag_to_mon,mon_moments) +# lag_nodes, lag_moments, mon_moments, lag_to_mon +# end + +function compute_monomial_domain_contribution(cut, + in_or_out, + b::MonomialBasis, + deg::Int) + compute_monomial_domain_contribution(cut,cut.geo,in_or_out,b,deg) +end + +function cut_facets(cut::EmbeddedDiscretization,geo::CSG.Geometry) + cut_facets(cut.bgmodel,geo) +end + +function compute_monomial_domain_contribution(cut, + geo::CSG.Geometry, + in_or_out::Integer, + b::MonomialBasis, + deg::Int) + + cut_io = CutInOrOut(in_or_out) + dir_Γᵉ = (-1)^(in_or_out==OUT) + # Embedded facets + Γᵉ = EmbeddedBoundary(cut,geo) + # Interior fitted cut facets + Λ = GhostSkeleton(cut) + cutf = cut_facets(cut,geo) + Γᶠ = SkeletonTriangulation(Λ,cutf,cut_io,geo) + # Boundary fitted cut facets + Γᵒ = BoundaryTriangulation(cutf,cut_io) + # Interior non-cut facets + Γᵇ = SkeletonTriangulation(Λ,cutf,in_or_out,geo) + # Boundary non-cut facets + Λ = BoundaryTriangulation(cut.bgmodel) + Γᵖ = BoundaryTriangulation(Λ,cutf,in_or_out,geo) + + D = num_dims(cut.bgmodel) + @check num_cells(Γᵉ) > 0 + J = int_c_b(Γᵉ,b,deg*D)*dir_Γᵉ + + int_c_b(Γᶠ.⁺,b,deg*D) + int_c_b(Γᶠ.⁻,b,deg*D) + if num_cells(Γᵇ) > 0 + J += int_c_b(Γᵇ.⁺,b,deg) + int_c_b(Γᵇ.⁻,b,deg) + end + if num_cells(Γᵒ) > 0 + J += int_c_b(Γᵒ,b,deg*D) + end + if num_cells(Γᵖ) > 0 + J += int_c_b(Γᵖ,b,deg) + end + J + +end + +function int_c_b(t::Triangulation,b::MonomialBasis,deg::Int) + + Dm = num_dims(get_background_model(t)) + dt = CellQuadrature(t,deg) + x_gp_ref_1d = dt.cell_point + facet_map = get_glue(t,Val{Dm}()).tface_to_mface_map + x_gp_ref = lazy_map(evaluate,facet_map,x_gp_ref_1d) + + cell_map = get_cell_map(get_background_model(t)) + facet_cell = get_glue(t,Val{Dm}()).tface_to_mface + facet_cell_map = lazy_map(Reindex(cell_map),facet_cell) + facet_cell_Jt = lazy_map(∇,facet_cell_map) + facet_cell_Jtx = lazy_map(evaluate,facet_cell_Jt,x_gp_ref) + + facet_n = get_facet_normal(t) + facet_nx = lazy_map(evaluate,facet_n,x_gp_ref_1d) + facet_nx_r = lazy_map(Broadcasting(push_normal),facet_cell_Jtx,facet_nx) + c = lazy_map(Broadcasting(⋅),facet_nx_r,x_gp_ref) + + v = Fill(b,num_cells(t)) + v_gp_ref = lazy_map(evaluate,v,x_gp_ref) + c_v = map(Broadcasting(*),v_gp_ref,c) + + facet_Jt = lazy_map(∇,facet_map) + facet_Jtx = lazy_map(evaluate,facet_Jt,x_gp_ref_1d) + + I_c_v_in_t = lazy_map(IntegrationMap(),c_v,dt.cell_weight,facet_Jtx) + + cont = DomainContribution() + add_contribution!(cont,t,I_c_v_in_t) + cont + +end + +function compute_monomial_cut_cell_moments(model::Triangulation, + facet_moments::DomainContribution, + b::MonomialBasis{D,T}) where {D,T} + cut_cell_to_moments = CutCellMoments(model,facet_moments) + for (trian,array) in facet_moments.dict + add_facet_moments!(cut_cell_to_moments,trian,array) + end + o = get_terms_degrees(b) + q = 1 ./ ( D .+ o ) + [ q .* d for d in cut_cell_to_moments.data ] +end + +function add_facet_moments!(ccm::CutCellMoments,trian,array::AbstractArray) + @abstractmethod +end + +function add_facet_moments!(ccm::CutCellMoments, + trian::SubFacetTriangulation, + array::AbstractArray) + add_facet_moments!(ccm,trian.subfacets,array) +end + +function add_facet_moments!(ccm::CutCellMoments, + sfd::SubFacetData, + array::AbstractArray) + facet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),sfd.facet_to_bgcell) + for i = 1:length(facet_to_cut_cell) + ccm.data[facet_to_cut_cell[i]] += array[i] + end +end + +function add_facet_moments!(ccm::CutCellMoments, + trian::SubFacetBoundaryTriangulation, + array::AbstractArray) + if length(trian.subfacet_to_facet) > 0 + subfacet_to_bgcell = lazy_map(Reindex(trian.facets.glue.face_to_cell),trian.subfacet_to_facet) + subfacet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),subfacet_to_bgcell) + l = length(subfacet_to_cut_cell) + for i = 1:l + ccm.data[subfacet_to_cut_cell[i]] += array[i] + end + else + add_facet_moments!(ccm,trian.facets,array) + end +end + +function add_facet_moments!(ccm::CutCellMoments, + trian::BoundaryTriangulation, + array::AbstractArray) + add_facet_moments!(ccm,trian.glue,array) +end + +function add_facet_moments!(ccm::CutCellMoments, + glue::FaceToCellGlue, + array::AbstractArray) + facet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),glue.face_to_cell) + cell_to_is_cut = findall(lazy_map(i->(i>0),facet_to_cut_cell)) + facet_to_cut_cell = lazy_map(Reindex(facet_to_cut_cell),cell_to_is_cut) + l = length(facet_to_cut_cell) + for i = 1:l + ccm.data[facet_to_cut_cell[i]] += array[cell_to_is_cut[i]] + end +end + +function add_facet_moments!(ccm::CutCellMoments, + trian::Triangulation, + array::AbstractArray) + + Dp = num_point_dims(trian) + Dc = num_cell_dims(trian) + @assert Dc == Dp-1 + @assert Dp == num_dims(get_background_model(trian)) + facet_to_bgcell = get_glue(trian,Val(Dp)).tface_to_mface + facet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),facet_to_bgcell) + l = length(facet_to_cut_cell) + for i = 1:l + ccm.data[facet_to_cut_cell[i]] += array[i] + end +end + +function get_nodes_and_change_of_basis(model::Triangulation, + cut, + b, + degree::Int) + D = num_dims(model) + T = eltype(eltype(get_node_coordinates(model))) + p = check_and_get_polytope(cut) + orders = tfill(degree,Val{D}()) + nodes, _ = compute_nodes(p,orders) + dofs = LagrangianDofBasis(T,nodes) + change = evaluate(dofs,b) + change = transpose(inv(change)) + change = Fill(change,num_cells(model)) + nodes, change +end + +function map_to_ref_space!(moments::AbstractArray, + nodes::Vector{<:Point}, + model::Triangulation) + cell_map = get_cell_map(model) + cell_Jt = lazy_map(∇,cell_map) + cell_detJt = lazy_map(Operation(det),cell_Jt) + cell_nodes = Fill(nodes,num_cells(model)) + detJt = lazy_map(evaluate,cell_detJt,cell_nodes) + moments = lazy_map(Broadcasting(/),moments,detJt) +end + +@inline function check_and_get_polytope(cut) + _check_and_get_polytope(cut.bgmodel.grid) +end + +@inline function get_terms_degrees(b::MonomialBasis) + [ _get_terms_degrees(c) for c in b.terms ] +end + +function _get_terms_degrees(c::CartesianIndex) + d = 0 + for i in 1:length(c) + d += (c[i]-1) + end + d +end + diff --git a/src/MomentFitting/JacobiPolynomialBases.jl b/src/MomentFitting/JacobiPolynomialBases.jl new file mode 100644 index 0000000..4455052 --- /dev/null +++ b/src/MomentFitting/JacobiPolynomialBases.jl @@ -0,0 +1,360 @@ + +struct JacobiPolynomial <: Field end + +struct JacobiPolynomialBasis{D,T} <: AbstractVector{JacobiPolynomial} + orders::NTuple{D,Int} + terms::Vector{CartesianIndex{D}} + function JacobiPolynomialBasis{D}( + ::Type{T}, orders::NTuple{D,Int}, terms::Vector{CartesianIndex{D}}) where {D,T} + new{D,T}(orders,terms) + end +end + +@inline Base.size(a::JacobiPolynomialBasis{D,T}) where {D,T} = (length(a.terms)*num_components(T),) +@inline Base.getindex(a::JacobiPolynomialBasis,i::Integer) = JacobiPolynomial() +@inline Base.IndexStyle(::JacobiPolynomialBasis) = IndexLinear() + +function JacobiPolynomialBasis{D}( + ::Type{T}, orders::NTuple{D,Int}, filter::Function=_q_filter) where {D,T} + + terms = _define_terms(filter, orders) + JacobiPolynomialBasis{D}(T,orders,terms) +end + +function JacobiPolynomialBasis{D}( + ::Type{T}, order::Int, filter::Function=_q_filter) where {D,T} + + orders = tfill(order,Val{D}()) + JacobiPolynomialBasis{D}(T,orders,filter) +end + +# API + +function get_exponents(b::JacobiPolynomialBasis) + indexbase = 1 + [Tuple(t) .- indexbase for t in b.terms] +end + +function get_order(b::JacobiPolynomialBasis) + maximum(b.orders) +end + +function get_orders(b::JacobiPolynomialBasis) + b.orders +end + +return_type(::JacobiPolynomialBasis{D,T}) where {D,T} = T + +# Field implementation + +function return_cache(f::JacobiPolynomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T} + @assert D == length(eltype(x)) "Incorrect number of point components" + np = length(x) + ndof = length(f.terms)*num_components(T) + n = 1 + _maximum(f.orders) + r = CachedArray(zeros(T,(np,ndof))) + v = CachedArray(zeros(T,(ndof,))) + c = CachedArray(zeros(eltype(T),(D,n))) + (r, v, c) +end + +function evaluate!(cache,f::JacobiPolynomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T} + r, v, c = cache + np = length(x) + ndof = length(f.terms)*num_components(T) + n = 1 + _maximum(f.orders) + setsize!(r,(np,ndof)) + setsize!(v,(ndof,)) + setsize!(c,(D,n)) + for i in 1:np + @inbounds xi = x[i] + _evaluate_nd_jp!(v,xi,f.orders,f.terms,c) + for j in 1:ndof + @inbounds r[i,j] = v[j] + end + end + r.array +end + +function return_cache( + fg::FieldGradientArray{1,JacobiPolynomialBasis{D,V}}, + x::AbstractVector{<:Point}) where {D,V} + + f = fg.fa + @assert D == length(eltype(x)) "Incorrect number of point components" + np = length(x) + ndof = length(f.terms)*num_components(V) + xi = testitem(x) + T = gradient_type(V,xi) + n = 1 + _maximum(f.orders) + r = CachedArray(zeros(T,(np,ndof))) + v = CachedArray(zeros(T,(ndof,))) + c = CachedArray(zeros(eltype(T),(D,n))) + g = CachedArray(zeros(eltype(T),(D,n))) + (r, v, c, g) +end + +function evaluate!( + cache, + fg::FieldGradientArray{1,JacobiPolynomialBasis{D,T}}, + x::AbstractVector{<:Point}) where {D,T} + + f = fg.fa + r, v, c, g = cache + np = length(x) + ndof = length(f.terms) * num_components(T) + n = 1 + _maximum(f.orders) + setsize!(r,(np,ndof)) + setsize!(v,(ndof,)) + setsize!(c,(D,n)) + setsize!(g,(D,n)) + for i in 1:np + @inbounds xi = x[i] + _gradient_nd_jp!(v,xi,f.orders,f.terms,c,g,T) + for j in 1:ndof + @inbounds r[i,j] = v[j] + end + end + r.array +end + +function return_cache( + fg::FieldGradientArray{2,JacobiPolynomialBasis{D,V}}, + x::AbstractVector{<:Point}) where {D,V} + + f = fg.fa + @assert D == length(eltype(x)) "Incorrect number of point components" + np = length(x) + ndof = length(f.terms)*num_components(V) + xi = testitem(x) + T = gradient_type(gradient_type(V,xi),xi) + n = 1 + _maximum(f.orders) + r = CachedArray(zeros(T,(np,ndof))) + v = CachedArray(zeros(T,(ndof,))) + c = CachedArray(zeros(eltype(T),(D,n))) + g = CachedArray(zeros(eltype(T),(D,n))) + h = CachedArray(zeros(eltype(T),(D,n))) + (r, v, c, g, h) +end + +function evaluate!( + cache, + fg::FieldGradientArray{2,JacobiPolynomialBasis{D,T}}, + x::AbstractVector{<:Point}) where {D,T} + + f = fg.fa + r, v, c, g, h = cache + np = length(x) + ndof = length(f.terms) * num_components(T) + n = 1 + _maximum(f.orders) + setsize!(r,(np,ndof)) + setsize!(v,(ndof,)) + setsize!(c,(D,n)) + setsize!(g,(D,n)) + setsize!(h,(D,n)) + for i in 1:np + @inbounds xi = x[i] + _hessian_nd_jp!(v,xi,f.orders,f.terms,c,g,h,T) + for j in 1:ndof + @inbounds r[i,j] = v[j] + end + end + r.array +end + +# Optimizing evaluation at a single point + +function return_cache(f::JacobiPolynomialBasis{D,T},x::Point) where {D,T} + ndof = length(f.terms)*num_components(T) + r = CachedArray(zeros(T,(ndof,))) + xs = [x] + cf = return_cache(f,xs) + r, cf, xs +end + +function evaluate!(cache,f::JacobiPolynomialBasis{D,T},x::Point) where {D,T} + r, cf, xs = cache + xs[1] = x + v = evaluate!(cf,f,xs) + ndof = size(v,2) + setsize!(r,(ndof,)) + a = r.array + copyto!(a,v) + a +end + +function return_cache( + f::FieldGradientArray{N,JacobiPolynomialBasis{D,V}}, x::Point) where {N,D,V} + xs = [x] + cf = return_cache(f,xs) + v = evaluate!(cf,f,xs) + r = CachedArray(zeros(eltype(v),(size(v,2),))) + r, cf, xs +end + +function evaluate!( + cache, f::FieldGradientArray{N,JacobiPolynomialBasis{D,V}}, x::Point) where {N,D,V} + r, cf, xs = cache + xs[1] = x + v = evaluate!(cf,f,xs) + ndof = size(v,2) + setsize!(r,(ndof,)) + a = r.array + copyto!(a,v) + a +end + +# Helpers + +function _evaluate_1d_jp!(v::AbstractMatrix{T},x,order,d) where T + n = order + 1 + z = one(T) + @inbounds v[d,1] = z + if n > 1 + ξ = ( 2*x[d] - 1 ) + for i in 2:n + @inbounds v[d,i] = sqrt(2*i-1)*jacobi(ξ,i-1,0,0) + end + end +end + +function _gradient_1d_jp!(v::AbstractMatrix{T},x,order,d) where T + n = order + 1 + z = zero(T) + @inbounds v[d,1] = z + if n > 1 + ξ = ( 2*x[d] - 1 ) + for i in 2:n + @inbounds v[d,i] = sqrt(2*i-1)*i*jacobi(ξ,i-2,1,1) + end + end +end + +function _hessian_1d_jp!(v::AbstractMatrix{T},x,order,d) where T + n = order + 1 + z = zero(T) + @inbounds v[d,1] = z + if n > 1 + @inbounds v[d,2] = z + ξ = ( 2*x[d] - 1 ) + for i in 3:n + @inbounds v[d,i] = sqrt(2*i-1)*(i*(i+1)/2)*jacobi(ξ,i-3,2,2) + end + end +end + +function _evaluate_nd_jp!( + v::AbstractVector{V}, + x, + orders, + terms::AbstractVector{CartesianIndex{D}}, + c::AbstractMatrix{T}) where {V,T,D} + + dim = D + for d in 1:dim + _evaluate_1d_jp!(c,x,orders[d],d) + end + + o = one(T) + k = 1 + + for ci in terms + + s = o + for d in 1:dim + @inbounds s *= c[d,ci[d]] + end + + k = _set_value!(v,s,k) + + end + +end + +function _gradient_nd_jp!( + v::AbstractVector{G}, + x, + orders, + terms::AbstractVector{CartesianIndex{D}}, + c::AbstractMatrix{T}, + g::AbstractMatrix{T}, + ::Type{V}) where {G,T,D,V} + + dim = D + for d in 1:dim + _evaluate_1d_jp!(c,x,orders[d],d) + _gradient_1d_jp!(g,x,orders[d],d) + end + + z = zero(Mutable(VectorValue{D,T})) + o = one(T) + k = 1 + + for ci in terms + + s = z + for i in eachindex(s) + @inbounds s[i] = o + end + for q in 1:dim + for d in 1:dim + if d != q + @inbounds s[q] *= c[d,ci[d]] + else + @inbounds s[q] *= g[d,ci[d]] + end + end + end + + k = _set_gradient!(v,s,k,V) + + end + +end + +function _hessian_nd_jp!( + v::AbstractVector{G}, + x, + orders, + terms::AbstractVector{CartesianIndex{D}}, + c::AbstractMatrix{T}, + g::AbstractMatrix{T}, + h::AbstractMatrix{T}, + ::Type{V}) where {G,T,D,V} + + dim = D + for d in 1:dim + _evaluate_1d_jp!(c,x,orders[d],d) + _gradient_1d_jp!(g,x,orders[d],d) + _hessian_1d_jp!(h,x,orders[d],d) + end + + z = zero(Mutable(TensorValue{D,D,T})) + o = one(T) + k = 1 + + for ci in terms + + s = z + for i in eachindex(s) + @inbounds s[i] = o + end + for r in 1:dim + for q in 1:dim + for d in 1:dim + if d != q && d != r + @inbounds s[r,q] *= c[d,ci[d]] + elseif d == q && d ==r + @inbounds s[r,q] *= h[d,ci[d]] + else + @inbounds s[r,q] *= g[d,ci[d]] + end + end + end + end + + k = _set_gradient!(v,s,k,V) + + end + +end diff --git a/src/MomentFitting/MomentFitting.jl b/src/MomentFitting/MomentFitting.jl new file mode 100644 index 0000000..d15387f --- /dev/null +++ b/src/MomentFitting/MomentFitting.jl @@ -0,0 +1,63 @@ +module MomentFitting + +using Gridap +using Gridap.Helpers +using Gridap.Arrays +using Gridap.Fields + +using Gridap.CellData + +using Gridap.Geometry +using Gridap.Geometry: get_reffes +using Gridap.Geometry: get_cell_to_parent_cell +using Gridap.Geometry: FaceToCellGlue +using Gridap.Geometry: push_normal + +using Gridap.Polynomials +using Gridap.Polynomials: jacobi +using Gridap.Polynomials: _q_filter +using Gridap.Polynomials: _define_terms +using Gridap.Polynomials: _maximum +using Gridap.Polynomials: _set_value! + +using Gridap.ReferenceFEs: get_polytope +using Gridap.ReferenceFEs: compute_nodes +using Gridap.ReferenceFEs: LagrangianDofBasis + +using Gridap.Arrays: CompressedArray + +using Gridap.ReferenceFEs +using Gridap.ReferenceFEs: Quadrature +using Gridap.ReferenceFEs: GenericQuadrature + +using GridapEmbedded.Interfaces +using GridapEmbedded.Interfaces: SubFacetTriangulation +using GridapEmbedded.Interfaces: SubFacetBoundaryTriangulation +using GridapEmbedded.Interfaces: CutInOrOut + +using GridapEmbedded.LevelSetCutters: _check_and_get_polytope + +using GridapEmbedded.CSG + +import Gridap.Polynomials: get_exponents +import Gridap.Polynomials: get_order +import Gridap.Polynomials: get_orders + +import Gridap.Arrays: return_type +import Gridap.Arrays: return_cache +import Gridap.Arrays: evaluate! + +import GridapEmbedded.Interfaces: cut_facets + +import LinearAlgebra: pinv, cond +using FillArrays + +export MomentFittingMeasures +export MomentFittingQuad + +include("JacobiPolynomialBases.jl") + +include("CutCellMoments.jl") + +end #module + diff --git a/test/MomentFitting/MomentFittingTests.jl b/test/MomentFitting/MomentFittingTests.jl new file mode 100644 index 0000000..bb8c2a3 --- /dev/null +++ b/test/MomentFitting/MomentFittingTests.jl @@ -0,0 +1,94 @@ +module MomentFittingTest + + using Test + using Gridap + using GridapEmbedded + using GridapEmbedded.MomentFitting + + function run_test(domain,partition,geom,degree) + bgmodel = CartesianDiscreteModel(domain,partition) + cutgeo = cut(bgmodel,geom) + Ωᵃ = Triangulation(cutgeo,ACTIVE_IN,geom) + dΩᵃ = Measure(MomentFittingQuad(Ωᵃ,cutgeo,degree)) + Ωᶜ = Triangulation(cutgeo) + dΩᶜ = Measure(Ωᶜ,num_dims(bgmodel)*degree,degree) + @test sum(∫(1)dΩᵃ) - sum(∫(1)dΩᶜ) + 1 ≈ 1 + if num_dims(bgmodel) == 2 + uᵃ = CellField(x->(x[1]+x[2])^degree,Ωᵃ) + uᶜ = CellField(x->(x[1]+x[2])^degree,Ωᶜ) + else + uᵃ = CellField(x->(x[1]+x[2]+x[3])^degree,Ωᵃ) + uᶜ = CellField(x->(x[1]+x[2]+x[3])^degree,Ωᶜ) + end + @test sum(∫(uᵃ)dΩᵃ) - sum(∫(uᶜ)dΩᶜ) + 1 ≈ 1 + end + + # CASE 2D 1 + n = 6 + partition = (n,n) + domain = (0,1,0,1) + geom = disk(0.42,x0=Point(0.5,0.5)) + # println("CASE 2D 1") + for deg in 1:10 # 19 + # println(deg) + run_test(domain,partition,geom,deg) + end + + # CASE 2D 2 + n = 6 + partition = (n,n) + domain = (0,1,0,1) + geom = square(L=0.53,x0=Point(0.5,0.5)) + # println("CASE 2D 2") + for deg in 1:10 # 19 + # println(deg) + run_test(domain,partition,geom,deg) + end + + # CASE 2D 3 + n = 6 + m = 7 + partition = (n,m) + domain = (0,1,0,1) + geom = disk(0.42,x0=Point(0.5,0.5)) + # println("CASE 2D 3") + for deg in 1:10 # 19 + # println(deg) + run_test(domain,partition,geom,deg) + end + + # CASE 3D 1 + n = 6 + partition = (n,n,n) + domain = (0,1,0,1,0,1) + geom = sphere(0.42,x0=Point(0.5,0.5,0.5)) + # println("CASE 3D 1") + for deg in 1:5 # 5 + # println(deg) + run_test(domain,partition,geom,deg) + end + + # CASE 3D 2 + n = 6 + partition = (n,n,n) + domain = (0,1,0,1,0,1) + geom = cube(L=0.53,x0=Point(0.5,0.5,0.5)) + # println("CASE 3D 2") + for deg in 1:5 # 7 + # println(deg) + run_test(domain,partition,geom,deg) + end + + # CASE 3D 3 + n = 2 + partition = (n,n,n) + domain = (0,1,0,1,0,1) + geom = plane(x0=Point(0.3,0.5,0.5),v=VectorValue(-1.0,0.0,0.0)) + # println("CASE 3D 3") + for deg in 1:5 # 11 + # println(deg) + run_test(domain,partition,geom,deg) + end + +end # module + diff --git a/test/MomentFitting/runtests.jl b/test/MomentFitting/runtests.jl new file mode 100644 index 0000000..a3372a8 --- /dev/null +++ b/test/MomentFitting/runtests.jl @@ -0,0 +1,7 @@ +module MomentFittingTests + +using Test + +@testset "MomentFitting" begin include("MomentFittingTests.jl") end + +end # module diff --git a/test/runtests.jl b/test/runtests.jl index 4ede043..4ba7497 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -14,9 +14,12 @@ if MiniQhull.QHULL_WRAPPER_LOADED[] @time @testset "AgFEM" begin include("AgFEMTests/runtests.jl") end + @time @testset "MomentFitting" begin include("MomentFittingTests/runtests.jl") end + include(joinpath(@__DIR__,"..","examples","runexamples.jl")) else + @warn "MiniQhull not properly installed. Tests are not executed." end From 7e46b0ad334eea764b3e62dfcdc6ca1c722db1e1 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Wed, 22 Feb 2023 13:41:02 +0100 Subject: [PATCH 02/19] fix testing --- test/{MomentFitting => MomentFittingTests}/MomentFittingTests.jl | 0 test/{MomentFitting => MomentFittingTests}/runtests.jl | 0 2 files changed, 0 insertions(+), 0 deletions(-) rename test/{MomentFitting => MomentFittingTests}/MomentFittingTests.jl (100%) rename test/{MomentFitting => MomentFittingTests}/runtests.jl (100%) diff --git a/test/MomentFitting/MomentFittingTests.jl b/test/MomentFittingTests/MomentFittingTests.jl similarity index 100% rename from test/MomentFitting/MomentFittingTests.jl rename to test/MomentFittingTests/MomentFittingTests.jl diff --git a/test/MomentFitting/runtests.jl b/test/MomentFittingTests/runtests.jl similarity index 100% rename from test/MomentFitting/runtests.jl rename to test/MomentFittingTests/runtests.jl From 050aa8302fbf01fa6ef887b2d9f386962ebe5743 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Wed, 22 Feb 2023 13:43:31 +0100 Subject: [PATCH 03/19] clean dependencies in MomentFitting module --- src/LevelSetCutters/LevelSetCutters.jl | 4 +++ src/MomentFitting/CutCellMoments.jl | 4 --- src/MomentFitting/MomentFitting.jl | 39 ++++++-------------------- 3 files changed, 13 insertions(+), 34 deletions(-) diff --git a/src/LevelSetCutters/LevelSetCutters.jl b/src/LevelSetCutters/LevelSetCutters.jl index f12bf89..7e1b31c 100644 --- a/src/LevelSetCutters/LevelSetCutters.jl +++ b/src/LevelSetCutters/LevelSetCutters.jl @@ -98,6 +98,10 @@ function cut_facets(background::DiscreteModel,geom::DiscreteGeometry) cut_facets(cutter,background,geom) end +function cut_facets(cut::EmbeddedDiscretization,geo::CSG.Geometry) + cut_facets(cut.bgmodel,geo) +end + function _cut_ls_facets(model::DiscreteModel,geom) D = num_cell_dims(model) grid = Grid(ReferenceFE{D-1},model) diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFitting/CutCellMoments.jl index 27bce17..e490ab8 100644 --- a/src/MomentFitting/CutCellMoments.jl +++ b/src/MomentFitting/CutCellMoments.jl @@ -217,10 +217,6 @@ function compute_monomial_domain_contribution(cut, compute_monomial_domain_contribution(cut,cut.geo,in_or_out,b,deg) end -function cut_facets(cut::EmbeddedDiscretization,geo::CSG.Geometry) - cut_facets(cut.bgmodel,geo) -end - function compute_monomial_domain_contribution(cut, geo::CSG.Geometry, in_or_out::Integer, diff --git a/src/MomentFitting/MomentFitting.jl b/src/MomentFitting/MomentFitting.jl index d15387f..e00f01c 100644 --- a/src/MomentFitting/MomentFitting.jl +++ b/src/MomentFitting/MomentFitting.jl @@ -1,57 +1,36 @@ module MomentFitting using Gridap -using Gridap.Helpers using Gridap.Arrays +using Gridap.CellData using Gridap.Fields +using Gridap.Geometry +using Gridap.Helpers +using Gridap.Polynomials +using Gridap.ReferenceFEs -using Gridap.CellData +using GridapEmbedded.Interfaces +using GridapEmbedded.CSG + +using FillArrays -using Gridap.Geometry -using Gridap.Geometry: get_reffes -using Gridap.Geometry: get_cell_to_parent_cell using Gridap.Geometry: FaceToCellGlue using Gridap.Geometry: push_normal - -using Gridap.Polynomials using Gridap.Polynomials: jacobi using Gridap.Polynomials: _q_filter using Gridap.Polynomials: _define_terms using Gridap.Polynomials: _maximum using Gridap.Polynomials: _set_value! -using Gridap.ReferenceFEs: get_polytope -using Gridap.ReferenceFEs: compute_nodes -using Gridap.ReferenceFEs: LagrangianDofBasis - -using Gridap.Arrays: CompressedArray - -using Gridap.ReferenceFEs -using Gridap.ReferenceFEs: Quadrature -using Gridap.ReferenceFEs: GenericQuadrature - -using GridapEmbedded.Interfaces using GridapEmbedded.Interfaces: SubFacetTriangulation using GridapEmbedded.Interfaces: SubFacetBoundaryTriangulation using GridapEmbedded.Interfaces: CutInOrOut - using GridapEmbedded.LevelSetCutters: _check_and_get_polytope -using GridapEmbedded.CSG - -import Gridap.Polynomials: get_exponents -import Gridap.Polynomials: get_order -import Gridap.Polynomials: get_orders - import Gridap.Arrays: return_type import Gridap.Arrays: return_cache import Gridap.Arrays: evaluate! -import GridapEmbedded.Interfaces: cut_facets - -import LinearAlgebra: pinv, cond -using FillArrays - export MomentFittingMeasures export MomentFittingQuad From 131b80c7e0a2741e136554e50cfda6e5badd1ec3 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Wed, 22 Feb 2023 13:48:29 +0100 Subject: [PATCH 04/19] ignoring vim files --- .gitignore | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/.gitignore b/.gitignore index 684a8e2..97dd38a 100644 --- a/.gitignore +++ b/.gitignore @@ -6,4 +6,5 @@ /docs/build/ /docs/site/ Manifest.toml -.vscode/ \ No newline at end of file +.vscode/ +*.swp From a10dc055e9a27faa510291ef75f360f81e15b6e6 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Wed, 22 Feb 2023 14:40:58 +0100 Subject: [PATCH 05/19] code cleanings --- src/MomentFitting/CutCellMoments.jl | 85 ------------------- src/MomentFitting/MomentFitting.jl | 3 + ...FittingTests.jl => CutCellMomentsTests.jl} | 2 +- .../JacobiPolynomialBasesTests.jl | 5 ++ test/MomentFittingTests/runtests.jl | 4 +- 5 files changed, 12 insertions(+), 87 deletions(-) rename test/MomentFittingTests/{MomentFittingTests.jl => CutCellMomentsTests.jl} (98%) create mode 100644 test/MomentFittingTests/JacobiPolynomialBasesTests.jl diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFitting/CutCellMoments.jl index e490ab8..5edc19c 100644 --- a/src/MomentFitting/CutCellMoments.jl +++ b/src/MomentFitting/CutCellMoments.jl @@ -18,45 +18,6 @@ function CutCellMoments(trian::Triangulation, CutCellMoments(data,bgcell_to_cut_cell) end -function MomentFittingMeasures(cut,degree::Int) - MomentFittingMeasures(cut,cut.geo,degree) -end - -function MomentFittingMeasures(cut,in_or_out,degree::Int) - MomentFittingMeasures(cut,cut.geo,in_or_out,degree) -end - -function MomentFittingMeasures(cut,geo::CSG.Geometry,degree::Int) - MomentFittingMeasures(cut,cut.geo,IN,degree) -end - -function MomentFittingMeasures(cut, - geo::CSG.Geometry, - in_or_out, - degree::Int) - - Ωᶜ = Triangulation(cut,CUT,geo) - Ωⁱ = Triangulation(cut,in_or_out,geo) - - ccell_to_point_vals, ccell_to_weight_vals = # - compute_lag_moments_from_leg(cut,cut.geo,in_or_out,degree) - ccell_to_weight_vals = collect(get_array(ccell_to_weight_vals)) - - nq = num_cells(Ωᶜ) - ptrs = collect(1:nq) - ccell_to_point = Fill(ccell_to_point_vals,nq) - ccell_to_weight = CompressedArray(ccell_to_weight_vals,ptrs) - ccell_to_quad = map(1:nq) do i - GenericQuadrature(ccell_to_point[i],ccell_to_weight[i]) - end - - dΩᶜ = CellQuadrature( # - ccell_to_quad,ccell_to_point,ccell_to_weight,Ωᶜ,ReferenceDomain()) - dΩⁱ = CellQuadrature(Ωⁱ,degree) - Measure(dΩᶜ), Measure(dΩⁱ), Measure(lazy_append(dΩᶜ,dΩⁱ)) - -end - function MomentFittingQuad(active_mesh::Triangulation, cut, degree::Int) @@ -117,40 +78,6 @@ function MomentFittingQuad(active_mesh::Triangulation, acell_to_quad,acell_to_point,acell_to_weight,active_mesh,ReferenceDomain()) end -# function compute_lag_moments(cut::EmbeddedDiscretization{D,T}, -# deg::Int) where{D,T} -# t = Triangulation(cut,cut.geo,CUT_IN) -# b = JacobiPolynomialBasis{D}(T,deg) -# p = check_and_get_polytope(cut) -# orders = tfill(deg,Val{D}()) -# nodes, _ = compute_nodes(p,orders) -# dofs = LagrangianDofBasis(T,nodes) -# change = evaluate(dofs,b) -# println(cond(change)) -# rtol = sqrt(eps(real(float(one(eltype(change)))))) -# change = pinv(change,rtol=rtol) -# l = linear_combination(change,b) -# v = Fill(l,num_cells(t)) -# dt = CellQuadrature(t,deg*D) -# x_gp_ref_1d = dt.cell_point -# cell_map = get_cell_ref_map(t) -# x_gp_ref = lazy_map(evaluate,cell_map,x_gp_ref_1d) -# v_gp_ref = lazy_map(evaluate,v,x_gp_ref) -# cell_Jt = lazy_map(∇,cell_map) -# cell_Jtx = lazy_map(evaluate,cell_Jt,x_gp_ref_1d) -# I_v_in_t = lazy_map(IntegrationMap(),v_gp_ref,dt.cell_weight,cell_Jtx) -# cbgm = DiscreteModel(cut,cut.geo,CUT) -# moments = [ zero(first(I_v_in_t)) for i in 1:num_cells(cbgm) ] -# cell_to_bgcell = get_cell_to_bgcell(t) -# cell_to_parent_cell = get_cell_to_parent_cell(cbgm) -# bgcell_to_cell = zeros(Int32,num_cells(get_parent_model(cbgm))) -# bgcell_to_cell[cell_to_parent_cell] .= 1:length(cell_to_parent_cell) -# for i in 1:num_cells(t) -# moments[bgcell_to_cell[cell_to_bgcell[i]]] += I_v_in_t[i] -# end -# nodes, moments -# end - function Pᵢ(i::Int) P = [] a = (-1)^i @@ -198,18 +125,6 @@ function compute_lag_moments_from_leg(cut, lag_nodes, lag_moments end -# function compute_cell_moments(cut::EmbeddedDiscretization{D,T}, -# degree::Int) where{D,T} -# bgtrian = Triangulation(cut.bgmodel) -# b = MonomialBasis{D}(T,degree) -# cut_bgmodel = DiscreteModel(cut,cut.geo,CUT) -# mon_contribs = compute_monomial_domain_contribution(cut,b,degree) -# mon_moments = compute_monomial_cut_cell_moments(cut_bgmodel,mon_contribs,b) -# lag_nodes, lag_to_mon = get_nodes_and_change_of_basis(cut_bgmodel,cut,b,degree) -# lag_moments = lazy_map(*,lag_to_mon,mon_moments) -# lag_nodes, lag_moments, mon_moments, lag_to_mon -# end - function compute_monomial_domain_contribution(cut, in_or_out, b::MonomialBasis, diff --git a/src/MomentFitting/MomentFitting.jl b/src/MomentFitting/MomentFitting.jl index e00f01c..555b864 100644 --- a/src/MomentFitting/MomentFitting.jl +++ b/src/MomentFitting/MomentFitting.jl @@ -30,6 +30,9 @@ using GridapEmbedded.LevelSetCutters: _check_and_get_polytope import Gridap.Arrays: return_type import Gridap.Arrays: return_cache import Gridap.Arrays: evaluate! +import Gridap.Polynomials: get_order +import Gridap.Polynomials: get_orders +import Gridap.Polynomials: get_exponents export MomentFittingMeasures export MomentFittingQuad diff --git a/test/MomentFittingTests/MomentFittingTests.jl b/test/MomentFittingTests/CutCellMomentsTests.jl similarity index 98% rename from test/MomentFittingTests/MomentFittingTests.jl rename to test/MomentFittingTests/CutCellMomentsTests.jl index bb8c2a3..1d4ae5a 100644 --- a/test/MomentFittingTests/MomentFittingTests.jl +++ b/test/MomentFittingTests/CutCellMomentsTests.jl @@ -1,4 +1,4 @@ -module MomentFittingTest +module CutCellMomentsTests using Test using Gridap diff --git a/test/MomentFittingTests/JacobiPolynomialBasesTests.jl b/test/MomentFittingTests/JacobiPolynomialBasesTests.jl new file mode 100644 index 0000000..d5b8ee9 --- /dev/null +++ b/test/MomentFittingTests/JacobiPolynomialBasesTests.jl @@ -0,0 +1,5 @@ +module JacobiPolynomialBasisTests + +# TBD + +end # module diff --git a/test/MomentFittingTests/runtests.jl b/test/MomentFittingTests/runtests.jl index a3372a8..2d8d551 100644 --- a/test/MomentFittingTests/runtests.jl +++ b/test/MomentFittingTests/runtests.jl @@ -2,6 +2,8 @@ module MomentFittingTests using Test -@testset "MomentFitting" begin include("MomentFittingTests.jl") end +@testset "CutCellMoments" begin include("CutCellMomentsTests.jl") end + +@testset "JacobiPolynomialBases" begin include("JacobiPolynomialBasesTests.jl") end end # module From f3b027ad70f19fb0cc7064d25d9a97979ee60f73 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Fri, 24 Feb 2023 18:06:16 +0100 Subject: [PATCH 06/19] bugfix in MomemtFittingQuadrature for OUT domains --- src/MomentFitting/CutCellMoments.jl | 19 +++++++++---------- 1 file changed, 9 insertions(+), 10 deletions(-) diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFitting/CutCellMoments.jl index 5edc19c..91e2b07 100644 --- a/src/MomentFitting/CutCellMoments.jl +++ b/src/MomentFitting/CutCellMoments.jl @@ -140,26 +140,26 @@ function compute_monomial_domain_contribution(cut, cut_io = CutInOrOut(in_or_out) dir_Γᵉ = (-1)^(in_or_out==OUT) + cutf = cut_facets(cut,geo) # Embedded facets Γᵉ = EmbeddedBoundary(cut,geo) # Interior fitted cut facets - Λ = GhostSkeleton(cut) - cutf = cut_facets(cut,geo) - Γᶠ = SkeletonTriangulation(Λ,cutf,cut_io,geo) + Γᶠ = SkeletonTriangulation(cutf,cut_io,geo) # Boundary fitted cut facets - Γᵒ = BoundaryTriangulation(cutf,cut_io) + Γᵒ = BoundaryTriangulation(cutf,cut_io,geo) # Interior non-cut facets - Γᵇ = SkeletonTriangulation(Λ,cutf,in_or_out,geo) + Γᵇ = SkeletonTriangulation(cutf,in_or_out,geo) # Boundary non-cut facets - Λ = BoundaryTriangulation(cut.bgmodel) - Γᵖ = BoundaryTriangulation(Λ,cutf,in_or_out,geo) + Γᵖ = BoundaryTriangulation(cutf,in_or_out,geo) D = num_dims(cut.bgmodel) @check num_cells(Γᵉ) > 0 J = int_c_b(Γᵉ,b,deg*D)*dir_Γᵉ + - int_c_b(Γᶠ.⁺,b,deg*D) + int_c_b(Γᶠ.⁻,b,deg*D) + int_c_b(Γᶠ.⁺,b,deg*D) + + int_c_b(Γᶠ.⁻,b,deg*D) if num_cells(Γᵇ) > 0 - J += int_c_b(Γᵇ.⁺,b,deg) + int_c_b(Γᵇ.⁻,b,deg) + J += int_c_b(Γᵇ.⁺,b,deg) + + int_c_b(Γᵇ.⁻,b,deg) end if num_cells(Γᵒ) > 0 J += int_c_b(Γᵒ,b,deg*D) @@ -327,4 +327,3 @@ function _get_terms_degrees(c::CartesianIndex) end d end - From b39ffb62018ff447f50a2bb7bd7440f94a9f7d3d Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Fri, 24 Feb 2023 18:31:31 +0100 Subject: [PATCH 07/19] testing MomentFittingQuad(..., OUT ) --- .../MomentFittingTests/CutCellMomentsTests.jl | 21 +++++++++++++++---- 1 file changed, 17 insertions(+), 4 deletions(-) diff --git a/test/MomentFittingTests/CutCellMomentsTests.jl b/test/MomentFittingTests/CutCellMomentsTests.jl index 1d4ae5a..50f3d5f 100644 --- a/test/MomentFittingTests/CutCellMomentsTests.jl +++ b/test/MomentFittingTests/CutCellMomentsTests.jl @@ -5,12 +5,14 @@ module CutCellMomentsTests using GridapEmbedded using GridapEmbedded.MomentFitting - function run_test(domain,partition,geom,degree) + function run_test(domain,partition,geom,degree,in_or_out=IN) + physical = (in_or_out==IN) ? PHYSICAL_IN : PHYSICAL_OUT + active = (in_or_out==IN) ? ACTIVE_IN : ACTIVE_OUT bgmodel = CartesianDiscreteModel(domain,partition) cutgeo = cut(bgmodel,geom) - Ωᵃ = Triangulation(cutgeo,ACTIVE_IN,geom) - dΩᵃ = Measure(MomentFittingQuad(Ωᵃ,cutgeo,degree)) - Ωᶜ = Triangulation(cutgeo) + Ωᵃ = Triangulation(cutgeo,active,geom) + dΩᵃ = Measure(MomentFittingQuad(Ωᵃ,cutgeo,in_or_out,degree)) + Ωᶜ = Triangulation(cutgeo,physical) dΩᶜ = Measure(Ωᶜ,num_dims(bgmodel)*degree,degree) @test sum(∫(1)dΩᵃ) - sum(∫(1)dΩᶜ) + 1 ≈ 1 if num_dims(bgmodel) == 2 @@ -57,6 +59,17 @@ module CutCellMomentsTests run_test(domain,partition,geom,deg) end + # CASE 2D 4 + n = 6 + partition = (n,n) + domain = (0,1,0,1) + geom = disk(0.42,x0=Point(0.5,0.5)) + # println("CASE 2D 4") + for deg in 1:10 # 19 + # println(deg) + run_test(domain,partition,geom,deg,OUT) + end + # CASE 3D 1 n = 6 partition = (n,n,n) From c13f42d8a9da98775d2cf8d0983173a1bc9ae022 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Mon, 27 Feb 2023 16:32:57 +0100 Subject: [PATCH 08/19] cut/split +80 char lines --- src/AgFEM/AggregateBoundingBoxes.jl | 24 ++++++++++++++---------- 1 file changed, 14 insertions(+), 10 deletions(-) diff --git a/src/AgFEM/AggregateBoundingBoxes.jl b/src/AgFEM/AggregateBoundingBoxes.jl index ffcd7b0..20da176 100644 --- a/src/AgFEM/AggregateBoundingBoxes.jl +++ b/src/AgFEM/AggregateBoundingBoxes.jl @@ -1,6 +1,6 @@ function init_bboxes(cell_to_coords) # RMK: Assuming first node is min and last node is max of BBox - [ [cell_to_coords[c][1],cell_to_coords[c][end]] for c in 1:length(cell_to_coords) ] + [ map(x-> [ x[1], x[2] ], cell_to_coords ) ] end function init_bboxes(cell_to_coords,cut::EmbeddedDiscretization) @@ -29,12 +29,13 @@ function init_cut_bboxes(cut,ccell_to_bgcell) end function init_subcell_bboxes(cut,inscell_to_subcell) - subcell_to_points = lazy_map(Reindex(cut.subcells.cell_to_points),inscell_to_subcell) - point_to_coords = cut.subcells.point_to_coords - subcell_to_coords = lazy_map(subcell_to_points) do points + subcell_to_point = cut.subcells.cell_to_points + point_to_coords = cut.subcells.point_to_coords + inscell_to_points = lazy_map(Reindex(subcell_to_points),inscell_to_subcell) + inscell_to_coords = lazy_map(inscell_to_points) do points point_to_coords[points] end - lazy_map(compute_subcell_bbox,subcell_to_coords) + lazy_map(compute_subcell_bbox,insncell_to_coords) end function compute_subcell_bbox(subcell_to_coords) @@ -111,8 +112,9 @@ function compute_bbox_dfaces(model::DiscreteModel,cell_to_agg_bbox) D = num_dims(gt) for d = 1:D-1 dface_to_Dfaces = get_faces(gt,d,D) - d_bboxes = - [ _compute_bbox_dface(dface_to_Dfaces,cell_to_agg_bbox,face) for face in 1:num_faces(gt,d) ] + d_bboxes = map(1:num_faces(gt,d)) do face + _compute_bbox_dface(dface_to_Dfaces,cell_to_agg_bbox,face) + end bboxes = push!(bboxes,d_bboxes) end bboxes = push!(bboxes,cell_to_agg_bbox) @@ -121,8 +123,8 @@ end function _compute_bbox_dface(dface_to_Dfaces,cell_to_agg_bbox,i) cells_around_dface_i = getindex(dface_to_Dfaces,i) bboxes_around_dface_i = cell_to_agg_bbox[cells_around_dface_i] - bbmins = [ bboxes_around_dface_i[i][1].data for i in 1:length(bboxes_around_dface_i) ] - bbmaxs = [ bboxes_around_dface_i[i][2].data for i in 1:length(bboxes_around_dface_i) ] + bbmins = map( i->i[1].data, bboxes_around_dface_i ) + bbmaxs = map( i->i[2].data, bboxes_around_dface_i ) [min.(bbmins...),max.(bbmaxs...)] end @@ -130,7 +132,9 @@ function _compute_cell_to_dface_bboxes(model::DiscreteModel,dbboxes) gt = get_grid_topology(model) trian = Triangulation(model) ctc = get_cell_coordinates(trian) - bboxes = [ __compute_cell_to_dface_bboxes(gt,ctc,dbboxes,cell) for cell in 1:num_cells(model) ] + bboxes = map(1:num_cells(model)) do cell + __compute_cell_to_dface_bboxes(gt,ctc,dbboxes,cell) + end CellPoint(bboxes,trian,PhysicalDomain()).cell_ref_point end From 25103f1968403034fa8c28902b1d97821a89f46b Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Mon, 27 Feb 2023 16:34:21 +0100 Subject: [PATCH 09/19] aggregate bounding boxes for IN or OUT --- src/AgFEM/AggregateBoundingBoxes.jl | 68 +++++++++++++++++++---------- 1 file changed, 45 insertions(+), 23 deletions(-) diff --git a/src/AgFEM/AggregateBoundingBoxes.jl b/src/AgFEM/AggregateBoundingBoxes.jl index 20da176..27eb10f 100644 --- a/src/AgFEM/AggregateBoundingBoxes.jl +++ b/src/AgFEM/AggregateBoundingBoxes.jl @@ -1,25 +1,26 @@ function init_bboxes(cell_to_coords) # RMK: Assuming first node is min and last node is max of BBox - [ map(x-> [ x[1], x[2] ], cell_to_coords ) ] + map(x-> [ x[1], x[end] ], vec(cell_to_coords) ) end -function init_bboxes(cell_to_coords,cut::EmbeddedDiscretization) +function init_bboxes(cell_to_coords,cut::EmbeddedDiscretization,in_or_out) bgcell_to_cbboxes = init_bboxes(cell_to_coords) cut_bgtrian = Triangulation(cut,CUT,cut.geo) cut_bgmodel = get_active_model(cut_bgtrian) ccell_to_bgcell = get_cell_to_parent_cell(cut_bgmodel) - ccell_to_cbboxes = init_cut_bboxes(cut,ccell_to_bgcell) + ccell_to_cbboxes = init_cut_bboxes(cut,ccell_to_bgcell,in_or_out) for (cc,cbb) in enumerate(ccell_to_cbboxes) bgcell_to_cbboxes[ccell_to_bgcell[cc]] = cbb end bgcell_to_cbboxes end -function init_cut_bboxes(cut,ccell_to_bgcell) +function init_cut_bboxes(cut,ccell_to_bgcell,in_or_out) subcell_to_inout = compute_subcell_to_inout(cut,cut.geo) # WARNING! This function has a bug - subcell_is_in = lazy_map(i->i==IN,subcell_to_inout) + subcell_is_in = lazy_map(i->i==in_or_out,subcell_to_inout) + subcell_to_bgcell = cut.subcells.cell_to_bgcell inscell_to_subcell = findall(subcell_is_in) - inscell_to_bgcell = lazy_map(Reindex(cut.subcells.cell_to_bgcell),inscell_to_subcell) + inscell_to_bgcell = lazy_map(Reindex(subcell_to_bgcell),inscell_to_subcell) inscell_to_bboxes = init_subcell_bboxes(cut,inscell_to_subcell) lazy_map(ccell_to_bgcell) do bg bg_to_scs = findall(map(x->x==bg,inscell_to_bgcell)) @@ -29,13 +30,12 @@ function init_cut_bboxes(cut,ccell_to_bgcell) end function init_subcell_bboxes(cut,inscell_to_subcell) - subcell_to_point = cut.subcells.cell_to_points - point_to_coords = cut.subcells.point_to_coords + subcell_to_points = cut.subcells.cell_to_points + point_to_coords = cut.subcells.point_to_coords inscell_to_points = lazy_map(Reindex(subcell_to_points),inscell_to_subcell) - inscell_to_coords = lazy_map(inscell_to_points) do points - point_to_coords[points] - end - lazy_map(compute_subcell_bbox,insncell_to_coords) + inscell_to_coords = lazy_map( # + Broadcasting(Reindex(point_to_coords)),inscell_to_points) + lazy_map(compute_subcell_bbox,inscell_to_coords) end function compute_subcell_bbox(subcell_to_coords) @@ -89,18 +89,27 @@ function compute_cell_bboxes(trian::Triangulation,cell_to_root) root_to_agg_bbox end -function compute_cell_bboxes(model::DiscreteModelPortion,cut::EmbeddedDiscretization,cell_to_root) - compute_cell_bboxes(get_parent_model(model),cut,cell_to_root) +function compute_cell_bboxes(model::DiscreteModelPortion, + cut::EmbeddedDiscretization, + cell_to_root, + in_or_out) + compute_cell_bboxes(get_parent_model(model),cut,cell_to_root,in_or_out) end -function compute_cell_bboxes(model::DiscreteModel,cut::EmbeddedDiscretization,cell_to_root) +function compute_cell_bboxes(model::DiscreteModel, + cut::EmbeddedDiscretization, + cell_to_root, + in_or_out) trian = Triangulation(model) - compute_cell_bboxes(trian,cut,cell_to_root) + compute_cell_bboxes(trian,cut,cell_to_root,in_or_out) end -function compute_cell_bboxes(trian::Triangulation,cut::EmbeddedDiscretization,cell_to_root) +function compute_cell_bboxes(trian::Triangulation, + cut::EmbeddedDiscretization, + cell_to_root, + in_or_out) cell_to_coords = get_cell_coordinates(trian) - root_to_agg_bbox = init_bboxes(cell_to_coords,cut) + root_to_agg_bbox = init_bboxes(cell_to_coords,cut,in_or_out) compute_bboxes!(root_to_agg_bbox,cell_to_root) reset_bboxes_at_cut_cells!(root_to_agg_bbox,cell_to_root,cell_to_coords) root_to_agg_bbox @@ -165,13 +174,26 @@ function compute_cell_to_dface_bboxes(model::DiscreteModelPortion,cell_to_root) bboxes[get_cell_to_parent_cell(model)] end -function compute_cell_to_dface_bboxes(model::DiscreteModel,cut::EmbeddedDiscretization,cell_to_root) - cbboxes = compute_cell_bboxes(model,cut,cell_to_root) +function compute_cell_to_dface_bboxes(model::DiscreteModel, + cut::EmbeddedDiscretization, + cell_to_root) + compute_cell_to_dface_bboxes(model,cut,cell_to_root,IN) +end + +function compute_cell_to_dface_bboxes(model::DiscreteModel, + cut::EmbeddedDiscretization, + cell_to_root, + in_or_out) + cbboxes = compute_cell_bboxes(model,cut,cell_to_root,in_or_out) dbboxes = compute_bbox_dfaces(model,cbboxes) _compute_cell_to_dface_bboxes(model,dbboxes) end -function compute_cell_to_dface_bboxes(model::DiscreteModelPortion,cut::EmbeddedDiscretization,cell_to_root) - bboxes = compute_cell_to_dface_bboxes(get_parent_model(model),cut,cell_to_root) +function compute_cell_to_dface_bboxes(model::DiscreteModelPortion, + cut::EmbeddedDiscretization, + cell_to_root, + in_or_out) + pmodel = get_parent_model(model) + bboxes = compute_cell_to_dface_bboxes(pmodel,cut,cell_to_root,in_or_out) bboxes[get_cell_to_parent_cell(model)] -end \ No newline at end of file +end From c98da54fffb6ca8e3b57be38d57d783ba2fb6110 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Thu, 23 Mar 2023 18:50:24 +0100 Subject: [PATCH 10/19] definition of AbstractEmbeddedDiscretization --- src/Interfaces/EmbeddedDiscretizations.jl | 12 ++++++- .../EmbeddedFacetDiscretizations.jl | 12 +++++-- src/Interfaces/Interfaces.jl | 1 + src/LevelSetCutters/LevelSetCutters.jl | 4 +++ src/MomentFitting/CutCellMoments.jl | 33 +++++++++++-------- src/MomentFitting/MomentFitting.jl | 2 ++ 6 files changed, 47 insertions(+), 17 deletions(-) diff --git a/src/Interfaces/EmbeddedDiscretizations.jl b/src/Interfaces/EmbeddedDiscretizations.jl index b482905..b739078 100644 --- a/src/Interfaces/EmbeddedDiscretizations.jl +++ b/src/Interfaces/EmbeddedDiscretizations.jl @@ -1,5 +1,7 @@ -struct EmbeddedDiscretization{Dp,T} <: GridapType +abstract type AbstractEmbeddedDiscretization <: GridapType end + +struct EmbeddedDiscretization{Dp,T} <: AbstractEmbeddedDiscretization bgmodel::DiscreteModel ls_to_bgcell_to_inoutcut::Vector{Vector{Int8}} subcells::SubCellData{Dp,Dp,T} @@ -46,6 +48,14 @@ function DiscreteModel(cut::EmbeddedDiscretization,geo::CSG.Geometry,in_or_out) #DiscreteModel(cut.bgmodel,cell_list) end +function get_background_model(cut::EmbeddedDiscretization) + cut.bgmodel +end + +function get_geometry(cut::EmbeddedDiscretization) + cut.geo +end + function compute_bgcell_to_inoutcut(cut::EmbeddedDiscretization,geo::CSG.Geometry) tree = get_tree(geo) diff --git a/src/Interfaces/EmbeddedFacetDiscretizations.jl b/src/Interfaces/EmbeddedFacetDiscretizations.jl index f28179a..b8b6caa 100644 --- a/src/Interfaces/EmbeddedFacetDiscretizations.jl +++ b/src/Interfaces/EmbeddedFacetDiscretizations.jl @@ -1,5 +1,5 @@ -struct EmbeddedFacetDiscretization{Dc,Dp,T} <: GridapType +struct EmbeddedFacetDiscretization{Dc,Dp,T} <: AbstractEmbeddedDiscretization bgmodel::DiscreteModel{Dp,Dp} ls_to_facet_to_inoutcut::Vector{Vector{Int8}} subfacets::SubCellData{Dc,Dp,T} @@ -8,6 +8,14 @@ struct EmbeddedFacetDiscretization{Dc,Dp,T} <: GridapType geo::CSG.Geometry end +function get_background_model(cut::EmbeddedFacetDiscretization) + cut.bgmodel +end + +function get_geometry(cut::EmbeddedFacetDiscretization) + cut.geo +end + function SkeletonTriangulation(cut::EmbeddedFacetDiscretization) SkeletonTriangulation(cut,PHYSICAL_IN) end @@ -243,7 +251,7 @@ struct SubFacetBoundaryTriangulation{Dc,Dp,T} <: Triangulation{Dc,Dp} subfacet_to_facet_map = _setup_cell_ref_map(subfacets,subgrid) face_ref_map = lazy_map(Reindex(glue.tface_to_mface_map),subfacet_to_facet) cell_ref_map = lazy_map(∘,face_ref_map,subfacet_to_facet_map) - + new{Dc,Dp,T}( facets, subfacets, diff --git a/src/Interfaces/Interfaces.jl b/src/Interfaces/Interfaces.jl index 4d1a84a..8f351d6 100644 --- a/src/Interfaces/Interfaces.jl +++ b/src/Interfaces/Interfaces.jl @@ -8,6 +8,7 @@ using Gridap.Geometry using Gridap.CellData using Gridap.Visualization +import GridapEmbedded.CSG: get_geometry import Gridap.Geometry: UnstructuredGrid import Gridap.Geometry: BoundaryTriangulation import Gridap.Geometry: SkeletonTriangulation diff --git a/src/LevelSetCutters/LevelSetCutters.jl b/src/LevelSetCutters/LevelSetCutters.jl index 7e1b31c..aa5b147 100644 --- a/src/LevelSetCutters/LevelSetCutters.jl +++ b/src/LevelSetCutters/LevelSetCutters.jl @@ -102,6 +102,10 @@ function cut_facets(cut::EmbeddedDiscretization,geo::CSG.Geometry) cut_facets(cut.bgmodel,geo) end +function cut_facets(cut::EmbeddedFacetDiscretization,args...) + cut +end + function _cut_ls_facets(model::DiscreteModel,geom) D = num_cell_dims(model) grid = Grid(ReferenceFE{D-1},model) diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFitting/CutCellMoments.jl index 91e2b07..e7292da 100644 --- a/src/MomentFitting/CutCellMoments.jl +++ b/src/MomentFitting/CutCellMoments.jl @@ -19,20 +19,22 @@ function CutCellMoments(trian::Triangulation, end function MomentFittingQuad(active_mesh::Triangulation, - cut, + cut::AbstractEmbeddedDiscretization, degree::Int) - MomentFittingQuad(active_mesh,cut,cut.geo,degree) + geo = get_geometry(cut) + MomentFittingQuad(active_mesh,cut,geo,degree) end function MomentFittingQuad(active_mesh::Triangulation, - cut, + cut::AbstractEmbeddedDiscretization, in_or_out, degree::Int) - MomentFittingQuad(active_mesh,cut,cut.geo,in_or_out,degree) + geo = get_geometry(cut) + MomentFittingQuad(active_mesh,cut,geo,in_or_out,degree) end function MomentFittingQuad(active_mesh::Triangulation, - cut, + cut::AbstractEmbeddedDiscretization, geo::CSG.Geometry, degree::Int) @@ -40,7 +42,7 @@ function MomentFittingQuad(active_mesh::Triangulation, end function MomentFittingQuad(active_mesh::Triangulation, - cut, + cut::AbstractEmbeddedDiscretization, geo::CSG.Geometry, in_or_out, degree::Int) @@ -107,7 +109,7 @@ function legendreToMonomial(n::Int,d::Int) B end -function compute_lag_moments_from_leg(cut, +function compute_lag_moments_from_leg(cut::AbstractEmbeddedDiscretization, geo::CSG.Geometry, in_or_out, degree::Int) @@ -125,14 +127,15 @@ function compute_lag_moments_from_leg(cut, lag_nodes, lag_moments end -function compute_monomial_domain_contribution(cut, +function compute_monomial_domain_contribution(cut::AbstractEmbeddedDiscretization, in_or_out, b::MonomialBasis, deg::Int) - compute_monomial_domain_contribution(cut,cut.geo,in_or_out,b,deg) + geo = get_geometry(cut) + compute_monomial_domain_contribution(cut,geo,in_or_out,b,deg) end -function compute_monomial_domain_contribution(cut, +function compute_monomial_domain_contribution(cut::AbstractEmbeddedDiscretization, geo::CSG.Geometry, in_or_out::Integer, b::MonomialBasis, @@ -152,7 +155,7 @@ function compute_monomial_domain_contribution(cut, # Boundary non-cut facets Γᵖ = BoundaryTriangulation(cutf,in_or_out,geo) - D = num_dims(cut.bgmodel) + D = num_dims(get_background_model(cut)) @check num_cells(Γᵉ) > 0 J = int_c_b(Γᵉ,b,deg*D)*dir_Γᵉ + int_c_b(Γᶠ.⁺,b,deg*D) + @@ -286,7 +289,7 @@ function add_facet_moments!(ccm::CutCellMoments, end function get_nodes_and_change_of_basis(model::Triangulation, - cut, + cut::AbstractEmbeddedDiscretization, b, degree::Int) D = num_dims(model) @@ -312,8 +315,10 @@ function map_to_ref_space!(moments::AbstractArray, moments = lazy_map(Broadcasting(/),moments,detJt) end -@inline function check_and_get_polytope(cut) - _check_and_get_polytope(cut.bgmodel.grid) +@inline function check_and_get_polytope(cut::AbstractEmbeddedDiscretization) + bgmodel = get_background_model(cut) + grid = get_grid(bgmodel) + _check_and_get_polytope(grid) end @inline function get_terms_degrees(b::MonomialBasis) diff --git a/src/MomentFitting/MomentFitting.jl b/src/MomentFitting/MomentFitting.jl index 555b864..29d675c 100644 --- a/src/MomentFitting/MomentFitting.jl +++ b/src/MomentFitting/MomentFitting.jl @@ -25,6 +25,8 @@ using Gridap.Polynomials: _set_value! using GridapEmbedded.Interfaces: SubFacetTriangulation using GridapEmbedded.Interfaces: SubFacetBoundaryTriangulation using GridapEmbedded.Interfaces: CutInOrOut +using GridapEmbedded.Interfaces: AbstractEmbeddedDiscretization +using GridapEmbedded.Interfaces: get_geometry using GridapEmbedded.LevelSetCutters: _check_and_get_polytope import Gridap.Arrays: return_type From d16484b74cc50179f1a3627632c335537cbaa6b3 Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Fri, 31 Mar 2023 15:35:58 +0200 Subject: [PATCH 11/19] [AggregateBoundingBoxes] Remove duplicate code / Simplify extensibility --- src/AgFEM/AggregateBoundingBoxes.jl | 74 ++++++------------- .../PoissonModalC0AgFEMTests.jl | 6 +- 2 files changed, 26 insertions(+), 54 deletions(-) diff --git a/src/AgFEM/AggregateBoundingBoxes.jl b/src/AgFEM/AggregateBoundingBoxes.jl index 27eb10f..4604144 100644 --- a/src/AgFEM/AggregateBoundingBoxes.jl +++ b/src/AgFEM/AggregateBoundingBoxes.jl @@ -1,9 +1,13 @@ +function init_bboxes(cell_to_coords,args...;kwargs...) + @abstractmethod +end + function init_bboxes(cell_to_coords) # RMK: Assuming first node is min and last node is max of BBox map(x-> [ x[1], x[end] ], vec(cell_to_coords) ) end -function init_bboxes(cell_to_coords,cut::EmbeddedDiscretization,in_or_out) +function init_bboxes(cell_to_coords,cut::EmbeddedDiscretization;in_or_out=IN) bgcell_to_cbboxes = init_bboxes(cell_to_coords) cut_bgtrian = Triangulation(cut,CUT,cut.geo) cut_bgmodel = get_active_model(cut_bgtrian) @@ -16,7 +20,7 @@ function init_bboxes(cell_to_coords,cut::EmbeddedDiscretization,in_or_out) end function init_cut_bboxes(cut,ccell_to_bgcell,in_or_out) - subcell_to_inout = compute_subcell_to_inout(cut,cut.geo) # WARNING! This function has a bug + subcell_to_inout = compute_subcell_to_inout(cut,cut.geo) subcell_is_in = lazy_map(i->i==in_or_out,subcell_to_inout) subcell_to_bgcell = cut.subcells.cell_to_bgcell inscell_to_subcell = findall(subcell_is_in) @@ -72,44 +76,25 @@ function reset_bboxes_at_cut_cells!(root_to_agg_bbox, end end -function compute_cell_bboxes(model::DiscreteModelPortion,cell_to_root) - compute_cell_bboxes(get_parent_model(model),cell_to_root) -end - -function compute_cell_bboxes(model::DiscreteModel,cell_to_root) - trian = Triangulation(model) - compute_cell_bboxes(trian,cell_to_root) -end - -function compute_cell_bboxes(trian::Triangulation,cell_to_root) - cell_to_coords = get_cell_coordinates(trian) - root_to_agg_bbox = init_bboxes(cell_to_coords) - compute_bboxes!(root_to_agg_bbox,cell_to_root) - reset_bboxes_at_cut_cells!(root_to_agg_bbox,cell_to_root,cell_to_coords) - root_to_agg_bbox +function compute_cell_bboxes(model,cell_to_root,args...;kwargs...) + @abstractmethod end function compute_cell_bboxes(model::DiscreteModelPortion, - cut::EmbeddedDiscretization, - cell_to_root, - in_or_out) - compute_cell_bboxes(get_parent_model(model),cut,cell_to_root,in_or_out) + cell_to_root,args...;kwargs...) + compute_cell_bboxes(get_parent_model(model),cell_to_root,args...;kwargs...) end function compute_cell_bboxes(model::DiscreteModel, - cut::EmbeddedDiscretization, - cell_to_root, - in_or_out) + cell_to_root,args...;kwargs...) trian = Triangulation(model) - compute_cell_bboxes(trian,cut,cell_to_root,in_or_out) + compute_cell_bboxes(trian,cell_to_root,args...;kwargs...) end function compute_cell_bboxes(trian::Triangulation, - cut::EmbeddedDiscretization, - cell_to_root, - in_or_out) + cell_to_root,args...;kwargs...) cell_to_coords = get_cell_coordinates(trian) - root_to_agg_bbox = init_bboxes(cell_to_coords,cut,in_or_out) + root_to_agg_bbox = init_bboxes(cell_to_coords,args...;kwargs...) compute_bboxes!(root_to_agg_bbox,cell_to_root) reset_bboxes_at_cut_cells!(root_to_agg_bbox,cell_to_root,cell_to_coords) root_to_agg_bbox @@ -163,37 +148,20 @@ function __compute_cell_to_dface_bboxes(gt::GridTopology,ctc,dbboxes,cell::Int) cdbboxes = vcat(cdbboxes,[dbboxes[D][cell][1],dbboxes[D][cell][end]]) end -function compute_cell_to_dface_bboxes(model::DiscreteModel,cell_to_root) - cbboxes = compute_cell_bboxes(model,cell_to_root) - dbboxes = compute_bbox_dfaces(model,cbboxes) - _compute_cell_to_dface_bboxes(model,dbboxes) -end - -function compute_cell_to_dface_bboxes(model::DiscreteModelPortion,cell_to_root) - bboxes = compute_cell_to_dface_bboxes(get_parent_model(model),cell_to_root) - bboxes[get_cell_to_parent_cell(model)] +function compute_cell_to_dface_bboxes(model,cell_to_root,args...;kwargs...) + @abstractmethod end function compute_cell_to_dface_bboxes(model::DiscreteModel, - cut::EmbeddedDiscretization, - cell_to_root) - compute_cell_to_dface_bboxes(model,cut,cell_to_root,IN) -end - -function compute_cell_to_dface_bboxes(model::DiscreteModel, - cut::EmbeddedDiscretization, - cell_to_root, - in_or_out) - cbboxes = compute_cell_bboxes(model,cut,cell_to_root,in_or_out) + cell_to_root,args...;kwargs...) + cbboxes = compute_cell_bboxes(model,cell_to_root,args...;kwargs...) dbboxes = compute_bbox_dfaces(model,cbboxes) _compute_cell_to_dface_bboxes(model,dbboxes) end function compute_cell_to_dface_bboxes(model::DiscreteModelPortion, - cut::EmbeddedDiscretization, - cell_to_root, - in_or_out) + cell_to_root,args...;kwargs...) pmodel = get_parent_model(model) - bboxes = compute_cell_to_dface_bboxes(pmodel,cut,cell_to_root,in_or_out) + bboxes = compute_cell_to_dface_bboxes(pmodel,cell_to_root,args...;kwargs...) bboxes[get_cell_to_parent_cell(model)] -end +end \ No newline at end of file diff --git a/test/GridapEmbeddedTests/PoissonModalC0AgFEMTests.jl b/test/GridapEmbeddedTests/PoissonModalC0AgFEMTests.jl index 5a6e00c..2ee8759 100644 --- a/test/GridapEmbeddedTests/PoissonModalC0AgFEMTests.jl +++ b/test/GridapEmbeddedTests/PoissonModalC0AgFEMTests.jl @@ -28,7 +28,11 @@ module PoissonModalC0AgFEMTests strategy = AggregateAllCutCells() aggregates = aggregate(strategy,cutgeo,geom) - bboxes = compute_cell_to_dface_bboxes(model,cutgeo,aggregates) + + # Only for testing purposes + bboxes = compute_cell_to_dface_bboxes(model,aggregates) + bboxes = compute_cell_to_dface_bboxes(model,aggregates,cutgeo) + bboxes = compute_cell_to_dface_bboxes(model,aggregates,cutgeo;in_or_out=IN) Ω_bg = Triangulation(bgmodel) Ω = Triangulation(cutgeo,PHYSICAL) From a42556cbbdbbc5c0f09e0a18d6b8cdf947067961 Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Fri, 31 Mar 2023 17:33:27 +0200 Subject: [PATCH 12/19] [MomentFitting] Instantiate the quads a la Gridap way --- src/Exports.jl | 4 +- src/MomentFitting/CutCellMoments.jl | 82 +++++++++---------- src/MomentFitting/MomentFitting.jl | 4 +- .../MomentFittingTests/CutCellMomentsTests.jl | 5 +- 4 files changed, 45 insertions(+), 50 deletions(-) diff --git a/src/Exports.jl b/src/Exports.jl index d7439de..523876d 100644 --- a/src/Exports.jl +++ b/src/Exports.jl @@ -47,4 +47,6 @@ end @publish AgFEM AggregateCutCellsByThreshold @publish AgFEM AggregateAllCutCells @publish AgFEM compute_cell_bboxes -@publish AgFEM compute_cell_to_dface_bboxes \ No newline at end of file +@publish AgFEM compute_cell_to_dface_bboxes + +@publish MomentFitting momentfitted \ No newline at end of file diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFitting/CutCellMoments.jl index e7292da..803c1b9 100644 --- a/src/MomentFitting/CutCellMoments.jl +++ b/src/MomentFitting/CutCellMoments.jl @@ -1,51 +1,25 @@ +struct MomentFitted <: QuadratureName end +const momentfitted = MomentFitted() -struct CutCellMoments - data::Vector{Vector{Float64}} - bgcell_to_cut_cell::Vector{Int32} -end - -function CutCellMoments(trian::Triangulation, - facet_moments::DomainContribution) - fi = [ testitem(array) for (trian,array) in facet_moments.dict ] - li = map(length,fi) - @assert all(li .== first(li)) - bgmodel = get_background_model(trian) - Dm = num_dims(bgmodel) - cell_to_parent_cell = get_glue(trian,Val{Dm}()).tface_to_mface - data = [ zero(first(fi)) for i in 1:length(cell_to_parent_cell) ] - bgcell_to_cut_cell = zeros(Int32,num_cells(bgmodel)) - bgcell_to_cut_cell[cell_to_parent_cell] .= 1:length(cell_to_parent_cell) - CutCellMoments(data,bgcell_to_cut_cell) -end - -function MomentFittingQuad(active_mesh::Triangulation, - cut::AbstractEmbeddedDiscretization, - degree::Int) - geo = get_geometry(cut) - MomentFittingQuad(active_mesh,cut,geo,degree) +function Quadrature(trian::Grid,::MomentFitted,args...;kwargs...) + @abstractmethod end -function MomentFittingQuad(active_mesh::Triangulation, - cut::AbstractEmbeddedDiscretization, - in_or_out, - degree::Int) +function Quadrature(trian::Grid, + ::MomentFitted, + cut::AbstractEmbeddedDiscretization, + degree::Int; + in_or_out=IN) geo = get_geometry(cut) - MomentFittingQuad(active_mesh,cut,geo,in_or_out,degree) + Quadrature(trian,momentfitted,cut,geo,degree;in_or_out=in_or_out) end -function MomentFittingQuad(active_mesh::Triangulation, - cut::AbstractEmbeddedDiscretization, - geo::CSG.Geometry, - degree::Int) - - MomentFittingQuad(active_mesh,cut,geo,IN,degree) -end - -function MomentFittingQuad(active_mesh::Triangulation, - cut::AbstractEmbeddedDiscretization, - geo::CSG.Geometry, - in_or_out, - degree::Int) +function Quadrature(active_mesh::Grid, + ::MomentFitted, + cut::AbstractEmbeddedDiscretization, + geo::CSG.Geometry, + degree::Int; + in_or_out=IN) acell_to_point_vals, acell_to_weight_vals = # compute_lag_moments_from_leg(cut,geo,in_or_out,degree) @@ -73,11 +47,29 @@ function MomentFittingQuad(active_mesh::Triangulation, acell_to_point = CompressedArray(acell_to_point_vals,acell_to_point_ptrs) acell_to_weight = CompressedArray(acell_to_weight_vals,acell_to_weight_ptrs) acell_to_quad = map(1:length(acell_to_point)) do i - GenericQuadrature(acell_to_point[i],acell_to_weight[i]) + GenericQuadrature(acell_to_point[i],acell_to_weight[i],"Moment-fitted quadrature of degree $degree") end + acell_to_quad = CompressedArray(acell_to_quad,1:length(acell_to_quad)) + +end - CellQuadrature( # - acell_to_quad,acell_to_point,acell_to_weight,active_mesh,ReferenceDomain()) +struct CutCellMoments + data::Vector{Vector{Float64}} + bgcell_to_cut_cell::Vector{Int32} +end + +function CutCellMoments(trian::Triangulation, + facet_moments::DomainContribution) + fi = [ testitem(array) for (trian,array) in facet_moments.dict ] + li = map(length,fi) + @assert all(li .== first(li)) + bgmodel = get_background_model(trian) + Dm = num_dims(bgmodel) + cell_to_parent_cell = get_glue(trian,Val{Dm}()).tface_to_mface + data = [ zero(first(fi)) for i in 1:length(cell_to_parent_cell) ] + bgcell_to_cut_cell = zeros(Int32,num_cells(bgmodel)) + bgcell_to_cut_cell[cell_to_parent_cell] .= 1:length(cell_to_parent_cell) + CutCellMoments(data,bgcell_to_cut_cell) end function Pᵢ(i::Int) diff --git a/src/MomentFitting/MomentFitting.jl b/src/MomentFitting/MomentFitting.jl index 29d675c..497b64b 100644 --- a/src/MomentFitting/MomentFitting.jl +++ b/src/MomentFitting/MomentFitting.jl @@ -8,6 +8,7 @@ using Gridap.Geometry using Gridap.Helpers using Gridap.Polynomials using Gridap.ReferenceFEs +import Gridap.ReferenceFEs: Quadrature using GridapEmbedded.Interfaces using GridapEmbedded.CSG @@ -36,8 +37,7 @@ import Gridap.Polynomials: get_order import Gridap.Polynomials: get_orders import Gridap.Polynomials: get_exponents -export MomentFittingMeasures -export MomentFittingQuad +export momentfitted include("JacobiPolynomialBases.jl") diff --git a/test/MomentFittingTests/CutCellMomentsTests.jl b/test/MomentFittingTests/CutCellMomentsTests.jl index 50f3d5f..b47e2a0 100644 --- a/test/MomentFittingTests/CutCellMomentsTests.jl +++ b/test/MomentFittingTests/CutCellMomentsTests.jl @@ -2,8 +2,8 @@ module CutCellMomentsTests using Test using Gridap + using Gridap.ReferenceFEs using GridapEmbedded - using GridapEmbedded.MomentFitting function run_test(domain,partition,geom,degree,in_or_out=IN) physical = (in_or_out==IN) ? PHYSICAL_IN : PHYSICAL_OUT @@ -11,7 +11,8 @@ module CutCellMomentsTests bgmodel = CartesianDiscreteModel(domain,partition) cutgeo = cut(bgmodel,geom) Ωᵃ = Triangulation(cutgeo,active,geom) - dΩᵃ = Measure(MomentFittingQuad(Ωᵃ,cutgeo,in_or_out,degree)) + quad = Quadrature(momentfitted,cutgeo,degree,in_or_out=in_or_out) + dΩᵃ = Measure(Ωᵃ,quad) Ωᶜ = Triangulation(cutgeo,physical) dΩᶜ = Measure(Ωᶜ,num_dims(bgmodel)*degree,degree) @test sum(∫(1)dΩᵃ) - sum(∫(1)dΩᶜ) + 1 ≈ 1 From 082eb567d593d0eab53ac5a5839bc165e002d877 Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Mon, 3 Apr 2023 10:12:21 +0200 Subject: [PATCH 13/19] Renaming + JacobiPolynomial bases to Gridap --- src/Exports.jl | 2 +- src/GridapEmbedded.jl | 2 +- .../MomentFittedQuadratures.jl} | 32 ++ src/MomentFitting/JacobiPolynomialBases.jl | 360 ------------------ src/MomentFitting/MomentFitting.jl | 47 --- .../MomentFittedQuadraturesTests.jl} | 0 test/MomentFittedQuadraturesTests/runtests.jl | 7 + .../JacobiPolynomialBasesTests.jl | 5 - test/MomentFittingTests/runtests.jl | 9 - test/runtests.jl | 2 +- 10 files changed, 42 insertions(+), 424 deletions(-) rename src/{MomentFitting/CutCellMoments.jl => MomentFittedQuadratures/MomentFittedQuadratures.jl} (93%) delete mode 100644 src/MomentFitting/JacobiPolynomialBases.jl delete mode 100644 src/MomentFitting/MomentFitting.jl rename test/{MomentFittingTests/CutCellMomentsTests.jl => MomentFittedQuadraturesTests/MomentFittedQuadraturesTests.jl} (100%) create mode 100644 test/MomentFittedQuadraturesTests/runtests.jl delete mode 100644 test/MomentFittingTests/JacobiPolynomialBasesTests.jl delete mode 100644 test/MomentFittingTests/runtests.jl diff --git a/src/Exports.jl b/src/Exports.jl index 523876d..8a6f2a5 100644 --- a/src/Exports.jl +++ b/src/Exports.jl @@ -49,4 +49,4 @@ end @publish AgFEM compute_cell_bboxes @publish AgFEM compute_cell_to_dface_bboxes -@publish MomentFitting momentfitted \ No newline at end of file +@publish MomentFittedQuadratures momentfitted \ No newline at end of file diff --git a/src/GridapEmbedded.jl b/src/GridapEmbedded.jl index d815f25..b84e57e 100644 --- a/src/GridapEmbedded.jl +++ b/src/GridapEmbedded.jl @@ -8,7 +8,7 @@ include("LevelSetCutters/LevelSetCutters.jl") include("AgFEM/AgFEM.jl") -include("MomentFitting/MomentFitting.jl") +include("MomentFittedQuadratures/MomentFittedQuadratures.jl") include("Exports.jl") diff --git a/src/MomentFitting/CutCellMoments.jl b/src/MomentFittedQuadratures/MomentFittedQuadratures.jl similarity index 93% rename from src/MomentFitting/CutCellMoments.jl rename to src/MomentFittedQuadratures/MomentFittedQuadratures.jl index 803c1b9..00fbb70 100644 --- a/src/MomentFitting/CutCellMoments.jl +++ b/src/MomentFittedQuadratures/MomentFittedQuadratures.jl @@ -1,3 +1,32 @@ +module MomentFittedQuadratures + +using Gridap +using Gridap.Arrays +using Gridap.CellData +using Gridap.Fields +using Gridap.Geometry +using Gridap.Helpers +using Gridap.Polynomials +using Gridap.ReferenceFEs +import Gridap.ReferenceFEs: Quadrature + +using GridapEmbedded.Interfaces +using GridapEmbedded.CSG + +using FillArrays + +using Gridap.Geometry: FaceToCellGlue +using Gridap.Geometry: push_normal + +using GridapEmbedded.Interfaces: SubFacetTriangulation +using GridapEmbedded.Interfaces: SubFacetBoundaryTriangulation +using GridapEmbedded.Interfaces: CutInOrOut +using GridapEmbedded.Interfaces: AbstractEmbeddedDiscretization +using GridapEmbedded.Interfaces: get_geometry +using GridapEmbedded.LevelSetCutters: _check_and_get_polytope + +export momentfitted + struct MomentFitted <: QuadratureName end const momentfitted = MomentFitted() @@ -324,3 +353,6 @@ function _get_terms_degrees(c::CartesianIndex) end d end + +end #module + diff --git a/src/MomentFitting/JacobiPolynomialBases.jl b/src/MomentFitting/JacobiPolynomialBases.jl deleted file mode 100644 index 4455052..0000000 --- a/src/MomentFitting/JacobiPolynomialBases.jl +++ /dev/null @@ -1,360 +0,0 @@ - -struct JacobiPolynomial <: Field end - -struct JacobiPolynomialBasis{D,T} <: AbstractVector{JacobiPolynomial} - orders::NTuple{D,Int} - terms::Vector{CartesianIndex{D}} - function JacobiPolynomialBasis{D}( - ::Type{T}, orders::NTuple{D,Int}, terms::Vector{CartesianIndex{D}}) where {D,T} - new{D,T}(orders,terms) - end -end - -@inline Base.size(a::JacobiPolynomialBasis{D,T}) where {D,T} = (length(a.terms)*num_components(T),) -@inline Base.getindex(a::JacobiPolynomialBasis,i::Integer) = JacobiPolynomial() -@inline Base.IndexStyle(::JacobiPolynomialBasis) = IndexLinear() - -function JacobiPolynomialBasis{D}( - ::Type{T}, orders::NTuple{D,Int}, filter::Function=_q_filter) where {D,T} - - terms = _define_terms(filter, orders) - JacobiPolynomialBasis{D}(T,orders,terms) -end - -function JacobiPolynomialBasis{D}( - ::Type{T}, order::Int, filter::Function=_q_filter) where {D,T} - - orders = tfill(order,Val{D}()) - JacobiPolynomialBasis{D}(T,orders,filter) -end - -# API - -function get_exponents(b::JacobiPolynomialBasis) - indexbase = 1 - [Tuple(t) .- indexbase for t in b.terms] -end - -function get_order(b::JacobiPolynomialBasis) - maximum(b.orders) -end - -function get_orders(b::JacobiPolynomialBasis) - b.orders -end - -return_type(::JacobiPolynomialBasis{D,T}) where {D,T} = T - -# Field implementation - -function return_cache(f::JacobiPolynomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T} - @assert D == length(eltype(x)) "Incorrect number of point components" - np = length(x) - ndof = length(f.terms)*num_components(T) - n = 1 + _maximum(f.orders) - r = CachedArray(zeros(T,(np,ndof))) - v = CachedArray(zeros(T,(ndof,))) - c = CachedArray(zeros(eltype(T),(D,n))) - (r, v, c) -end - -function evaluate!(cache,f::JacobiPolynomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T} - r, v, c = cache - np = length(x) - ndof = length(f.terms)*num_components(T) - n = 1 + _maximum(f.orders) - setsize!(r,(np,ndof)) - setsize!(v,(ndof,)) - setsize!(c,(D,n)) - for i in 1:np - @inbounds xi = x[i] - _evaluate_nd_jp!(v,xi,f.orders,f.terms,c) - for j in 1:ndof - @inbounds r[i,j] = v[j] - end - end - r.array -end - -function return_cache( - fg::FieldGradientArray{1,JacobiPolynomialBasis{D,V}}, - x::AbstractVector{<:Point}) where {D,V} - - f = fg.fa - @assert D == length(eltype(x)) "Incorrect number of point components" - np = length(x) - ndof = length(f.terms)*num_components(V) - xi = testitem(x) - T = gradient_type(V,xi) - n = 1 + _maximum(f.orders) - r = CachedArray(zeros(T,(np,ndof))) - v = CachedArray(zeros(T,(ndof,))) - c = CachedArray(zeros(eltype(T),(D,n))) - g = CachedArray(zeros(eltype(T),(D,n))) - (r, v, c, g) -end - -function evaluate!( - cache, - fg::FieldGradientArray{1,JacobiPolynomialBasis{D,T}}, - x::AbstractVector{<:Point}) where {D,T} - - f = fg.fa - r, v, c, g = cache - np = length(x) - ndof = length(f.terms) * num_components(T) - n = 1 + _maximum(f.orders) - setsize!(r,(np,ndof)) - setsize!(v,(ndof,)) - setsize!(c,(D,n)) - setsize!(g,(D,n)) - for i in 1:np - @inbounds xi = x[i] - _gradient_nd_jp!(v,xi,f.orders,f.terms,c,g,T) - for j in 1:ndof - @inbounds r[i,j] = v[j] - end - end - r.array -end - -function return_cache( - fg::FieldGradientArray{2,JacobiPolynomialBasis{D,V}}, - x::AbstractVector{<:Point}) where {D,V} - - f = fg.fa - @assert D == length(eltype(x)) "Incorrect number of point components" - np = length(x) - ndof = length(f.terms)*num_components(V) - xi = testitem(x) - T = gradient_type(gradient_type(V,xi),xi) - n = 1 + _maximum(f.orders) - r = CachedArray(zeros(T,(np,ndof))) - v = CachedArray(zeros(T,(ndof,))) - c = CachedArray(zeros(eltype(T),(D,n))) - g = CachedArray(zeros(eltype(T),(D,n))) - h = CachedArray(zeros(eltype(T),(D,n))) - (r, v, c, g, h) -end - -function evaluate!( - cache, - fg::FieldGradientArray{2,JacobiPolynomialBasis{D,T}}, - x::AbstractVector{<:Point}) where {D,T} - - f = fg.fa - r, v, c, g, h = cache - np = length(x) - ndof = length(f.terms) * num_components(T) - n = 1 + _maximum(f.orders) - setsize!(r,(np,ndof)) - setsize!(v,(ndof,)) - setsize!(c,(D,n)) - setsize!(g,(D,n)) - setsize!(h,(D,n)) - for i in 1:np - @inbounds xi = x[i] - _hessian_nd_jp!(v,xi,f.orders,f.terms,c,g,h,T) - for j in 1:ndof - @inbounds r[i,j] = v[j] - end - end - r.array -end - -# Optimizing evaluation at a single point - -function return_cache(f::JacobiPolynomialBasis{D,T},x::Point) where {D,T} - ndof = length(f.terms)*num_components(T) - r = CachedArray(zeros(T,(ndof,))) - xs = [x] - cf = return_cache(f,xs) - r, cf, xs -end - -function evaluate!(cache,f::JacobiPolynomialBasis{D,T},x::Point) where {D,T} - r, cf, xs = cache - xs[1] = x - v = evaluate!(cf,f,xs) - ndof = size(v,2) - setsize!(r,(ndof,)) - a = r.array - copyto!(a,v) - a -end - -function return_cache( - f::FieldGradientArray{N,JacobiPolynomialBasis{D,V}}, x::Point) where {N,D,V} - xs = [x] - cf = return_cache(f,xs) - v = evaluate!(cf,f,xs) - r = CachedArray(zeros(eltype(v),(size(v,2),))) - r, cf, xs -end - -function evaluate!( - cache, f::FieldGradientArray{N,JacobiPolynomialBasis{D,V}}, x::Point) where {N,D,V} - r, cf, xs = cache - xs[1] = x - v = evaluate!(cf,f,xs) - ndof = size(v,2) - setsize!(r,(ndof,)) - a = r.array - copyto!(a,v) - a -end - -# Helpers - -function _evaluate_1d_jp!(v::AbstractMatrix{T},x,order,d) where T - n = order + 1 - z = one(T) - @inbounds v[d,1] = z - if n > 1 - ξ = ( 2*x[d] - 1 ) - for i in 2:n - @inbounds v[d,i] = sqrt(2*i-1)*jacobi(ξ,i-1,0,0) - end - end -end - -function _gradient_1d_jp!(v::AbstractMatrix{T},x,order,d) where T - n = order + 1 - z = zero(T) - @inbounds v[d,1] = z - if n > 1 - ξ = ( 2*x[d] - 1 ) - for i in 2:n - @inbounds v[d,i] = sqrt(2*i-1)*i*jacobi(ξ,i-2,1,1) - end - end -end - -function _hessian_1d_jp!(v::AbstractMatrix{T},x,order,d) where T - n = order + 1 - z = zero(T) - @inbounds v[d,1] = z - if n > 1 - @inbounds v[d,2] = z - ξ = ( 2*x[d] - 1 ) - for i in 3:n - @inbounds v[d,i] = sqrt(2*i-1)*(i*(i+1)/2)*jacobi(ξ,i-3,2,2) - end - end -end - -function _evaluate_nd_jp!( - v::AbstractVector{V}, - x, - orders, - terms::AbstractVector{CartesianIndex{D}}, - c::AbstractMatrix{T}) where {V,T,D} - - dim = D - for d in 1:dim - _evaluate_1d_jp!(c,x,orders[d],d) - end - - o = one(T) - k = 1 - - for ci in terms - - s = o - for d in 1:dim - @inbounds s *= c[d,ci[d]] - end - - k = _set_value!(v,s,k) - - end - -end - -function _gradient_nd_jp!( - v::AbstractVector{G}, - x, - orders, - terms::AbstractVector{CartesianIndex{D}}, - c::AbstractMatrix{T}, - g::AbstractMatrix{T}, - ::Type{V}) where {G,T,D,V} - - dim = D - for d in 1:dim - _evaluate_1d_jp!(c,x,orders[d],d) - _gradient_1d_jp!(g,x,orders[d],d) - end - - z = zero(Mutable(VectorValue{D,T})) - o = one(T) - k = 1 - - for ci in terms - - s = z - for i in eachindex(s) - @inbounds s[i] = o - end - for q in 1:dim - for d in 1:dim - if d != q - @inbounds s[q] *= c[d,ci[d]] - else - @inbounds s[q] *= g[d,ci[d]] - end - end - end - - k = _set_gradient!(v,s,k,V) - - end - -end - -function _hessian_nd_jp!( - v::AbstractVector{G}, - x, - orders, - terms::AbstractVector{CartesianIndex{D}}, - c::AbstractMatrix{T}, - g::AbstractMatrix{T}, - h::AbstractMatrix{T}, - ::Type{V}) where {G,T,D,V} - - dim = D - for d in 1:dim - _evaluate_1d_jp!(c,x,orders[d],d) - _gradient_1d_jp!(g,x,orders[d],d) - _hessian_1d_jp!(h,x,orders[d],d) - end - - z = zero(Mutable(TensorValue{D,D,T})) - o = one(T) - k = 1 - - for ci in terms - - s = z - for i in eachindex(s) - @inbounds s[i] = o - end - for r in 1:dim - for q in 1:dim - for d in 1:dim - if d != q && d != r - @inbounds s[r,q] *= c[d,ci[d]] - elseif d == q && d ==r - @inbounds s[r,q] *= h[d,ci[d]] - else - @inbounds s[r,q] *= g[d,ci[d]] - end - end - end - end - - k = _set_gradient!(v,s,k,V) - - end - -end diff --git a/src/MomentFitting/MomentFitting.jl b/src/MomentFitting/MomentFitting.jl deleted file mode 100644 index 497b64b..0000000 --- a/src/MomentFitting/MomentFitting.jl +++ /dev/null @@ -1,47 +0,0 @@ -module MomentFitting - -using Gridap -using Gridap.Arrays -using Gridap.CellData -using Gridap.Fields -using Gridap.Geometry -using Gridap.Helpers -using Gridap.Polynomials -using Gridap.ReferenceFEs -import Gridap.ReferenceFEs: Quadrature - -using GridapEmbedded.Interfaces -using GridapEmbedded.CSG - -using FillArrays - -using Gridap.Geometry: FaceToCellGlue -using Gridap.Geometry: push_normal -using Gridap.Polynomials: jacobi -using Gridap.Polynomials: _q_filter -using Gridap.Polynomials: _define_terms -using Gridap.Polynomials: _maximum -using Gridap.Polynomials: _set_value! - -using GridapEmbedded.Interfaces: SubFacetTriangulation -using GridapEmbedded.Interfaces: SubFacetBoundaryTriangulation -using GridapEmbedded.Interfaces: CutInOrOut -using GridapEmbedded.Interfaces: AbstractEmbeddedDiscretization -using GridapEmbedded.Interfaces: get_geometry -using GridapEmbedded.LevelSetCutters: _check_and_get_polytope - -import Gridap.Arrays: return_type -import Gridap.Arrays: return_cache -import Gridap.Arrays: evaluate! -import Gridap.Polynomials: get_order -import Gridap.Polynomials: get_orders -import Gridap.Polynomials: get_exponents - -export momentfitted - -include("JacobiPolynomialBases.jl") - -include("CutCellMoments.jl") - -end #module - diff --git a/test/MomentFittingTests/CutCellMomentsTests.jl b/test/MomentFittedQuadraturesTests/MomentFittedQuadraturesTests.jl similarity index 100% rename from test/MomentFittingTests/CutCellMomentsTests.jl rename to test/MomentFittedQuadraturesTests/MomentFittedQuadraturesTests.jl diff --git a/test/MomentFittedQuadraturesTests/runtests.jl b/test/MomentFittedQuadraturesTests/runtests.jl new file mode 100644 index 0000000..2cc8d84 --- /dev/null +++ b/test/MomentFittedQuadraturesTests/runtests.jl @@ -0,0 +1,7 @@ +module MomentFittedQuadraturesTests + +using Test + +@testset "MomentFittedQuadraturesTests" begin include("MomentFittedQuadraturesTests.jl") end + +end # module diff --git a/test/MomentFittingTests/JacobiPolynomialBasesTests.jl b/test/MomentFittingTests/JacobiPolynomialBasesTests.jl deleted file mode 100644 index d5b8ee9..0000000 --- a/test/MomentFittingTests/JacobiPolynomialBasesTests.jl +++ /dev/null @@ -1,5 +0,0 @@ -module JacobiPolynomialBasisTests - -# TBD - -end # module diff --git a/test/MomentFittingTests/runtests.jl b/test/MomentFittingTests/runtests.jl deleted file mode 100644 index 2d8d551..0000000 --- a/test/MomentFittingTests/runtests.jl +++ /dev/null @@ -1,9 +0,0 @@ -module MomentFittingTests - -using Test - -@testset "CutCellMoments" begin include("CutCellMomentsTests.jl") end - -@testset "JacobiPolynomialBases" begin include("JacobiPolynomialBasesTests.jl") end - -end # module diff --git a/test/runtests.jl b/test/runtests.jl index 765e5a4..8b1363f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -19,7 +19,7 @@ if QHULL_LOADED @time @testset "AgFEM" begin include("AgFEMTests/runtests.jl") end - @time @testset "MomentFitting" begin include("MomentFittingTests/runtests.jl") end + @time @testset "MomentFittedQuadratures" begin include("MomentFittedQuadraturesTests/runtests.jl") end include(joinpath(@__DIR__,"..","examples","runexamples.jl")) From c967f101135a825e5a1672c9d67fabd57dff6b4c Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Mon, 3 Apr 2023 10:12:53 +0200 Subject: [PATCH 14/19] Temporarily adding manifest file --- Manifest.toml | 619 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 619 insertions(+) create mode 100644 Manifest.toml diff --git a/Manifest.toml b/Manifest.toml new file mode 100644 index 0000000..881bf2d --- /dev/null +++ b/Manifest.toml @@ -0,0 +1,619 @@ +# This file is machine-generated - editing it directly is not advised + +julia_version = "1.8.3" +manifest_format = "2.0" +project_hash = "3d6e007212983dcf6fe7381d5b73ff895d20ab4b" + +[[deps.AbstractFFTs]] +deps = ["ChainRulesCore", "LinearAlgebra"] +git-tree-sha1 = "16b6dbc4cf7caee4e1e75c49485ec67b667098a0" +uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c" +version = "1.3.1" + +[[deps.AbstractTrees]] +git-tree-sha1 = "03e0550477d86222521d254b741d470ba17ea0b5" +uuid = "1520ce14-60c1-5f80-bbc7-55ef81b5835c" +version = "0.3.4" + +[[deps.Adapt]] +deps = ["LinearAlgebra", "Requires"] +git-tree-sha1 = "cc37d689f599e8df4f464b2fa3870ff7db7492ef" +uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" +version = "3.6.1" + +[[deps.ArgCheck]] +git-tree-sha1 = "a3a402a35a2f7e0b87828ccabbd5ebfbebe356b4" +uuid = "dce04be8-c92d-5529-be00-80e4d2c0e197" +version = "2.3.0" + +[[deps.ArgTools]] +uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" +version = "1.1.1" + +[[deps.ArnoldiMethod]] +deps = ["LinearAlgebra", "Random", "StaticArrays"] +git-tree-sha1 = "f87e559f87a45bece9c9ed97458d3afe98b1ebb9" +uuid = "ec485272-7323-5ecc-a04f-4719b315124d" +version = "0.1.0" + +[[deps.ArrayInterface]] +deps = ["Adapt", "LinearAlgebra", "Requires", "SnoopPrecompile", "SparseArrays", "SuiteSparse"] +git-tree-sha1 = "53bb1691fb8560633ed8c0fa11d8b0900aaa805c" +uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" +version = "7.4.2" + +[[deps.ArrayLayouts]] +deps = ["FillArrays", "LinearAlgebra", "SparseArrays"] +git-tree-sha1 = "4aff5fa660eb95c2e0deb6bcdabe4d9a96bc4667" +uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" +version = "0.8.18" + +[[deps.Artifacts]] +uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" + +[[deps.AutoHashEquals]] +git-tree-sha1 = "45bb6705d93be619b81451bb2006b7ee5d4e4453" +uuid = "15f4f7f2-30c1-5605-9d31-71845cf9641f" +version = "0.2.0" + +[[deps.BSON]] +git-tree-sha1 = "2208958832d6e1b59e49f53697483a84ca8d664e" +uuid = "fbb218c0-5317-5bc6-957e-2ee96dd4b1f0" +version = "0.3.7" + +[[deps.Base64]] +uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" + +[[deps.BlockArrays]] +deps = ["ArrayLayouts", "FillArrays", "LinearAlgebra"] +git-tree-sha1 = "3b15c61bcece7c426ea641d143c808ace3661973" +uuid = "8e7c35d0-a365-5155-bbbb-fb81a777f24e" +version = "0.16.25" + +[[deps.ChainRulesCore]] +deps = ["Compat", "LinearAlgebra", "SparseArrays"] +git-tree-sha1 = "c6d890a52d2c4d55d326439580c3b8d0875a77d9" +uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" +version = "1.15.7" + +[[deps.ChangesOfVariables]] +deps = ["ChainRulesCore", "LinearAlgebra", "Test"] +git-tree-sha1 = "485193efd2176b88e6622a39a246f8c5b600e74e" +uuid = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" +version = "0.1.6" + +[[deps.CodecZlib]] +deps = ["TranscodingStreams", "Zlib_jll"] +git-tree-sha1 = "9c209fb7536406834aa938fb149964b985de6c83" +uuid = "944b1d66-785c-5afd-91f1-9de20f533193" +version = "0.7.1" + +[[deps.Combinatorics]] +git-tree-sha1 = "08c8b6831dc00bfea825826be0bc8336fc369860" +uuid = "861a8166-3701-5b0c-9a16-15d98fcdc6aa" +version = "1.0.2" + +[[deps.CommonSubexpressions]] +deps = ["MacroTools", "Test"] +git-tree-sha1 = "7b8a93dba8af7e3b42fecabf646260105ac373f7" +uuid = "bbf7d656-a473-5ed7-a52c-81e309532950" +version = "0.3.0" + +[[deps.Compat]] +deps = ["Dates", "LinearAlgebra", "UUIDs"] +git-tree-sha1 = "7a60c856b9fa189eb34f5f8a6f6b5529b7942957" +uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" +version = "4.6.1" + +[[deps.CompilerSupportLibraries_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" +version = "0.5.2+0" + +[[deps.ConstructionBase]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "89a9db8d28102b094992472d333674bd1a83ce2a" +uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" +version = "1.5.1" + +[[deps.DataStructures]] +deps = ["Compat", "InteractiveUtils", "OrderedCollections"] +git-tree-sha1 = "d1fff3a548102f48987a52a2e0d114fa97d730f0" +uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" +version = "0.18.13" + +[[deps.Dates]] +deps = ["Printf"] +uuid = "ade2ca70-3891-5945-98fb-dc099432e06a" + +[[deps.DiffResults]] +deps = ["StaticArraysCore"] +git-tree-sha1 = "782dd5f4561f5d267313f23853baaaa4c52ea621" +uuid = "163ba53b-c6d8-5494-b064-1a9d43ac40c5" +version = "1.1.0" + +[[deps.DiffRules]] +deps = ["IrrationalConstants", "LogExpFunctions", "NaNMath", "Random", "SpecialFunctions"] +git-tree-sha1 = "a4ad7ef19d2cdc2eff57abbbe68032b1cd0bd8f8" +uuid = "b552c78f-8df3-52c6-915a-8e097449b14b" +version = "1.13.0" + +[[deps.Distances]] +deps = ["LinearAlgebra", "SparseArrays", "Statistics", "StatsAPI"] +git-tree-sha1 = "49eba9ad9f7ead780bfb7ee319f962c811c6d3b2" +uuid = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7" +version = "0.10.8" + +[[deps.Distributed]] +deps = ["Random", "Serialization", "Sockets"] +uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" + +[[deps.DocStringExtensions]] +deps = ["LibGit2"] +git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d" +uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" +version = "0.9.3" + +[[deps.Downloads]] +deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] +uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" +version = "1.6.0" + +[[deps.FFTW]] +deps = ["AbstractFFTs", "FFTW_jll", "LinearAlgebra", "MKL_jll", "Preferences", "Reexport"] +git-tree-sha1 = "f9818144ce7c8c41edf5c4c179c684d92aa4d9fe" +uuid = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341" +version = "1.6.0" + +[[deps.FFTW_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "c6033cc3892d0ef5bb9cd29b7f2f0331ea5184ea" +uuid = "f5851436-0d7a-5f13-b9de-f02708fd171a" +version = "3.3.10+0" + +[[deps.FastGaussQuadrature]] +deps = ["LinearAlgebra", "SpecialFunctions", "StaticArrays"] +git-tree-sha1 = "58d83dd5a78a36205bdfddb82b1bb67682e64487" +uuid = "442a2c76-b920-505d-bb47-c5924d526838" +version = "0.4.9" + +[[deps.FileIO]] +deps = ["Pkg", "Requires", "UUIDs"] +git-tree-sha1 = "7be5f99f7d15578798f338f5433b6c432ea8037b" +uuid = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549" +version = "1.16.0" + +[[deps.FileWatching]] +uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" + +[[deps.FillArrays]] +deps = ["LinearAlgebra", "Random", "SparseArrays", "Statistics"] +git-tree-sha1 = "7072f1e3e5a8be51d525d64f63d3ec1287ff2790" +uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" +version = "0.13.11" + +[[deps.FiniteDiff]] +deps = ["ArrayInterface", "LinearAlgebra", "Requires", "Setfield", "SparseArrays", "StaticArrays"] +git-tree-sha1 = "03fcb1c42ec905d15b305359603888ec3e65f886" +uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" +version = "2.19.0" + +[[deps.ForwardDiff]] +deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions", "StaticArrays"] +git-tree-sha1 = "00e252f4d706b3d55a8863432e742bf5717b498d" +uuid = "f6369f11-7733-5829-9624-2563aa707210" +version = "0.10.35" + +[[deps.Future]] +deps = ["Random"] +uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" + +[[deps.Gridap]] +deps = ["AbstractTrees", "BSON", "BlockArrays", "Combinatorics", "DataStructures", "DocStringExtensions", "FastGaussQuadrature", "FileIO", "FillArrays", "ForwardDiff", "JLD2", "JSON", "LineSearches", "LinearAlgebra", "NLsolve", "NearestNeighbors", "PolynomialBases", "QuadGK", "Random", "SparseArrays", "SparseMatricesCSR", "StaticArrays", "Test", "WriteVTK"] +git-tree-sha1 = "9e0dd5beb345775ee9d4e9d579ff123999d3f381" +repo-rev = "master" +repo-url = "https://github.com/gridap/Gridap.jl.git" +uuid = "56d4f2e9-7ea1-5844-9cf6-b9c51ca7ce8e" +version = "0.17.17" + +[[deps.Inflate]] +git-tree-sha1 = "5cd07aab533df5170988219191dfad0519391428" +uuid = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9" +version = "0.1.3" + +[[deps.IntelOpenMP_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "d979e54b71da82f3a65b62553da4fc3d18c9004c" +uuid = "1d5cc7b8-4909-519e-a0f8-d0f5ad9712d0" +version = "2018.0.3+2" + +[[deps.InteractiveUtils]] +deps = ["Markdown"] +uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" + +[[deps.InverseFunctions]] +deps = ["Test"] +git-tree-sha1 = "49510dfcb407e572524ba94aeae2fced1f3feb0f" +uuid = "3587e190-3f89-42d0-90ee-14403ec27112" +version = "0.1.8" + +[[deps.IrrationalConstants]] +git-tree-sha1 = "630b497eafcc20001bba38a4651b327dcfc491d2" +uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" +version = "0.2.2" + +[[deps.JLD2]] +deps = ["FileIO", "MacroTools", "Mmap", "OrderedCollections", "Pkg", "Printf", "Reexport", "Requires", "TranscodingStreams", "UUIDs"] +git-tree-sha1 = "42c17b18ced77ff0be65957a591d34f4ed57c631" +uuid = "033835bb-8acc-5ee8-8aae-3f567f8a3819" +version = "0.4.31" + +[[deps.JLLWrappers]] +deps = ["Preferences"] +git-tree-sha1 = "abc9885a7ca2052a736a600f7fa66209f96506e1" +uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210" +version = "1.4.1" + +[[deps.JSON]] +deps = ["Dates", "Mmap", "Parsers", "Unicode"] +git-tree-sha1 = "3c837543ddb02250ef42f4738347454f95079d4e" +uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6" +version = "0.21.3" + +[[deps.LazyArtifacts]] +deps = ["Artifacts", "Pkg"] +uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3" + +[[deps.LibCURL]] +deps = ["LibCURL_jll", "MozillaCACerts_jll"] +uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" +version = "0.6.3" + +[[deps.LibCURL_jll]] +deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"] +uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" +version = "7.84.0+0" + +[[deps.LibGit2]] +deps = ["Base64", "NetworkOptions", "Printf", "SHA"] +uuid = "76f85450-5226-5b5a-8eaa-529ad045b433" + +[[deps.LibSSH2_jll]] +deps = ["Artifacts", "Libdl", "MbedTLS_jll"] +uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" +version = "1.10.2+0" + +[[deps.Libdl]] +uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" + +[[deps.Libiconv_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "c7cb1f5d892775ba13767a87c7ada0b980ea0a71" +uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531" +version = "1.16.1+2" + +[[deps.LightGraphs]] +deps = ["ArnoldiMethod", "DataStructures", "Distributed", "Inflate", "LinearAlgebra", "Random", "SharedArrays", "SimpleTraits", "SparseArrays", "Statistics"] +git-tree-sha1 = "432428df5f360964040ed60418dd5601ecd240b6" +uuid = "093fc24a-ae57-5d10-9952-331d41423f4d" +version = "1.3.5" + +[[deps.LightXML]] +deps = ["Libdl", "XML2_jll"] +git-tree-sha1 = "e129d9391168c677cd4800f5c0abb1ed8cb3794f" +uuid = "9c8b4983-aa76-5018-a973-4c85ecc9e179" +version = "0.9.0" + +[[deps.LineSearches]] +deps = ["LinearAlgebra", "NLSolversBase", "NaNMath", "Parameters", "Printf"] +git-tree-sha1 = "7bbea35cec17305fc70a0e5b4641477dc0789d9d" +uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255" +version = "7.2.0" + +[[deps.LinearAlgebra]] +deps = ["Libdl", "libblastrampoline_jll"] +uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" + +[[deps.LogExpFunctions]] +deps = ["ChainRulesCore", "ChangesOfVariables", "DocStringExtensions", "InverseFunctions", "IrrationalConstants", "LinearAlgebra"] +git-tree-sha1 = "0a1b7c2863e44523180fdb3146534e265a91870b" +uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" +version = "0.3.23" + +[[deps.Logging]] +uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" + +[[deps.MKL_jll]] +deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "Pkg"] +git-tree-sha1 = "2ce8695e1e699b68702c03402672a69f54b8aca9" +uuid = "856f044c-d86e-5d09-b602-aeab76dc8ba7" +version = "2022.2.0+0" + +[[deps.MacroTools]] +deps = ["Markdown", "Random"] +git-tree-sha1 = "42324d08725e200c23d4dfb549e0d5d89dede2d2" +uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" +version = "0.5.10" + +[[deps.Markdown]] +deps = ["Base64"] +uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" + +[[deps.MbedTLS_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" +version = "2.28.0+0" + +[[deps.MiniQhull]] +deps = ["QhullMiniWrapper_jll"] +git-tree-sha1 = "9dc837d180ee49eeb7c8b77bb1c860452634b0d1" +uuid = "978d7f02-9e05-4691-894f-ae31a51d76ca" +version = "0.4.0" + +[[deps.Mmap]] +uuid = "a63ad114-7e13-5084-954f-fe012c677804" + +[[deps.MozillaCACerts_jll]] +uuid = "14a3606d-f60d-562e-9121-12d972cd8159" +version = "2022.2.1" + +[[deps.NLSolversBase]] +deps = ["DiffResults", "Distributed", "FiniteDiff", "ForwardDiff"] +git-tree-sha1 = "a0b464d183da839699f4c79e7606d9d186ec172c" +uuid = "d41bc354-129a-5804-8e4c-c37616107c6c" +version = "7.8.3" + +[[deps.NLsolve]] +deps = ["Distances", "LineSearches", "LinearAlgebra", "NLSolversBase", "Printf", "Reexport"] +git-tree-sha1 = "019f12e9a1a7880459d0173c182e6a99365d7ac1" +uuid = "2774e3e8-f4cf-5e23-947b-6d7e65073b56" +version = "4.5.1" + +[[deps.NaNMath]] +deps = ["OpenLibm_jll"] +git-tree-sha1 = "0877504529a3e5c3343c6f8b4c0381e57e4387e4" +uuid = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" +version = "1.0.2" + +[[deps.NearestNeighbors]] +deps = ["Distances", "StaticArrays"] +git-tree-sha1 = "2c3726ceb3388917602169bed973dbc97f1b51a8" +uuid = "b8a86587-4115-5ab1-83bc-aa920d37bbce" +version = "0.4.13" + +[[deps.NetworkOptions]] +uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" +version = "1.2.0" + +[[deps.OpenBLAS_jll]] +deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] +uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" +version = "0.3.20+0" + +[[deps.OpenLibm_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "05823500-19ac-5b8b-9628-191a04bc5112" +version = "0.8.1+0" + +[[deps.OpenSpecFun_jll]] +deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1" +uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e" +version = "0.5.5+0" + +[[deps.OrderedCollections]] +git-tree-sha1 = "d321bf2de576bf25ec4d3e4360faca399afca282" +uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" +version = "1.6.0" + +[[deps.Parameters]] +deps = ["OrderedCollections", "UnPack"] +git-tree-sha1 = "34c0e9ad262e5f7fc75b10a9952ca7692cfc5fbe" +uuid = "d96e819e-fc66-5662-9728-84c9c7592b0a" +version = "0.12.3" + +[[deps.Parsers]] +deps = ["Dates", "SnoopPrecompile"] +git-tree-sha1 = "478ac6c952fddd4399e71d4779797c538d0ff2bf" +uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0" +version = "2.5.8" + +[[deps.Pkg]] +deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] +uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" +version = "1.8.0" + +[[deps.PolynomialBases]] +deps = ["ArgCheck", "AutoHashEquals", "FFTW", "FastGaussQuadrature", "LinearAlgebra", "Requires", "SimpleUnPack", "SpecialFunctions"] +git-tree-sha1 = "cfbe38cc661a5a418e198cb21cdc42258a586301" +uuid = "c74db56a-226d-5e98-8bb0-a6049094aeea" +version = "0.4.17" + +[[deps.Preferences]] +deps = ["TOML"] +git-tree-sha1 = "47e5f437cc0e7ef2ce8406ce1e7e24d44915f88d" +uuid = "21216c6a-2e73-6563-6e65-726566657250" +version = "1.3.0" + +[[deps.Printf]] +deps = ["Unicode"] +uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7" + +[[deps.QhullMiniWrapper_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Qhull_jll"] +git-tree-sha1 = "607cf73c03f8a9f83b36db0b86a3a9c14179621f" +uuid = "460c41e3-6112-5d7f-b78c-b6823adb3f2d" +version = "1.0.0+1" + +[[deps.Qhull_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] +git-tree-sha1 = "238dd7e2cc577281976b9681702174850f8d4cbc" +uuid = "784f63db-0788-585a-bace-daefebcd302b" +version = "8.0.1001+0" + +[[deps.QuadGK]] +deps = ["DataStructures", "LinearAlgebra"] +git-tree-sha1 = "6ec7ac8412e83d57e313393220879ede1740f9ee" +uuid = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" +version = "2.8.2" + +[[deps.REPL]] +deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] +uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" + +[[deps.Random]] +deps = ["SHA", "Serialization"] +uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" + +[[deps.Reexport]] +git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b" +uuid = "189a3867-3050-52da-a836-e630ba90ab69" +version = "1.2.2" + +[[deps.Requires]] +deps = ["UUIDs"] +git-tree-sha1 = "838a3a4188e2ded87a4f9f184b4b0d78a1e91cb7" +uuid = "ae029012-a4dd-5104-9daa-d747884805df" +version = "1.3.0" + +[[deps.SHA]] +uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" +version = "0.7.0" + +[[deps.Serialization]] +uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b" + +[[deps.Setfield]] +deps = ["ConstructionBase", "Future", "MacroTools", "StaticArraysCore"] +git-tree-sha1 = "e2cc6d8c88613c05e1defb55170bf5ff211fbeac" +uuid = "efcf1570-3423-57d1-acb7-fd33fddbac46" +version = "1.1.1" + +[[deps.SharedArrays]] +deps = ["Distributed", "Mmap", "Random", "Serialization"] +uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383" + +[[deps.SimpleTraits]] +deps = ["InteractiveUtils", "MacroTools"] +git-tree-sha1 = "5d7e3f4e11935503d3ecaf7186eac40602e7d231" +uuid = "699a6c99-e7fa-54fc-8d76-47d257e15c1d" +version = "0.9.4" + +[[deps.SimpleUnPack]] +git-tree-sha1 = "58e6353e72cde29b90a69527e56df1b5c3d8c437" +uuid = "ce78b400-467f-4804-87d8-8f486da07d0a" +version = "1.1.0" + +[[deps.SnoopPrecompile]] +deps = ["Preferences"] +git-tree-sha1 = "e760a70afdcd461cf01a575947738d359234665c" +uuid = "66db9d55-30c0-4569-8b51-7e840670fc0c" +version = "1.0.3" + +[[deps.Sockets]] +uuid = "6462fe0b-24de-5631-8697-dd941f90decc" + +[[deps.SparseArrays]] +deps = ["LinearAlgebra", "Random"] +uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" + +[[deps.SparseMatricesCSR]] +deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"] +git-tree-sha1 = "38677ca58e80b5cad2382e5a1848f93b054ad28d" +uuid = "a0a7dd2c-ebf4-11e9-1f05-cf50bc540ca1" +version = "0.6.7" + +[[deps.SpecialFunctions]] +deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] +git-tree-sha1 = "ef28127915f4229c971eb43f3fc075dd3fe91880" +uuid = "276daf66-3868-5448-9aa4-cd146d93841b" +version = "2.2.0" + +[[deps.StaticArrays]] +deps = ["LinearAlgebra", "Random", "StaticArraysCore", "Statistics"] +git-tree-sha1 = "b8d897fe7fa688e93aef573711cb207c08c9e11e" +uuid = "90137ffa-7385-5640-81b9-e52037218182" +version = "1.5.19" + +[[deps.StaticArraysCore]] +git-tree-sha1 = "6b7ba252635a5eff6a0b0664a41ee140a1c9e72a" +uuid = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" +version = "1.4.0" + +[[deps.Statistics]] +deps = ["LinearAlgebra", "SparseArrays"] +uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" + +[[deps.StatsAPI]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "45a7769a04a3cf80da1c1c7c60caf932e6f4c9f7" +uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0" +version = "1.6.0" + +[[deps.SuiteSparse]] +deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"] +uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" + +[[deps.TOML]] +deps = ["Dates"] +uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" +version = "1.0.0" + +[[deps.Tar]] +deps = ["ArgTools", "SHA"] +uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" +version = "1.10.1" + +[[deps.Test]] +deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] +uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" + +[[deps.TranscodingStreams]] +deps = ["Random", "Test"] +git-tree-sha1 = "94f38103c984f89cf77c402f2a68dbd870f8165f" +uuid = "3bb67fe8-82b1-5028-8e26-92a6c54297fa" +version = "0.9.11" + +[[deps.UUIDs]] +deps = ["Random", "SHA"] +uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" + +[[deps.UnPack]] +git-tree-sha1 = "387c1f73762231e86e0c9c5443ce3b4a0a9a0c2b" +uuid = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" +version = "1.0.2" + +[[deps.Unicode]] +uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" + +[[deps.WriteVTK]] +deps = ["Base64", "CodecZlib", "FillArrays", "LightXML", "TranscodingStreams"] +git-tree-sha1 = "49353f30da65f377cff0f934bb9f562a2c0441b9" +uuid = "64499a7a-5c06-52f2-abe2-ccb03c286192" +version = "1.17.1" + +[[deps.XML2_jll]] +deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"] +git-tree-sha1 = "93c41695bc1c08c46c5899f4fe06d6ead504bb73" +uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a" +version = "2.10.3+0" + +[[deps.Zlib_jll]] +deps = ["Libdl"] +uuid = "83775a58-1f1d-513f-b197-d71354ab007a" +version = "1.2.12+3" + +[[deps.libblastrampoline_jll]] +deps = ["Artifacts", "Libdl", "OpenBLAS_jll"] +uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" +version = "5.1.1+0" + +[[deps.nghttp2_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" +version = "1.48.0+0" + +[[deps.p7zip_jll]] +deps = ["Artifacts", "Libdl"] +uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0" +version = "17.4.0+0" From 6c64247ca6ad59a6a8809af680aa6adf8befd527 Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Tue, 4 Apr 2023 21:27:40 +0200 Subject: [PATCH 15/19] simpler interface for int_c_b() (moment fitted contributions) --- .../MomentFittedQuadratures.jl | 81 +++++++++---------- 1 file changed, 38 insertions(+), 43 deletions(-) diff --git a/src/MomentFittedQuadratures/MomentFittedQuadratures.jl b/src/MomentFittedQuadratures/MomentFittedQuadratures.jl index 00fbb70..c1b9580 100644 --- a/src/MomentFittedQuadratures/MomentFittedQuadratures.jl +++ b/src/MomentFittedQuadratures/MomentFittedQuadratures.jl @@ -177,56 +177,51 @@ function compute_monomial_domain_contribution(cut::AbstractEmbeddedDiscretizatio Γᵖ = BoundaryTriangulation(cutf,in_or_out,geo) D = num_dims(get_background_model(cut)) - @check num_cells(Γᵉ) > 0 - J = int_c_b(Γᵉ,b,deg*D)*dir_Γᵉ + - int_c_b(Γᶠ.⁺,b,deg*D) + - int_c_b(Γᶠ.⁻,b,deg*D) - if num_cells(Γᵇ) > 0 - J += int_c_b(Γᵇ.⁺,b,deg) + - int_c_b(Γᵇ.⁻,b,deg) - end - if num_cells(Γᵒ) > 0 - J += int_c_b(Γᵒ,b,deg*D) - end - if num_cells(Γᵖ) > 0 - J += int_c_b(Γᵖ,b,deg) - end - J + int_c_b(Γᵉ,b,deg*D)*dir_Γᵉ + int_c_b(Γᶠ,b,deg*D) + int_c_b(Γᵒ,b,deg*D) + + int_c_b(Γᵇ,b,deg) + int_c_b(Γᵖ,b,deg) end -function int_c_b(t::Triangulation,b::MonomialBasis,deg::Int) - - Dm = num_dims(get_background_model(t)) - dt = CellQuadrature(t,deg) - x_gp_ref_1d = dt.cell_point - facet_map = get_glue(t,Val{Dm}()).tface_to_mface_map - x_gp_ref = lazy_map(evaluate,facet_map,x_gp_ref_1d) - - cell_map = get_cell_map(get_background_model(t)) - facet_cell = get_glue(t,Val{Dm}()).tface_to_mface - facet_cell_map = lazy_map(Reindex(cell_map),facet_cell) - facet_cell_Jt = lazy_map(∇,facet_cell_map) - facet_cell_Jtx = lazy_map(evaluate,facet_cell_Jt,x_gp_ref) - - facet_n = get_facet_normal(t) - facet_nx = lazy_map(evaluate,facet_n,x_gp_ref_1d) - facet_nx_r = lazy_map(Broadcasting(push_normal),facet_cell_Jtx,facet_nx) - c = lazy_map(Broadcasting(⋅),facet_nx_r,x_gp_ref) - - v = Fill(b,num_cells(t)) - v_gp_ref = lazy_map(evaluate,v,x_gp_ref) - c_v = map(Broadcasting(*),v_gp_ref,c) - - facet_Jt = lazy_map(∇,facet_map) - facet_Jtx = lazy_map(evaluate,facet_Jt,x_gp_ref_1d) +function int_c_b(t::SkeletonTriangulation,b::MonomialBasis,deg::Int) + int_c_b(t.plus,b,deg) + int_c_b(t.minus,b,deg) +end - I_c_v_in_t = lazy_map(IntegrationMap(),c_v,dt.cell_weight,facet_Jtx) +function int_c_b(t::AppendedTriangulation,b::MonomialBasis,deg::Int) + int_c_b(t.a,b,deg) + int_c_b(t.b,b,deg) +end +function int_c_b(t::Triangulation,b::MonomialBasis,deg::Int) cont = DomainContribution() - add_contribution!(cont,t,I_c_v_in_t) + if num_cells(t) > 0 + Dm = num_dims(get_background_model(t)) + dt = CellQuadrature(t,deg) + x_gp_ref_1d = dt.cell_point + facet_map = get_glue(t,Val{Dm}()).tface_to_mface_map + x_gp_ref = lazy_map(evaluate,facet_map,x_gp_ref_1d) + + cell_map = get_cell_map(get_background_model(t)) + facet_cell = get_glue(t,Val{Dm}()).tface_to_mface + facet_cell_map = lazy_map(Reindex(cell_map),facet_cell) + facet_cell_Jt = lazy_map(∇,facet_cell_map) + facet_cell_Jtx = lazy_map(evaluate,facet_cell_Jt,x_gp_ref) + + facet_n = get_facet_normal(t) + facet_nx = lazy_map(evaluate,facet_n,x_gp_ref_1d) + facet_nx_r = lazy_map(Broadcasting(push_normal),facet_cell_Jtx,facet_nx) + c = lazy_map(Broadcasting(⋅),facet_nx_r,x_gp_ref) + + v = Fill(b,num_cells(t)) + v_gp_ref = lazy_map(evaluate,v,x_gp_ref) + c_v = map(Broadcasting(*),v_gp_ref,c) + + facet_Jt = lazy_map(∇,facet_map) + facet_Jtx = lazy_map(evaluate,facet_Jt,x_gp_ref_1d) + + I_c_v_in_t = lazy_map(IntegrationMap(),c_v,dt.cell_weight,facet_Jtx) + + add_contribution!(cont,t,I_c_v_in_t) + end cont - end function compute_monomial_cut_cell_moments(model::Triangulation, From 1c2dc52fde059e114a9f20f67cf98db271366714 Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Fri, 7 Apr 2023 11:29:08 +0200 Subject: [PATCH 16/19] Upgrading gh actions to 1.8 --- .github/workflows/ci.yml | 6 +++--- .github/workflows/ci_x86.yml | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index e84ec91..b047a8d 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -8,7 +8,7 @@ jobs: fail-fast: false matrix: version: - - '1.6' + - '1.8' os: - ubuntu-latest arch: @@ -42,7 +42,7 @@ jobs: fail-fast: false matrix: version: - - '1.6' + - '1.8' os: - ubuntu-latest arch: @@ -77,7 +77,7 @@ jobs: - uses: actions/checkout@v2 - uses: julia-actions/setup-julia@v1 with: - version: '1.6' + version: '1.8' - run: | julia --project=docs -e ' using Pkg diff --git a/.github/workflows/ci_x86.yml b/.github/workflows/ci_x86.yml index 49a3ef7..03ca2ab 100644 --- a/.github/workflows/ci_x86.yml +++ b/.github/workflows/ci_x86.yml @@ -14,7 +14,7 @@ jobs: fail-fast: false matrix: version: - - '1.6' + - '1.8' os: - ubuntu-latest arch: From 9093e028734241c7805867721e0a4f44a8f160de Mon Sep 17 00:00:00 2001 From: Pere Antoni Martorell Date: Thu, 11 May 2023 13:07:16 +0200 Subject: [PATCH 17/19] [MomentFitting] Optimizations --- .../MomentFittedQuadratures.jl | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/src/MomentFittedQuadratures/MomentFittedQuadratures.jl b/src/MomentFittedQuadratures/MomentFittedQuadratures.jl index c1b9580..5c6f114 100644 --- a/src/MomentFittedQuadratures/MomentFittedQuadratures.jl +++ b/src/MomentFittedQuadratures/MomentFittedQuadratures.jl @@ -212,7 +212,7 @@ function int_c_b(t::Triangulation,b::MonomialBasis,deg::Int) v = Fill(b,num_cells(t)) v_gp_ref = lazy_map(evaluate,v,x_gp_ref) - c_v = map(Broadcasting(*),v_gp_ref,c) + c_v = lazy_map(Broadcasting(*),v_gp_ref,c) facet_Jt = lazy_map(∇,facet_map) facet_Jtx = lazy_map(evaluate,facet_Jt,x_gp_ref_1d) @@ -250,8 +250,9 @@ function add_facet_moments!(ccm::CutCellMoments, sfd::SubFacetData, array::AbstractArray) facet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),sfd.facet_to_bgcell) + c = array_cache(array) for i = 1:length(facet_to_cut_cell) - ccm.data[facet_to_cut_cell[i]] += array[i] + ccm.data[facet_to_cut_cell[i]] += getindex!(c,array,i) end end @@ -262,8 +263,9 @@ function add_facet_moments!(ccm::CutCellMoments, subfacet_to_bgcell = lazy_map(Reindex(trian.facets.glue.face_to_cell),trian.subfacet_to_facet) subfacet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),subfacet_to_bgcell) l = length(subfacet_to_cut_cell) + c = array_cache(array) for i = 1:l - ccm.data[subfacet_to_cut_cell[i]] += array[i] + ccm.data[subfacet_to_cut_cell[i]] += getindex!(c,array,i) end else add_facet_moments!(ccm,trian.facets,array) @@ -283,8 +285,9 @@ function add_facet_moments!(ccm::CutCellMoments, cell_to_is_cut = findall(lazy_map(i->(i>0),facet_to_cut_cell)) facet_to_cut_cell = lazy_map(Reindex(facet_to_cut_cell),cell_to_is_cut) l = length(facet_to_cut_cell) + c = array_cache(array) for i = 1:l - ccm.data[facet_to_cut_cell[i]] += array[cell_to_is_cut[i]] + ccm.data[facet_to_cut_cell[i]] += getindex!(c,array,cell_to_is_cut[i]) end end @@ -299,8 +302,9 @@ function add_facet_moments!(ccm::CutCellMoments, facet_to_bgcell = get_glue(trian,Val(Dp)).tface_to_mface facet_to_cut_cell = lazy_map(Reindex(ccm.bgcell_to_cut_cell),facet_to_bgcell) l = length(facet_to_cut_cell) + c = array_cache(array) for i = 1:l - ccm.data[facet_to_cut_cell[i]] += array[i] + ccm.data[facet_to_cut_cell[i]] += getindex!(c,array,i) end end From 1d2b3a7e9b7982b6a5b8ffe91aaa62674caa1362 Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Wed, 16 Aug 2023 15:25:21 +0200 Subject: [PATCH 18/19] deleted Manifest file --- Manifest.toml | 619 -------------------------------------------------- 1 file changed, 619 deletions(-) delete mode 100644 Manifest.toml diff --git a/Manifest.toml b/Manifest.toml deleted file mode 100644 index 881bf2d..0000000 --- a/Manifest.toml +++ /dev/null @@ -1,619 +0,0 @@ -# This file is machine-generated - editing it directly is not advised - -julia_version = "1.8.3" -manifest_format = "2.0" -project_hash = "3d6e007212983dcf6fe7381d5b73ff895d20ab4b" - -[[deps.AbstractFFTs]] -deps = ["ChainRulesCore", "LinearAlgebra"] -git-tree-sha1 = "16b6dbc4cf7caee4e1e75c49485ec67b667098a0" -uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c" -version = "1.3.1" - -[[deps.AbstractTrees]] -git-tree-sha1 = "03e0550477d86222521d254b741d470ba17ea0b5" -uuid = "1520ce14-60c1-5f80-bbc7-55ef81b5835c" -version = "0.3.4" - -[[deps.Adapt]] -deps = ["LinearAlgebra", "Requires"] -git-tree-sha1 = "cc37d689f599e8df4f464b2fa3870ff7db7492ef" -uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" -version = "3.6.1" - -[[deps.ArgCheck]] -git-tree-sha1 = "a3a402a35a2f7e0b87828ccabbd5ebfbebe356b4" -uuid = "dce04be8-c92d-5529-be00-80e4d2c0e197" -version = "2.3.0" - -[[deps.ArgTools]] -uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" -version = "1.1.1" - -[[deps.ArnoldiMethod]] -deps = ["LinearAlgebra", "Random", "StaticArrays"] -git-tree-sha1 = "f87e559f87a45bece9c9ed97458d3afe98b1ebb9" -uuid = "ec485272-7323-5ecc-a04f-4719b315124d" -version = "0.1.0" - -[[deps.ArrayInterface]] -deps = ["Adapt", "LinearAlgebra", "Requires", "SnoopPrecompile", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "53bb1691fb8560633ed8c0fa11d8b0900aaa805c" -uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" -version = "7.4.2" - -[[deps.ArrayLayouts]] -deps = ["FillArrays", "LinearAlgebra", "SparseArrays"] -git-tree-sha1 = "4aff5fa660eb95c2e0deb6bcdabe4d9a96bc4667" -uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" -version = "0.8.18" - -[[deps.Artifacts]] -uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33" - -[[deps.AutoHashEquals]] -git-tree-sha1 = "45bb6705d93be619b81451bb2006b7ee5d4e4453" -uuid = "15f4f7f2-30c1-5605-9d31-71845cf9641f" -version = "0.2.0" - -[[deps.BSON]] -git-tree-sha1 = "2208958832d6e1b59e49f53697483a84ca8d664e" -uuid = "fbb218c0-5317-5bc6-957e-2ee96dd4b1f0" -version = "0.3.7" - -[[deps.Base64]] -uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" - -[[deps.BlockArrays]] -deps = ["ArrayLayouts", "FillArrays", "LinearAlgebra"] -git-tree-sha1 = "3b15c61bcece7c426ea641d143c808ace3661973" -uuid = "8e7c35d0-a365-5155-bbbb-fb81a777f24e" -version = "0.16.25" - -[[deps.ChainRulesCore]] -deps = ["Compat", "LinearAlgebra", "SparseArrays"] -git-tree-sha1 = "c6d890a52d2c4d55d326439580c3b8d0875a77d9" -uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" -version = "1.15.7" - -[[deps.ChangesOfVariables]] -deps = ["ChainRulesCore", "LinearAlgebra", "Test"] -git-tree-sha1 = "485193efd2176b88e6622a39a246f8c5b600e74e" -uuid = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" -version = "0.1.6" - -[[deps.CodecZlib]] -deps = ["TranscodingStreams", "Zlib_jll"] -git-tree-sha1 = "9c209fb7536406834aa938fb149964b985de6c83" -uuid = "944b1d66-785c-5afd-91f1-9de20f533193" -version = "0.7.1" - -[[deps.Combinatorics]] -git-tree-sha1 = "08c8b6831dc00bfea825826be0bc8336fc369860" -uuid = "861a8166-3701-5b0c-9a16-15d98fcdc6aa" -version = "1.0.2" - -[[deps.CommonSubexpressions]] -deps = ["MacroTools", "Test"] -git-tree-sha1 = "7b8a93dba8af7e3b42fecabf646260105ac373f7" -uuid = "bbf7d656-a473-5ed7-a52c-81e309532950" -version = "0.3.0" - -[[deps.Compat]] -deps = ["Dates", "LinearAlgebra", "UUIDs"] -git-tree-sha1 = "7a60c856b9fa189eb34f5f8a6f6b5529b7942957" -uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" -version = "4.6.1" - -[[deps.CompilerSupportLibraries_jll]] -deps = ["Artifacts", "Libdl"] -uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" -version = "0.5.2+0" - -[[deps.ConstructionBase]] -deps = ["LinearAlgebra"] -git-tree-sha1 = "89a9db8d28102b094992472d333674bd1a83ce2a" -uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" -version = "1.5.1" - -[[deps.DataStructures]] -deps = ["Compat", "InteractiveUtils", "OrderedCollections"] -git-tree-sha1 = "d1fff3a548102f48987a52a2e0d114fa97d730f0" -uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8" -version = "0.18.13" - -[[deps.Dates]] -deps = ["Printf"] -uuid = "ade2ca70-3891-5945-98fb-dc099432e06a" - -[[deps.DiffResults]] -deps = ["StaticArraysCore"] -git-tree-sha1 = "782dd5f4561f5d267313f23853baaaa4c52ea621" -uuid = "163ba53b-c6d8-5494-b064-1a9d43ac40c5" -version = "1.1.0" - -[[deps.DiffRules]] -deps = ["IrrationalConstants", "LogExpFunctions", "NaNMath", "Random", "SpecialFunctions"] -git-tree-sha1 = "a4ad7ef19d2cdc2eff57abbbe68032b1cd0bd8f8" -uuid = "b552c78f-8df3-52c6-915a-8e097449b14b" -version = "1.13.0" - -[[deps.Distances]] -deps = ["LinearAlgebra", "SparseArrays", "Statistics", "StatsAPI"] -git-tree-sha1 = "49eba9ad9f7ead780bfb7ee319f962c811c6d3b2" -uuid = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7" -version = "0.10.8" - -[[deps.Distributed]] -deps = ["Random", "Serialization", "Sockets"] -uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" - -[[deps.DocStringExtensions]] -deps = ["LibGit2"] -git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d" -uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" -version = "0.9.3" - -[[deps.Downloads]] -deps = ["ArgTools", "FileWatching", "LibCURL", "NetworkOptions"] -uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6" -version = "1.6.0" - -[[deps.FFTW]] -deps = ["AbstractFFTs", "FFTW_jll", "LinearAlgebra", "MKL_jll", "Preferences", "Reexport"] -git-tree-sha1 = "f9818144ce7c8c41edf5c4c179c684d92aa4d9fe" -uuid = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341" -version = "1.6.0" - -[[deps.FFTW_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "c6033cc3892d0ef5bb9cd29b7f2f0331ea5184ea" -uuid = "f5851436-0d7a-5f13-b9de-f02708fd171a" -version = "3.3.10+0" - -[[deps.FastGaussQuadrature]] -deps = ["LinearAlgebra", "SpecialFunctions", "StaticArrays"] -git-tree-sha1 = "58d83dd5a78a36205bdfddb82b1bb67682e64487" -uuid = "442a2c76-b920-505d-bb47-c5924d526838" -version = "0.4.9" - -[[deps.FileIO]] -deps = ["Pkg", "Requires", "UUIDs"] -git-tree-sha1 = "7be5f99f7d15578798f338f5433b6c432ea8037b" -uuid = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549" -version = "1.16.0" - -[[deps.FileWatching]] -uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee" - -[[deps.FillArrays]] -deps = ["LinearAlgebra", "Random", "SparseArrays", "Statistics"] -git-tree-sha1 = "7072f1e3e5a8be51d525d64f63d3ec1287ff2790" -uuid = "1a297f60-69ca-5386-bcde-b61e274b549b" -version = "0.13.11" - -[[deps.FiniteDiff]] -deps = ["ArrayInterface", "LinearAlgebra", "Requires", "Setfield", "SparseArrays", "StaticArrays"] -git-tree-sha1 = "03fcb1c42ec905d15b305359603888ec3e65f886" -uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" -version = "2.19.0" - -[[deps.ForwardDiff]] -deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions", "StaticArrays"] -git-tree-sha1 = "00e252f4d706b3d55a8863432e742bf5717b498d" -uuid = "f6369f11-7733-5829-9624-2563aa707210" -version = "0.10.35" - -[[deps.Future]] -deps = ["Random"] -uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820" - -[[deps.Gridap]] -deps = ["AbstractTrees", "BSON", "BlockArrays", "Combinatorics", "DataStructures", "DocStringExtensions", "FastGaussQuadrature", "FileIO", "FillArrays", "ForwardDiff", "JLD2", "JSON", "LineSearches", "LinearAlgebra", "NLsolve", "NearestNeighbors", "PolynomialBases", "QuadGK", "Random", "SparseArrays", "SparseMatricesCSR", "StaticArrays", "Test", "WriteVTK"] -git-tree-sha1 = "9e0dd5beb345775ee9d4e9d579ff123999d3f381" -repo-rev = "master" -repo-url = "https://github.com/gridap/Gridap.jl.git" -uuid = "56d4f2e9-7ea1-5844-9cf6-b9c51ca7ce8e" -version = "0.17.17" - -[[deps.Inflate]] -git-tree-sha1 = "5cd07aab533df5170988219191dfad0519391428" -uuid = "d25df0c9-e2be-5dd7-82c8-3ad0b3e990b9" -version = "0.1.3" - -[[deps.IntelOpenMP_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "d979e54b71da82f3a65b62553da4fc3d18c9004c" -uuid = "1d5cc7b8-4909-519e-a0f8-d0f5ad9712d0" -version = "2018.0.3+2" - -[[deps.InteractiveUtils]] -deps = ["Markdown"] -uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" - -[[deps.InverseFunctions]] -deps = ["Test"] -git-tree-sha1 = "49510dfcb407e572524ba94aeae2fced1f3feb0f" -uuid = "3587e190-3f89-42d0-90ee-14403ec27112" -version = "0.1.8" - -[[deps.IrrationalConstants]] -git-tree-sha1 = "630b497eafcc20001bba38a4651b327dcfc491d2" -uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" -version = "0.2.2" - -[[deps.JLD2]] -deps = ["FileIO", "MacroTools", "Mmap", "OrderedCollections", "Pkg", "Printf", "Reexport", "Requires", "TranscodingStreams", "UUIDs"] -git-tree-sha1 = "42c17b18ced77ff0be65957a591d34f4ed57c631" -uuid = "033835bb-8acc-5ee8-8aae-3f567f8a3819" -version = "0.4.31" - -[[deps.JLLWrappers]] -deps = ["Preferences"] -git-tree-sha1 = "abc9885a7ca2052a736a600f7fa66209f96506e1" -uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210" -version = "1.4.1" - -[[deps.JSON]] -deps = ["Dates", "Mmap", "Parsers", "Unicode"] -git-tree-sha1 = "3c837543ddb02250ef42f4738347454f95079d4e" -uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6" -version = "0.21.3" - -[[deps.LazyArtifacts]] -deps = ["Artifacts", "Pkg"] -uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3" - -[[deps.LibCURL]] -deps = ["LibCURL_jll", "MozillaCACerts_jll"] -uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21" -version = "0.6.3" - -[[deps.LibCURL_jll]] -deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"] -uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0" -version = "7.84.0+0" - -[[deps.LibGit2]] -deps = ["Base64", "NetworkOptions", "Printf", "SHA"] -uuid = "76f85450-5226-5b5a-8eaa-529ad045b433" - -[[deps.LibSSH2_jll]] -deps = ["Artifacts", "Libdl", "MbedTLS_jll"] -uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8" -version = "1.10.2+0" - -[[deps.Libdl]] -uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb" - -[[deps.Libiconv_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "c7cb1f5d892775ba13767a87c7ada0b980ea0a71" -uuid = "94ce4f54-9a6c-5748-9c1c-f9c7231a4531" -version = "1.16.1+2" - -[[deps.LightGraphs]] -deps = ["ArnoldiMethod", "DataStructures", "Distributed", "Inflate", "LinearAlgebra", "Random", "SharedArrays", "SimpleTraits", "SparseArrays", "Statistics"] -git-tree-sha1 = "432428df5f360964040ed60418dd5601ecd240b6" -uuid = "093fc24a-ae57-5d10-9952-331d41423f4d" -version = "1.3.5" - -[[deps.LightXML]] -deps = ["Libdl", "XML2_jll"] -git-tree-sha1 = "e129d9391168c677cd4800f5c0abb1ed8cb3794f" -uuid = "9c8b4983-aa76-5018-a973-4c85ecc9e179" -version = "0.9.0" - -[[deps.LineSearches]] -deps = ["LinearAlgebra", "NLSolversBase", "NaNMath", "Parameters", "Printf"] -git-tree-sha1 = "7bbea35cec17305fc70a0e5b4641477dc0789d9d" -uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255" -version = "7.2.0" - -[[deps.LinearAlgebra]] -deps = ["Libdl", "libblastrampoline_jll"] -uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" - -[[deps.LogExpFunctions]] -deps = ["ChainRulesCore", "ChangesOfVariables", "DocStringExtensions", "InverseFunctions", "IrrationalConstants", "LinearAlgebra"] -git-tree-sha1 = "0a1b7c2863e44523180fdb3146534e265a91870b" -uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" -version = "0.3.23" - -[[deps.Logging]] -uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" - -[[deps.MKL_jll]] -deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "Pkg"] -git-tree-sha1 = "2ce8695e1e699b68702c03402672a69f54b8aca9" -uuid = "856f044c-d86e-5d09-b602-aeab76dc8ba7" -version = "2022.2.0+0" - -[[deps.MacroTools]] -deps = ["Markdown", "Random"] -git-tree-sha1 = "42324d08725e200c23d4dfb549e0d5d89dede2d2" -uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" -version = "0.5.10" - -[[deps.Markdown]] -deps = ["Base64"] -uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" - -[[deps.MbedTLS_jll]] -deps = ["Artifacts", "Libdl"] -uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" -version = "2.28.0+0" - -[[deps.MiniQhull]] -deps = ["QhullMiniWrapper_jll"] -git-tree-sha1 = "9dc837d180ee49eeb7c8b77bb1c860452634b0d1" -uuid = "978d7f02-9e05-4691-894f-ae31a51d76ca" -version = "0.4.0" - -[[deps.Mmap]] -uuid = "a63ad114-7e13-5084-954f-fe012c677804" - -[[deps.MozillaCACerts_jll]] -uuid = "14a3606d-f60d-562e-9121-12d972cd8159" -version = "2022.2.1" - -[[deps.NLSolversBase]] -deps = ["DiffResults", "Distributed", "FiniteDiff", "ForwardDiff"] -git-tree-sha1 = "a0b464d183da839699f4c79e7606d9d186ec172c" -uuid = "d41bc354-129a-5804-8e4c-c37616107c6c" -version = "7.8.3" - -[[deps.NLsolve]] -deps = ["Distances", "LineSearches", "LinearAlgebra", "NLSolversBase", "Printf", "Reexport"] -git-tree-sha1 = "019f12e9a1a7880459d0173c182e6a99365d7ac1" -uuid = "2774e3e8-f4cf-5e23-947b-6d7e65073b56" -version = "4.5.1" - -[[deps.NaNMath]] -deps = ["OpenLibm_jll"] -git-tree-sha1 = "0877504529a3e5c3343c6f8b4c0381e57e4387e4" -uuid = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" -version = "1.0.2" - -[[deps.NearestNeighbors]] -deps = ["Distances", "StaticArrays"] -git-tree-sha1 = "2c3726ceb3388917602169bed973dbc97f1b51a8" -uuid = "b8a86587-4115-5ab1-83bc-aa920d37bbce" -version = "0.4.13" - -[[deps.NetworkOptions]] -uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908" -version = "1.2.0" - -[[deps.OpenBLAS_jll]] -deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] -uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" -version = "0.3.20+0" - -[[deps.OpenLibm_jll]] -deps = ["Artifacts", "Libdl"] -uuid = "05823500-19ac-5b8b-9628-191a04bc5112" -version = "0.8.1+0" - -[[deps.OpenSpecFun_jll]] -deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1" -uuid = "efe28fd5-8261-553b-a9e1-b2916fc3738e" -version = "0.5.5+0" - -[[deps.OrderedCollections]] -git-tree-sha1 = "d321bf2de576bf25ec4d3e4360faca399afca282" -uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" -version = "1.6.0" - -[[deps.Parameters]] -deps = ["OrderedCollections", "UnPack"] -git-tree-sha1 = "34c0e9ad262e5f7fc75b10a9952ca7692cfc5fbe" -uuid = "d96e819e-fc66-5662-9728-84c9c7592b0a" -version = "0.12.3" - -[[deps.Parsers]] -deps = ["Dates", "SnoopPrecompile"] -git-tree-sha1 = "478ac6c952fddd4399e71d4779797c538d0ff2bf" -uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0" -version = "2.5.8" - -[[deps.Pkg]] -deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] -uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" -version = "1.8.0" - -[[deps.PolynomialBases]] -deps = ["ArgCheck", "AutoHashEquals", "FFTW", "FastGaussQuadrature", "LinearAlgebra", "Requires", "SimpleUnPack", "SpecialFunctions"] -git-tree-sha1 = "cfbe38cc661a5a418e198cb21cdc42258a586301" -uuid = "c74db56a-226d-5e98-8bb0-a6049094aeea" -version = "0.4.17" - -[[deps.Preferences]] -deps = ["TOML"] -git-tree-sha1 = "47e5f437cc0e7ef2ce8406ce1e7e24d44915f88d" -uuid = "21216c6a-2e73-6563-6e65-726566657250" -version = "1.3.0" - -[[deps.Printf]] -deps = ["Unicode"] -uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7" - -[[deps.QhullMiniWrapper_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg", "Qhull_jll"] -git-tree-sha1 = "607cf73c03f8a9f83b36db0b86a3a9c14179621f" -uuid = "460c41e3-6112-5d7f-b78c-b6823adb3f2d" -version = "1.0.0+1" - -[[deps.Qhull_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] -git-tree-sha1 = "238dd7e2cc577281976b9681702174850f8d4cbc" -uuid = "784f63db-0788-585a-bace-daefebcd302b" -version = "8.0.1001+0" - -[[deps.QuadGK]] -deps = ["DataStructures", "LinearAlgebra"] -git-tree-sha1 = "6ec7ac8412e83d57e313393220879ede1740f9ee" -uuid = "1fd47b50-473d-5c70-9696-f719f8f3bcdc" -version = "2.8.2" - -[[deps.REPL]] -deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] -uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" - -[[deps.Random]] -deps = ["SHA", "Serialization"] -uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" - -[[deps.Reexport]] -git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b" -uuid = "189a3867-3050-52da-a836-e630ba90ab69" -version = "1.2.2" - -[[deps.Requires]] -deps = ["UUIDs"] -git-tree-sha1 = "838a3a4188e2ded87a4f9f184b4b0d78a1e91cb7" -uuid = "ae029012-a4dd-5104-9daa-d747884805df" -version = "1.3.0" - -[[deps.SHA]] -uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce" -version = "0.7.0" - -[[deps.Serialization]] -uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b" - -[[deps.Setfield]] -deps = ["ConstructionBase", "Future", "MacroTools", "StaticArraysCore"] -git-tree-sha1 = "e2cc6d8c88613c05e1defb55170bf5ff211fbeac" -uuid = "efcf1570-3423-57d1-acb7-fd33fddbac46" -version = "1.1.1" - -[[deps.SharedArrays]] -deps = ["Distributed", "Mmap", "Random", "Serialization"] -uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383" - -[[deps.SimpleTraits]] -deps = ["InteractiveUtils", "MacroTools"] -git-tree-sha1 = "5d7e3f4e11935503d3ecaf7186eac40602e7d231" -uuid = "699a6c99-e7fa-54fc-8d76-47d257e15c1d" -version = "0.9.4" - -[[deps.SimpleUnPack]] -git-tree-sha1 = "58e6353e72cde29b90a69527e56df1b5c3d8c437" -uuid = "ce78b400-467f-4804-87d8-8f486da07d0a" -version = "1.1.0" - -[[deps.SnoopPrecompile]] -deps = ["Preferences"] -git-tree-sha1 = "e760a70afdcd461cf01a575947738d359234665c" -uuid = "66db9d55-30c0-4569-8b51-7e840670fc0c" -version = "1.0.3" - -[[deps.Sockets]] -uuid = "6462fe0b-24de-5631-8697-dd941f90decc" - -[[deps.SparseArrays]] -deps = ["LinearAlgebra", "Random"] -uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" - -[[deps.SparseMatricesCSR]] -deps = ["LinearAlgebra", "SparseArrays", "SuiteSparse"] -git-tree-sha1 = "38677ca58e80b5cad2382e5a1848f93b054ad28d" -uuid = "a0a7dd2c-ebf4-11e9-1f05-cf50bc540ca1" -version = "0.6.7" - -[[deps.SpecialFunctions]] -deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] -git-tree-sha1 = "ef28127915f4229c971eb43f3fc075dd3fe91880" -uuid = "276daf66-3868-5448-9aa4-cd146d93841b" -version = "2.2.0" - -[[deps.StaticArrays]] -deps = ["LinearAlgebra", "Random", "StaticArraysCore", "Statistics"] -git-tree-sha1 = "b8d897fe7fa688e93aef573711cb207c08c9e11e" -uuid = "90137ffa-7385-5640-81b9-e52037218182" -version = "1.5.19" - -[[deps.StaticArraysCore]] -git-tree-sha1 = "6b7ba252635a5eff6a0b0664a41ee140a1c9e72a" -uuid = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" -version = "1.4.0" - -[[deps.Statistics]] -deps = ["LinearAlgebra", "SparseArrays"] -uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" - -[[deps.StatsAPI]] -deps = ["LinearAlgebra"] -git-tree-sha1 = "45a7769a04a3cf80da1c1c7c60caf932e6f4c9f7" -uuid = "82ae8749-77ed-4fe6-ae5f-f523153014b0" -version = "1.6.0" - -[[deps.SuiteSparse]] -deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"] -uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" - -[[deps.TOML]] -deps = ["Dates"] -uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" -version = "1.0.0" - -[[deps.Tar]] -deps = ["ArgTools", "SHA"] -uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" -version = "1.10.1" - -[[deps.Test]] -deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] -uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40" - -[[deps.TranscodingStreams]] -deps = ["Random", "Test"] -git-tree-sha1 = "94f38103c984f89cf77c402f2a68dbd870f8165f" -uuid = "3bb67fe8-82b1-5028-8e26-92a6c54297fa" -version = "0.9.11" - -[[deps.UUIDs]] -deps = ["Random", "SHA"] -uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" - -[[deps.UnPack]] -git-tree-sha1 = "387c1f73762231e86e0c9c5443ce3b4a0a9a0c2b" -uuid = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" -version = "1.0.2" - -[[deps.Unicode]] -uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5" - -[[deps.WriteVTK]] -deps = ["Base64", "CodecZlib", "FillArrays", "LightXML", "TranscodingStreams"] -git-tree-sha1 = "49353f30da65f377cff0f934bb9f562a2c0441b9" -uuid = "64499a7a-5c06-52f2-abe2-ccb03c286192" -version = "1.17.1" - -[[deps.XML2_jll]] -deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "Zlib_jll"] -git-tree-sha1 = "93c41695bc1c08c46c5899f4fe06d6ead504bb73" -uuid = "02c8fc9c-b97f-50b9-bbe4-9be30ff0a78a" -version = "2.10.3+0" - -[[deps.Zlib_jll]] -deps = ["Libdl"] -uuid = "83775a58-1f1d-513f-b197-d71354ab007a" -version = "1.2.12+3" - -[[deps.libblastrampoline_jll]] -deps = ["Artifacts", "Libdl", "OpenBLAS_jll"] -uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" -version = "5.1.1+0" - -[[deps.nghttp2_jll]] -deps = ["Artifacts", "Libdl"] -uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" -version = "1.48.0+0" - -[[deps.p7zip_jll]] -deps = ["Artifacts", "Libdl"] -uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0" -version = "17.4.0+0" From 2a08ea5084aa428fb421b9726fb03a6d5b9bbf3e Mon Sep 17 00:00:00 2001 From: Eric Neiva Date: Mon, 21 Aug 2023 12:53:06 +0200 Subject: [PATCH 19/19] Update NEWS.md --- NEWS.md | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/NEWS.md b/NEWS.md index 6314fbe..ca8b257 100644 --- a/NEWS.md +++ b/NEWS.md @@ -4,6 +4,11 @@ All notable changes to this project will be documented in this file. The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/), and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). +## [0.8.2] - 2023-08-22 + +### Added +- Moment fitted machinery. Since PR [#68](https://github.com/gridap/GridapEmbedded.jl/pull/68). + ## [0.8.0] - 2021-11-03 ### Changed