# Command for third order derivative

**URL:** <https://community.freefem.org/t/command-for-third-order-derivative/4057>\
**Category:** General Discussion\
**Created:** [September 6, 2025, 7:50am UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057 "2025-09-06T07:50:23Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![debendra](https://avatars.discourse-cdn.com/v4/letter/d/8baadc/32.png) [@debendra](https://community.freefem.org/u/debendra)\
**Post date:** [September 6, 2025, 7:50am UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/1 "2025-09-06T07:50:23Z")

</div>

// Third derivative example in FreeFem++

load “msh3” // only if you want 3D later, otherwise optional

// Define mesh  
mesh Th = square(10, 10);

// Define finite element space (need degree \>= 3 for 3rd derivatives)  
fespace Vh(Th, P3);

// Define a test function  
Vh u = x^3 + 2_x^2_y + y^3;

// Compute third derivatives  
Vh d3x = dx(dxx(u)); // ∂³u / ∂x³  
Vh d3y = dy(dyy(u)); // ∂³u / ∂y³  
Vh d2x1y = dy(dxx(u)); // ∂³u / ∂x²∂y  
Vh d1x2y = dx(dyy(u)); // ∂³u / ∂x∂y²

// Print values at a point (0.5, 0.5)  
real px = 0.5, py = 0.5;  
cout \<\< "u(x,y) = " \<\< u(px,py) \<\< endl;  
cout \<\< "d³u/dx³ = " \<\< d3x(px,py) \<\< endl;  
cout \<\< "d³u/dy³ = " \<\< d3y(px,py) \<\< endl;  
cout \<\< "d³u/dx²dy = " \<\< d2x1y(px,py) \<\< endl;  
cout \<\< "d³u/dxdy² = " \<\< d1x2y(px,py) \<\< endl;

// Plot results  
plot(u, value=1, cmm=“u(x,y)”);  
plot(d3x, value=1, cmm=“∂³u/∂x³”);  
plot(d3y, value=1, cmm=“∂³u/∂y³”);  
plot(d2x1y, value=1, cmm=“∂³u/∂x²∂y”);  
plot(d1x2y, value=1, cmm=“∂³u/∂x∂y²”);

What is the command for third order derivative?

This is not working

---

<div class="post-metadata">

**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\
**Post date:** [September 6, 2025, 11:31am UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/2 "2025-09-06T11:31:46Z")

</div>

A command for 3rd order derivative would be useless. Already, using second order derivatives commands leads to mistakes. This is because most of the time, functions in finite element spaces are continuous, but their derivatives are discontinuous. It follows that computing a second derivative inside the element misses the jump part of the first-order derivative.

Instead of that, if you want to compute second or third derivatives, you have to interpolate discontinuous functions into spaces of continuous functions before applying derivative operators dx(), dy(). This gives

[der3.edp](https://community.freefem.org/uploads/short-url/d8wQZONbtyw2dccBOG8MQNtq1xL.edp) (1.5 KB)

---

<div class="post-metadata">

**Author:** ![debendra](https://avatars.discourse-cdn.com/v4/letter/d/8baadc/32.png) [@debendra](https://community.freefem.org/u/debendra)\
**Post date:** [December 27, 2025, 5:18pm UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/3 "2025-12-27T17:18:12Z")

</div>

Thank you!!, Will it work in the weak formulation?

---

<div class="post-metadata">

**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\
**Post date:** [December 29, 2025, 6:54pm UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/4 "2025-12-29T18:54:47Z")

</div>

Yes, in a weak formulation the interpolated derivatives in a continuous space have to be understood as independent variables.  
For example if you want to solve  
-\Delta u=f  
you can write the system  
\partial\_x u=u\_x,  
\partial\_y u=u\_y,  
-\partial\_x u\_x-\partial\_y u\_y=f  
with unknowns u,u\_x,u\_y in P1. Then you write a variational formulation with 3 unknowns and 3 test functions.  
This is the mixed finite element method.

---

<div class="post-metadata">

**Author:** ![debenswain](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/debenswain/32/3800_2.png) [@debenswain](https://community.freefem.org/u/debenswain)\
**Post date:** [August 21, 2026, 1:18pm UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/5 "2026-08-21T13:18:25Z")

</div>

Dear [François Bouchut](https://community.freefem.org/u/fb77),

When we solve \Delta^2 u=f, u=0, For DG method, we need third order derivative for the approximation solution. How can we write the term the form \sum\_{K \in T\_h} (\nabla \Delta u, \nabla \Delta v)\_{\partal K} inside the solver environment?, Where T\_h is the triangulations.

---

<div class="post-metadata">

**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\
**Post date:** [August 21, 2026, 6:02pm UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/6 "2026-08-21T18:02:53Z")

</div>

You probably mean for u in at least P3 or P3dc. As commented above, the third order derivative operators are not implemented in FreeFem, thus you need to do something else, as described above.

---

<div class="post-metadata">

**Author:** ![debenswain](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/debenswain/32/3800_2.png) [@debenswain](https://community.freefem.org/u/debenswain)\
**Post date:** [August 24, 2026, 6:33am UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/7 "2026-08-24T06:33:30Z")

</div>

// Third derivative example in FreeFem++  
load “Element\_P3dc”  
load “Element\_P4dc”  
load “Element\_P4”  
verbosity = 0;  
// Define mesh  
mesh Th = square(10, 10);

// Define finite element space (need degree \>= 3 for 3rd derivatives)  
fespace Vh(Th, P3dc);  
Vh u, v;  
real f=0.0;  
real ud = 0.0;

solve Demo(u, v) = - intalledges(Th)(  
( mean(dx(dxx(u)))\*jump(v) +mean(dx(dxx(v)))jump(u) ) / nTonEdge  
)

- int2d(Th)(  
f\*v)

- ind1d(1,2,3,4, Th)(udv)

cout \<\< “L^2 norm” \<\< int2d(Th)((u^2)) \<\< endl;“”

The output shows" 18 : solve Demo(u, v) = - intalledges(Th)(  
19 : ( mean(dx(dxx(u)) current line = 19  
Assertion fail : (0)  
line :418, in file ./../femlib/DOperator.hpp  
error Assertion fail : (0)  
line :418, in file ./../femlib/DOperator.hpp  
code = 6 mpirank: 0

---

<div class="post-metadata">

**Author:** ![fb77](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/fb77/32/3796_2.png) [@fb77](https://community.freefem.org/u/fb77)\
**Post date:** [August 27, 2026, 2:05pm UTC](https://community.freefem.org/t/command-for-third-order-derivative/4057/8 "2026-08-27T14:05:11Z")

</div>

You can replace `dxx(u)` by `u1` in P1dc, and `dxx(v)` by `q1` in P1dc as

```auto
// Third derivative example in FreeFem++
load "Element_P3dc"
load "Element_P4dc"
load "Element_P4"
verbosity = 0;
// Define mesh
mesh Th = square(10, 10);

// Define finite element space (need degree >= 3 for 3rd derivatives)
fespace Vh(Th, P3dc);
Vh u, v;
real f=0.0;
real ud = 0.0;
fespace Vh1(Th,P1dc);
Vh1 u1,v1,p1,q1;

solve Demo(u,u1,p1, v,v1,q1) =
int2d(Th)((u1-dxx(u))*v1)
+int2d(Th)(p1*(q1-dxx(v)))
- intalledges(Th)(
( mean(dx(u1))*jump(v) +mean(dx(q1))*jump(u) ) / nTonEdge
)
-int2d(Th)(f*v)
-int1d(Th,1,2,3,4)(ud*v);
cout << "L^2 norm " << int2d(Th)((u^2)) << endl;

```

It compiles, but does not give a solution because the problem is ill-posed (matrix is singular since the space of solutions is invariant by adding a second-order polynomial).