# Weak Form of Higher order derivatives **URL:** \ **Category:** General Discussion\ **Created:** [September 5, 2021, 1:22pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190 "2021-09-05T13:22:47Z")\ **Posts on this page:** 19\ **Page:** 1
**Author:** ![sumantkr](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/sumantkr/32/2681_2.png) [@sumantkr](https://community.freefem.org/u/sumantkr)\ **Post date:** [September 5, 2021, 1:22pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/1 "2021-09-05T13:22:47Z")
Dear all Can we write third or fourth-order derivatives in FreeFEM? However, it is showing **identifiers errors** while writing third-order derivatives ( like dxxx( u ) ) . Kindly suggest me in solving this issue. ---
**Author:** ![sumantkr](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/sumantkr/32/2681_2.png) [@sumantkr](https://community.freefem.org/u/sumantkr)\ **Post date:** [September 6, 2021, 10:02am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/2 "2021-09-06T10:02:41Z")
Applying the Greens theorem multiple times will reduce the degree of PDEs term. It solved my problem. ---
**Author:** ![VaishnaviGnanam](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/vaishnavignanam/32/2613_2.png) [@VaishnaviGnanam](https://community.freefem.org/u/VaishnaviGnanam)\ **Post date:** [August 28, 2024, 8:40am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/3 "2024-08-28T08:40:15Z")
What about 6th order? ---
**Author:** ![sumantkr](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/sumantkr/32/2681_2.png) [@sumantkr](https://community.freefem.org/u/sumantkr)\ **Post date:** [August 28, 2024, 9:00am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/4 "2024-08-28T09:00:34Z")
You may use macro, such as: macro grad(u) [dx(u),dy(u)]// macro div(u1,u2) (dx(u1)+dy(u2))// macro Vgradu(v1,v2,u) (grad(u)'\*[v1,v2])// macro Grad(u1,u2) [grad(u1),grad(u2)]// macro VgradU(v1,v2,u1,u2) [Vgradu(v1,v2,u1),Vgradu(v1,v2,u2)]// ---
**Author:** ![VaishnaviGnanam](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/vaishnavignanam/32/2613_2.png) [@VaishnaviGnanam](https://community.freefem.org/u/VaishnaviGnanam)\ **Post date:** [August 28, 2024, 10:36am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/5 "2024-08-28T10:36:37Z")
In 1d case suppose we want to write u\_xxxxxx means in weak formulation we can pass the the derivative and we can get (u\_xxx, v\_xxx) or more we can pass. Now how we can tackle this? ---
**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 28, 2024, 12:02pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/6 "2024-08-28T12:02:59Z")
You can write your equation in a system form. For example if you want to solve u\_{xxxxxx}=f, you can write u\_x=u\_1,\quad (u\_1)\_x=u\_2,\quad (u\_2)\_x=u\_3,\quad (u\_3)\_x=u\_4,\quad (u\_4)\_x=u\_5,\quad (u\_5)\_x=f. Then you have the unknowns u,u\_1,u\_2,u\_3,u\_4,u\_5 in P1 and you write the variational formulation for v,v\_1,v\_2,v\_3,v\_4,v\_5 in P1 \int \partial\_x(u\_i) v\_i=\int u\_{i+1}v\_i,\quad i=0,...,4, \int \partial\_x(u\_5)v\_5=\int fv\_5. In order to have a well-posed system you need to introduce boundary conditions, by adding boundary terms. A 2d example is to solve the Bilaplace problem \Delta^2 u=f=1 in the unit disc with boundary conditions u=0 and \nabla u=0 on the boundary. [bilaplacian.edp](https://community.freefem.org/uploads/short-url/wVaQSgyarC0E23T6Ckm2FIk2WZf.edp) (1.0 KB) In this file two methods are used. The first is by the Morley finite element space. The second is by the above method, where you have the unknowns ub, uxb=dx(ub) uyb=dy(ub) ulb=dx(uxb)+dy(uyb) ---
**Author:** ![VaishnaviGnanam](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/vaishnavignanam/32/2613_2.png) [@VaishnaviGnanam](https://community.freefem.org/u/VaishnaviGnanam)\ **Post date:** [August 29, 2024, 5:55am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/7 "2024-08-29T05:55:46Z")
Dear sir , sorry for the inconvenience, i couldn’t understand it properly. can you please explain it with 1d case ? Thank you. ---
**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 29, 2024, 12:26pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/8 "2024-08-29T12:26:27Z")
As an example we can consider the 1d problem \partial^4\_x u=1 in (0,1) with boundary condition u=0 and \partial\_x u=0 on the whole boundary, which is the bilaplacian problem in 1d. The solution is given by u=x^2(1-x)^2/24. You can solve it by introducing the unknowns u, u1=dx(u), u2=dx(u1) It gives ```auto solve bilap1d(u,u1,u2,v,v1,v2)= -int1d(Th)(dx(u2)*dx(v2)) +int0d(Th)(dx(u2)*v2*N.x) -int1d(Th)(v2) +int1d(Th)((dx(u)-u1)*v) +int1d(Th)(-u1*dx(v1)-u2*v1) +int0d(Th)(1.e30*u*v2) ; ``` Whole file is [bilap1d.edp](https://community.freefem.org/uploads/short-url/e6e93AB8o87B3Wk0qfFCvDujhrx.edp) (818 Bytes) It contains also another possible resolution with only two unknowns u,u2 ---
**Author:** ![VaishnaviGnanam](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/vaishnavignanam/32/2613_2.png) [@VaishnaviGnanam](https://community.freefem.org/u/VaishnaviGnanam)\ **Post date:** [August 30, 2024, 5:56am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/9 "2024-08-30T05:56:47Z")
Dear sir, I am very new to freefem . So in the context, bilap1d(u,u1,u2,v,v1,v2) , if we say u1 means does it consider it as dx(u) and so on? and can you please explain about the code you provided ? Thank you. ---
**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 30, 2024, 10:18am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/10 "2024-08-30T10:18:16Z")
The variational formulation that is solved by bilap1d() is: find u, u\_1, u\_2 in P1 such that for all v, v\_1, v\_2 in P1 -\int \partial\_x u\_2\partial\_x v\_2+\int\_\Gamma \partial\_x u\_2 v\_2 N-\int v\_2 +\int (v\partial\_x u - v u\_1 ) +\int(-u\_1\partial\_x v\_1 - u\_2 v\_1)+\int\_\Gamma \xi u v\_2\ =0, where \int stands for \int\_0^1, \int\_\Gamma f=f(0)+f(1), and N is the exterior normal i.e. N(0)=-1, N(1)=1, and \xi=10^{30} is a large coefficient. In this formulation, you can take successively v\_1=0, v\_2=0, v arbitrary, then v=0, v\_2=0, v\_1 arbitrary, and finally v=0, v\_1=0, v\_2 arbitrary. The result is that it is equivalent to write -\int \partial\_x u\_2\partial\_x v\_2+\int\_\Gamma \partial\_x u\_2 v\_2 N-\int v\_2+\int\_\Gamma \xi u v\_2\ =0 for any v\_2 in P1, \int (v\partial\_x u - v u\_1 )=0 for any v in P1, \int(-u\_1\partial\_x v\_1 - u\_2 v\_1)=0 for any v\_1 in P1. The second line is an approximation of \partial\_x u=u\_1. The third line is an approximation of \partial\_x u\_1=u\_2 and u\_1(0)=0, u\_1(1)=0, because \int -u\_1\partial\_x v\_1=\int v\_1\partial\_x u\_1-\int\_\Gamma u\_1 v\_1. The first line is an approximation of \partial^2\_{xx}u\_2=1 and u(0)=0, u(1)=0, because -\int \partial\_x u\_2\partial\_x v\_2+\int\_\Gamma \partial\_x u\_2 v\_2 N=\int v\_2\partial^2\_{xx}u\_2 (but notice that \partial^2\_{xx} u\_2 is not a function, it contains Dirac masses). * * * Another possibility which is closer to the formulation of my first message (and maybe easier to understand) is to use the unknowns u, u1=dx(u), u2=dx(u1), u3=dx(u2). This gives ```auto fespace V(Th,P1); V u,u1,u2,u3,v,v1,v2,v3; solve bilap1d(u,u1,u2,u3,v,v1,v2,v3)= int1d(Th)((dx(u2)-u3)*v2) +int1d(Th)(dx(u3)*v3) -int1d(Th)(v3) +int1d(Th)((dx(u)-u1)*v) +int1d(Th)(-u1*dx(v1)-u2*v1) +int0d(Th)(1.e30*u*v3) ; ``` ---
**Author:** ![VaishnaviGnanam](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/vaishnavignanam/32/2613_2.png) [@VaishnaviGnanam](https://community.freefem.org/u/VaishnaviGnanam)\ **Post date:** [August 31, 2024, 11:53am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/11 "2024-08-31T11:53:31Z")
Sir , Thank you so much . one doubt sir, for the given equation \partial^4\_xu=1 how you are getting that weak formulation? and \partial\_xu2=\partial^3\_xu? Thank you. ---
**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 31, 2024, 12:19pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/12 "2024-08-31T12:19:06Z")
At the continuous level you just replace \partial^4\_x u=1 by the set of equations \partial\_x u= u\_1, \ \partial\_x u\_1=u\_2, \ \partial\_x u\_2=u\_3, \ \partial\_x u\_3=1. The first equation gives u\_1=\partial\_x u, thus replacing u\_1 by this value in the second equation gives u\_2=\partial^2\_x u. Then replacing u\_2 by this value in the third equation gives u\_3=\partial^3\_x u. Finally replacing u\_3 by this value in the fourth equation gives \partial^4\_x u=1. ---
**Author:** ![VaishnaviGnanam](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/vaishnavignanam/32/2613_2.png) [@VaishnaviGnanam](https://community.freefem.org/u/VaishnaviGnanam)\ **Post date:** [September 6, 2024, 5:20am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/13 "2024-09-06T05:20:21Z")
Thank you so much sir. Suppose if we are having u\_t,u\_x and u may be some non linear term f(u) also means how to deal with that sir ? Thank you. ---
**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 9, 2024, 10:53am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/14 "2024-09-09T10:53:26Z")
For time dependence and nonlinearities in principle you can use standard approaches. In particular for an evolution problem with u\_t and a Cauchy initial data you need a time scheme. The most simple is to use implicit Euler. For nonlinearities you have to solve the problem by an iterative algorithm like the Newton method. The issue however is stability, and that depends strongly on the equation you want to solve and on the scheme you use. ---
**Author:** ![RYP](https://avatars.discourse-cdn.com/v4/letter/r/b9e5f3/32.png) [@RYP](https://community.freefem.org/u/RYP)\ **Post date:** [June 28, 2026, 10:40am UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/15 "2026-06-28T10:40:20Z")
Hello, sir, I would like to ask why we use`u*v2`\*？`u*v1`_or_`u*v`\*does not seem to work. Besides, can we use`on`to specify boundary conditions? Looking forward to your reply. Thank you. ---
**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:** [June 28, 2026, 12:59pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/16 "2026-06-28T12:59:31Z")
I don’t see to which code lines you refer to… ---
**Author:** ![RYP](https://avatars.discourse-cdn.com/v4/letter/r/b9e5f3/32.png) [@RYP](https://community.freefem.org/u/RYP)\ **Post date:** [June 28, 2026, 1:14pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/17 "2026-06-28T13:14:56Z")
```auto solve bilap1d(u,u1,u2,v,v1,v2)= -int1d(Th)(dx(u2)*dx(v2)) +int0d(Th)(dx(u2)*v2*N.x) -int1d(Th)(v2) +int1d(Th)((dx(u)-u1)*v) +int1d(Th)(-u1*dx(v1)-u2*v1) +int0d(Th)(1.e30*u*v2) ; ``` The penalty term in the last line, it does not seem to work when written as `1.e30*u*v1`or`1.e30*u*v`. I wonder why it should be `u*v2`. In addition, can we use`on`to specify boundary conditions in this problem? Thank you in advance. ---
**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:** [June 28, 2026, 2:57pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/18 "2026-06-28T14:57:24Z")
Indeed I don’t know why taking `u*v2` is stable, and not the other choices. If we use `on(u=...)` it puts a term on the diagonal of the matrix, hence uses the variable which is declared dual of `u` in `bilap1d(u,u1,u2,v,v1,v2)`. It is thus equivalent to taking `u*v`, which is not stable. You could use `on(u=...)` with the stable variable `v2` if you state `bilap1d(u,u1,u2,v2,v1,v)` ---
**Author:** ![RYP](https://avatars.discourse-cdn.com/v4/letter/r/b9e5f3/32.png) [@RYP](https://community.freefem.org/u/RYP)\ **Post date:** [June 28, 2026, 4:24pm UTC](https://community.freefem.org/t/weak-form-of-higher-order-derivatives/1190/19 "2026-06-28T16:24:44Z")
Thank you very much for your explanation，I will give it a try.