Skip to content

Commit

Permalink
init qp
Browse files Browse the repository at this point in the history
  • Loading branch information
Abdelrahman912 committed Feb 11, 2025
1 parent d50f7fa commit 32c5763
Show file tree
Hide file tree
Showing 7 changed files with 60 additions and 37 deletions.
26 changes: 21 additions & 5 deletions src/ferrite-addons/PR883.jl
Original file line number Diff line number Diff line change
@@ -1,3 +1,16 @@

"""
QuadraturePoint{dim, T}
A simple helper to carry quadrature point information.
"""
struct QuadraturePoint{dim, T,VEC}
i::Int
ξ::VEC
end

QuadraturePoint(i::Int, ξ::VEC) where {dim, T,VEC} = QuadraturePoint{dim, T,VEC}(i, ξ)

using Ferrite
using StaticArrays
import Base: @propagate_inbounds
Expand Down Expand Up @@ -152,6 +165,7 @@ struct StaticQuadratureValues{T, N_t, dNdx_t, M_t, NumN, NumM} <: AbstractQuadra
N::SVector{NumN, N_t}
dNdx::SVector{NumN, dNdx_t}
M::SVector{NumM, M_t}
qp::QuadraturePoint
end

@propagate_inbounds Ferrite.getngeobasefunctions(qv::StaticQuadratureValues) = length(qv.M)
Expand Down Expand Up @@ -200,6 +214,8 @@ end

Ferrite.shape_value(siv::StaticInterpolationValues, qp::Int, i::Int) = siv.Nξ[i, qp]
Ferrite.getnbasefunctions(siv::StaticInterpolationValues) = getnbasefunctions(siv.ip)
Ferrite.shape_value(cv::FerriteUtils.StaticInterpolationValues, qp::QuadraturePoint, base_fun_idx::Int) = Ferrite.shape_value(cv, qp.i, base_fun_idx)
Ferrite.function_value(cv::FerriteUtils.StaticInterpolationValues, qp::QuadraturePoint, ue) = Ferrite.function_value(cv, qp.i, ue)

# Dispatch on DiffOrder parameter?
# Reuse functions for GeometryMapping - same signature but need access functions
Expand Down Expand Up @@ -227,17 +243,16 @@ end
struct StaticCellValues{FV, GM, Nqp, T}
fv::FV # StaticInterpolationValues
gm::GM # StaticInterpolationValues
#x::Tx # AbstractVector{<:Vec} or Nothing
weights::NTuple{Nqp, T}
qps::NTuple{Nqp,QuadraturePoint} # quadrature points
end

function StaticCellValues(cv::CellValues)
fv = StaticInterpolationValues(cv.fun_values)
gm = StaticInterpolationValues(cv.geo_mapping)
#sdim = sdim_from_gradtype(shape_gradient_type(cv))
#x = SaveCoords ? fill(zero(Vec{sdim}), getngeobasefunctions(cv)) : nothing
weights = ntuple(i -> getweights(cv.qr)[i], getnquadpoints(cv))
return StaticCellValues(fv, gm, weights)
qps = ntuple(i -> QuadraturePoint(i, getpoints(cv.qr)[i]), getnquadpoints(qr))
return StaticCellValues(fv, gm,weights, qps)
end

