# Cannot find bug in very easy fourth order ODE

**URL:** <https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637>\
**Category:** Uncategorized\
**Created:** [June 27, 2022, 9:51pm UTC](https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637 "2022-06-27T21:51:22Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![KaiJ](https://avatars.discourse-cdn.com/v4/letter/k/f14d63/32.png) [@KaiJ](https://fenicsproject.discourse.group/u/KaiJ)\
**Post date:** [June 27, 2022, 9:51pm UTC](https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637/1 "2022-06-27T21:51:22Z")

</div>

Hi!  
I want to solve a very easy fourth order problem with Robin BC on the 1D interval (0,1) in FEniCS. The weak form says:  
\<u’‘,v’‘\> + 1/b\_1 \* (u(0)v(0)+u(1)v(1)) + 1/b\_2\*(u’(0)v’(0)+u’(1)v’(1)) = \<f,v\>

For the rhs f specified in the code below, I know the true solution true\_sol. Both Expressions should be right. Nevertheless, when I try to plot the numerical solution in comparison to the true solution, the numerical solution is completely false. All the code is attached below.  
As I am a complete FEniCS beginner, I have no clue what I am doing wrong. Probably some very dumb rookie mistake… Thanks very much in advance for any help.

```auto
from fenics import *
import numpy as np
import matplotlib.pyplot as plt

n = 500
mesh = IntervalMesh(n, 0, 1)
V = FunctionSpace(mesh, 'CG', 2) 

b1 = 1
b2 = 1
u = TrialFunction(V)
v = TestFunction(V)
f = Expression('pow((2*pi),4)*sin(2*pi*x[0])', degree = 2)

a = dot(div(grad(u)), div(grad(v)))*dx + pow(b1,-1)*(u*v*ds) + pow(b2,-1)*(dot(grad(u),grad(v))*ds)
L = f*v*dx

sol = Function(V)

solve(a==L, sol)

true_sol = Expression('sin(2*pi*x[0]) - 2*pi*(2*x[0]-1)*(4*b1*(pow(pi,2)*(6*b2-2*pow(x[0],2)+2*x[0]+1)+3)+x[0]*(x[0]-1))*pow(24*b1 + 6*b2 + 1, -1)', domain=mesh, degree = 4, b1 = b1, b2 = b2)

#Plots
dom = np.arange(0,1,0.01)
values_sol = []
values_true = []
for i in dom:
    values_sol.append(sol(i))
    values_true.append(true_sol(i))
plt.plot(dom, values_sol, c='r', label = 'Calculated')
plt.plot(dom, values_true,'b', label = 'True')
plt.legend()
plt.savefig('Figure')

```

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [June 27, 2022, 10:04pm UTC](https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637/2 "2022-06-27T22:04:14Z")

</div>

See for instance: [Biharmonic equation — DOLFIN documentation](https://fenicsproject.org/olddocs/dolfin/2019.1.0/python/demos/biharmonic/demo_biharmonic.py.html)

---

<div class="post-metadata">

**Author:** ![KaiJ](https://avatars.discourse-cdn.com/v4/letter/k/f14d63/32.png) [@KaiJ](https://fenicsproject.discourse.group/u/KaiJ)\
**Post date:** [June 28, 2022, 6:17am UTC](https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637/3 "2022-06-28T06:17:31Z")

</div>

Thanks for the quick answer. Sadly, that file was exactly my starting point: I tried to adjust this to 1D and to my bilinear form/boundary conditions and arrived at the non-working solution above.

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [June 28, 2022, 6:30am UTC](https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637/4 "2022-06-28T06:30:13Z")

</div>

But you have removed all coupling conditions from your form. (i.e. all `dS` integrals, those over internal facets).

---

<div class="post-metadata">

**Author:** ![KaiJ](https://avatars.discourse-cdn.com/v4/letter/k/f14d63/32.png) [@KaiJ](https://fenicsproject.discourse.group/u/KaiJ)\
**Post date:** [June 30, 2022, 8:15am UTC](https://fenicsproject.discourse.group/t/cannot-find-bug-in-very-easy-fourth-order-ode/8637/5 "2022-06-30T08:15:07Z")

</div>

Thanks, don’t know how I could forget that…