September 4, 2025

What is Gridap

Gridap is a finite element library written entirely in Julia. It allows you to describe partial differential equations in a high-level, mathematical way: you define meshes, trial and test spaces, weak forms, and material laws almost as you would write them on paper. The library then takes care of assembling and solving the discrete system. Because it is pure Julia, Gridap benefits from the language’s speed, composability, and ease of use; you can directly integrate it with Julia’s numerical and scientific computing ecosystem. One of its main advantages is flexibility: you can extend it to handle different kinds of PDEs, from simple Poisson problems to elasticity and multiphysics, without having to rely on external C++ backends. This makes it attractive for researchers and developers who want both performance and a concise, expressive way of implementing FEM.

Tutorial Code


using Gridap
using GridapGmsh

model = DiscreteModelFromFile("./models/solid.json")
path = "./results/script/"
files_path = path*"data/"
mkpath(files_path)

writevtk(model,"model")

order = 1

reffe = ReferenceFE(lagrangian,VectorValue{3,Float64},order)
V0 = TestFESpace(model,reffe;
  conformity=:H1,
  dirichlet_tags=["surface_1","surface_2"],
  dirichlet_masks=[(true,false,false), (true,true,true)])


g1(x) = VectorValue(0.005,0.0,0.0)
g2(x) = VectorValue(0.0,0.0,0.0)


U = TrialFESpace(V0,[g1,g2])


const E = 70.0e9
const ν = 0.33
const λ = (E*ν)/((1+ν)*(1-2*ν))
const μ = E/(2*(1+ν))
σ(ε) = λ*tr(ε)*one(ε) + 2*μ*ε


degree = 2*order
Ω = Triangulation(model)
dΩ = Measure(Ω,degree)


a(u,v) = ∫( ε(v) ⊙ (σ∘ε(u)) )*dΩ
l(v) = 0

op = AffineFEOperator(a,l,U,V0)
uh = solve(op)


writevtk(Ω,files_path*"Omega",cellfields=["uh"=>uh,"epsi"=>ε(uh),"sigma"=>σ∘ε(uh)])

  \lambda = \frac{E \, \nu}{(1+\nu)(1-2\nu)}, 
  \quad 
  \mu = \frac{E}{2(1+\nu)}


Hooke’s\space law (isotropic\space linear\space elasticity)\\
  \sigma(\varepsilon) = \lambda \, \mathrm{tr}(\varepsilon) \, I \;+\; 2 \mu \, \varepsilon
\text{Small strain tensor (symmetric gradient)}\\
  \varepsilon(u) = \tfrac{1}{2} \left( \nabla u + \nabla u^\top \right)

Result