# function StaticCellValues(cv::CellValues, ::Val{SaveCoords}=Val(true)) where SaveCoords
Expand Down Expand Up @@ -285,6 +300,7 @@ function _quadrature_point_values(fe_v::StaticCellValues, q_point::Int, cell_coo

Nx, dNdx = calculate_mapped_values(fe_v.fv, q_point, mapping)
M = fe_v.gm.Nξ[:, q_point]
qp = fe_v.qps[q_point]
end
return StaticQuadratureValues(detJdV, Nx, dNdx, M)
return StaticQuadratureValues(detJdV, Nx, dNdx, M,qp)
end
19 changes: 18 additions & 1 deletion src/ferrite-addons/gpu/adapt.jl
Original file line number Diff line number Diff line change
@@ -1,18 +1,35 @@
############
# adapt.jl #
############
import Tensors: Vec
Adapt.@adapt_structure QuadratureValuesIterator
Adapt.@adapt_structure StaticQuadratureValues
Adapt.@adapt_structure DeviceCellIterator
Adapt.@adapt_structure DeviceGrid
Adapt.@adapt_structure DeviceDofHandler
Adapt.@adapt_structure DeviceDofHandlerData
Adapt.@adapt_structure DeviceSubDofHandlerData
#Adapt.@adapt_structure QuadraturePoint
Adapt.@adapt_structure Vec

function Adapt.adapt_structure(to, cv::CellValues)
fv = Adapt.adapt(to, StaticInterpolationValues(cv.fun_values))
gm = Adapt.adapt(to, StaticInterpolationValues(cv.geo_mapping))
n_quadoints = cv.qr.weights |> length
weights = Adapt.adapt(to, ntuple(i -> cv.qr.weights[i], n_quadoints))
return StaticCellValues(fv, gm, weights)
#qps = Adapt.adapt(to,ntuple(i -> Adapt.adapt_structure(to,QuadraturePoint(i, cv.qr.points[i])), n_quadoints))
qps = Adapt.adapt(to,ntuple(i -> QuadraturePoint(i, Adapt.adapt_structure(to,Vec{1}((0.0f0,)))), n_quadoints))
return StaticCellValues(fv, gm, weights,qps)
end

function Adapt.adapt_structure(to,qp::QuadraturePoint)
@show "A7a"
i = Adapt.adapt(to, qp.i)
ξ = Adapt.adapt_structure(to, qp.ξ)
return QuadraturePoint(i, ξ)
end

function Adapt.adapt_structure(to,qp::Vec{dim,T}) where {dim,T}
@show "A7a2"
return Vec{dim,T}(Adapt.adapt(to, qp.data))
end
5 changes: 1 addition & 4 deletions src/ferrite-addons/gpu/device_memalloc.jl
Original file line number Diff line number Diff line change
Expand Up @@ -29,12 +29,9 @@ function try_allocate_shared_mem(::Type{AbstractMemShape{Tv}}, block_dim::Ti, n_
error("please provide concrete implementation for $(typeof(AbstractMemShape{Tv})).")
end



cellmem(shared_mem::AbstractDeviceSharedMem, ::Integer) = error("please provide concrete implementation for $(typeof(shared_mem)).")
cellmem(shared_mem::AbstractDeviceGlobalMem, ::Integer) = error("please provide concrete implementation for $(typeof(shared_mem)).")


function allocate_global_mem(::Type{AbstractMemShape{Tv}}, n_cells::Ti, n_basefuncs::Ti) where {Ti <: Integer, Tv <: Real}
error("please provide concrete implementation for $(typeof(AbstractMemShape{Tv})).")
end
Expand Down Expand Up @@ -62,4 +59,4 @@ struct KeFeCellMem{MatrixType, VectorType} <: AbstractCellMem
fe::VectorType
end

struct NoCellMem <: AbstractCellMem end # mainly for testing purposes
struct NoCellMem <: AbstractCellMem end # mainly for testing purposes
13 changes: 0 additions & 13 deletions src/ferrite-addons/quadrature_iterator.jl
Original file line number Diff line number Diff line change
@@ -1,13 +1,3 @@
"""
QuadraturePoint{dim, T}
A simple helper to carry quadrature point information.
"""
struct QuadraturePoint{dim, T}
i::Int
ξ::Vec{dim, T}
end

"""
QuadratureIterator(::QuadratureRule)
QuadratureIterator(::FacetQuadratureRule, local_face_idx::Int)
Expand Down Expand Up @@ -44,6 +34,3 @@ Ferrite.shape_gradient(fv::FacetValues, qp::QuadraturePoint, base_fun_idx::Int)
Ferrite.function_value(fv::FacetValues, qp::QuadraturePoint, ue) = Ferrite.function_value(fv, qp.i, ue)
Ferrite.function_gradient(fv::FacetValues, qp::QuadraturePoint, ue) = Ferrite.function_gradient(fv, qp.i, ue)
Ferrite.getnormal(fv::FacetValues, qp::QuadraturePoint) = Ferrite.getnormal(fv, qp.i)

Ferrite.shape_value(cv::FerriteUtils.StaticInterpolationValues, qp::QuadraturePoint, base_fun_idx::Int) = Ferrite.shape_value(cv, qp.i, base_fun_idx)
Ferrite.function_value(cv::FerriteUtils.StaticInterpolationValues, qp::QuadraturePoint, ue) = Ferrite.function_value(cv, qp.i, ue)
25 changes: 14 additions & 11 deletions src/modeling/core/coefficients.jl
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
Cell_Cache = Union{CellCache, FerriteUtils.DeviceCellCache}

"""
FieldCoefficient(data, interpolation)
Expand Down Expand Up @@ -36,7 +38,7 @@ function _create_field_coefficient_cache(coefficient::FieldCoefficient{<:Vec{<:A
return FieldCoefficientCache(coefficient.elementwise_data, FerriteUtils.StaticInterpolationValues(fv.ip, SMatrix{Nξs[1], Nξs[2]}(fv.Nξ), nothing))
end

function evaluate_coefficient(cache::FieldCoefficientCache{T}, geometry_cache, qp::QuadraturePoint, t) where T
function evaluate_coefficient(cache::FieldCoefficientCache{T}, geometry_cache::Cell_Cache, qp::QuadraturePoint, t) where T
@unpack elementwise_data, cv = cache
val = zero(T)
cellidx = cellid(geometry_cache)
Expand All @@ -46,17 +48,18 @@ function evaluate_coefficient(cache::FieldCoefficientCache{T}, geometry_cache, q
end
return val
end
# GPU coefficient evaluation!
function evaluate_coefficient(cache::FieldCoefficientCache{T}, geometry_cache::FerriteUtils.DeviceCellCache, qv::FerriteUtils.StaticQuadratureValues, t) where T
@unpack elementwise_data, cv = cache
val = zero(T)
cellidx = Ferrite.cellid(geometry_cache)

@inbounds for i in 1:getnbasefunctions(cv)
val += Ferrite.shape_value(qv, i) * elementwise_data[i, cellidx]
end
return val
end
# # GPU coefficient evaluation!
# function evaluate_coefficient(cache::FieldCoefficientCache{T}, geometry_cache::FerriteUtils.DeviceCellCache, qv::FerriteUtils.StaticQuadratureValues, t) where T
# @unpack elementwise_data, cv = cache
# val = zero(T)
# cellidx = Ferrite.cellid(geometry_cache)

# @inbounds for i in 1:getnbasefunctions(cv)
# val += Ferrite.shape_value(qv, i) * elementwise_data[i, cellidx]
# end
# return val
# end


"""
Expand Down
4 changes: 3 additions & 1 deletion src/utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -50,10 +50,12 @@ include("ferrite-addons/gpu/adapt.jl")

end

import .FerriteUtils: QuadraturePoint

include("ferrite-addons/collections.jl")
include("ferrite-addons/quadrature_iterator.jl")

function celldofsview(dh::Ferrite.AbstractDofHandler, i::Ti) where {Ti <: Integer}
function celldofsview(dh::Ferrite.AbstractDofHandler, i::Integer)
ndofs = ndofs_per_cell(dh, i)
offset = dh.cell_dofs_offset[i]
return @views dh.cell_dofs[offset:(offset+ndofs-1)]
Expand Down
5 changes: 3 additions & 2 deletions test/gpu/test_coefficients.jl
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,8 @@ function coeffs_kernel!(Vals, dh, coeff_cache, cv, t)
quadrature_iterator = QuadratureValuesIterator(cv, coords)
n_quad_points = length(quadrature_iterator) # number of quadrature points
for (i, qv) in pairs(quadrature_iterator) # pairs will fetch the current index (i) and the `StaticQuadratureValue` (qv)
fx = evaluate_coefficient(coeff_cache, cell, qv, t)
qp = qv.qp
fx = evaluate_coefficient(coeff_cache, cell, qp, t)
Vals[(cell_id-1)*n_quad_points+i] = fx
end
end
Expand Down Expand Up @@ -58,8 +59,8 @@ import Adapt: adapt_structure
n_cells = grid.cells |> length
sdh = first(dh.subdofhandlers)

cu_dh = adapt_structure(cuda_strategy, dh)
strategy = Thunderbolt.CudaAssemblyStrategy(Float32, Int32)
cu_dh = adapt_structure(strategy, dh)
device_dh = deep_adapt(strategy, cu_dh.gpudata)
GC.@preserve cu_dh begin
@testset "ConstantCoefficient($val" for val [1.0f0]
Expand Down

0 comments on commit 32c5763

Please sign in to comment.