# Finite Element Assembly using map over elements

**URL:** <https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576>\
**Category:** Numerics\
**Created:** [April 4, 2021, 8:56pm UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576 "2021-04-04T20:56:00Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![ziolai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ziolai/32/23422_2.png) [@ziolai](https://discourse.julialang.org/u/ziolai)\
**Post date:** [April 4, 2021, 8:56pm UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/1 "2021-04-04T20:56:00Z")

</div>

Hi,

I am wondering whether the loop over the elements in a finite element assembly can be replace by a flavor of map over the list of elements. I am happy to consider the linear elements on a one-dimensional mesh to start from.

Assume thus having on the interval (or computational domain) 0 \<= x \<= 1 three elements (or subintervals) 0 \<= x \<= 1/3, 1/3 \<= x \<= 2/3 and 2/3 \<= x \<= 1. We thus have as connectivity matrix

e = [1 2; 2 3; 3 4]

and a list of coordinates p = [0 1/3 2/3 1].

Is there a way to compute the length (area) of each subinterval by a map on the matrix e?

Thanks. Domenico.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [April 4, 2021, 8:59pm UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/2 "2021-04-04T20:59:42Z")

</div>

Have you seen [Gridap.jl](https://github.com/gridap/Gridap.jl)? There’s no need to implement FEM from scratch.

---

<div class="post-metadata">

**Author:** ![PetrKryslUCSD](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petrkryslucsd/32/215825_2.png) [@PetrKryslUCSD](https://discourse.julialang.org/u/PetrKryslUCSD)\
**Post date:** [April 4, 2021, 11:53pm UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/3 "2021-04-04T23:53:17Z")

</div>

It can be done. There may not be any advantage to re-doing something that works well.

---

<div class="post-metadata">

**Author:** ![ziolai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ziolai/32/23422_2.png) [@ziolai](https://discourse.julialang.org/u/ziolai)\
**Post date:** [April 5, 2021, 7:33am UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/4 "2021-04-05T07:33:34Z")

</div>

Dear Petr and Steven,

Sincere thanks to both of you for getting in touch.

Yes, I have seen GridApp.jl (and similar such as JuliaFEM.jl and JuAFEM.jl) and the overview given here [JuliaPDE · GitHub](https://github.com/juliapde) . And yes, I do agree that there might be no intrinsic need to replace the loop over elements by a map.

My interest, however, is in gaining a better understanding of Julia in general and of coding PDE solvers in Julia in particular. Replacing the for-loop by a map might constitute a fun coding exercise. Being able to go beyond a blind use of a package like GridApp.jl might put me in the position to contribute to the package in a next stage.

I therefore remains curious in the question as originally posted.

Greetings, Domenico.

---

<div class="post-metadata">

**Author:** ![PetrKryslUCSD](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petrkryslucsd/32/215825_2.png) [@PetrKryslUCSD](https://discourse.julialang.org/u/PetrKryslUCSD)\
**Post date:** [April 5, 2021, 8:25am UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/5 "2021-04-05T08:25:32Z")

</div>

Fair enough. Here is to get you started:

```julia
julia> map(i -> p[e[i, 2]]-p[e[i, 1]], 1:size(e, 1))
3-element Vector{Float64}:
 3.33333e-01
 3.33333e-01
 3.33333e-01

```

---

<div class="post-metadata">

**Author:** ![ziolai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ziolai/32/23422_2.png) [@ziolai](https://discourse.julialang.org/u/ziolai)\
**Post date:** [April 5, 2021, 8:50am UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/6 "2021-04-05T08:50:16Z")

</div>

Wonderful, very nice!

This, however, assumes that the variables p and e are global variables to the anonymous function in the input variables i. Is there a way to avoid this? I tried a map over both e and p, leading to the error that e and p are of different size.

Can mapreduce be exploited to assemble the stiffness matrix, that is, can reduce be instructed to take node-element connectivity into account?

Thanks again, Domenico.

---

<div class="post-metadata">

**Author:** ![ziolai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ziolai/32/23422_2.png) [@ziolai](https://discourse.julialang.org/u/ziolai)\
**Post date:** [April 5, 2021, 9:24am UTC](https://discourse.julialang.org/t/finite-element-assembly-using-map-over-elements/58576/7 "2021-04-05T09:24:19Z")

</div>

@petr: your advice has been instrumental. Thx!

Below how to remove dependence on global variables.

Below can be extended to elementary matrix contribution. The set of matrix contributions can subsequently be assembled into the global matrix using the sparse function.

Can mapreduce be used for the global assembly?

How do map and explicit for loop compare in assembling matrices?

Does map allow easy extension to pmap for parellel assembly?

```julia
function compute_element_area(elem_id, e, p)
  area_id = p[e[elem_id,2]] - p[e[elem_id,1]]
  return area_id
end 

elem_ids = [1:3] 
e = [1 2; 2 3; 3 4]
p = [0 1/3 2/3 1]

map(elem_ids) do elem_id 
    compute_element_area(elem_id, e, p)
end

```