### Write i-coordinate in FEniCS

217
views
0
7 months ago by
Dear FEniCS community. How to write the momentum equation with i-coordinate in FEniCS?

where u = (u,v).

Community: FEniCS Project

0
7 months ago by
Dear Raihan,

in general, in order to get the x-equation of a vector-variable equation, you have to multiply with the x unit vector. This means that the resulting x-equation will be scalar. So the trial and test functions of the Navier-Stokes problem must live on a scalar FunctionSpace from now on. This is the main difference.

So in order to write the x-momentum you will define first trial and test functions that are scalar and representing u_x. Also, the second difference is that you must define boundary conditions for u_x and not for u. That is all you need.

There is a cool feature in Fenics, called split that allows you to extract the scalar components of a vector. You might need it if you have any initial or boundary conditions for u instead of u_x. This is an example code of a simple vector PDE, solved for x-component, with vector initial condition that needs to be splitted in order to get the initial condition for u_x.

dt = 1e-3

mesh = UnitSquareMesh(100, 100)

V = VectorFunctionSpace(mesh, 'P', 3)

Q = FunctionSpace(mesh, 'P', 1)

u = TrialFunction(Q)

v = TestFunction(Q)

u_n = Function(V)

u_n = interpolate(Expression(('10.0', '1.0'), degree=2), V)

u_x_n, u_y_n = u_n.split(deepcopy = True)

bc = DirichletBC(Q, Constant((1.0)), 'on_boundary')

A = assemble(a1)

RHS_x = 1.0/dt*u_x_n*v*dx

b1 = assemble(RHS_x)

bc.apply(A, b1)

u = Function(Q)

t = 0.0

while error > 1e-3:
t += dt
start_time3 = time.time()
b1 = assemble(RHS_x)
bc.apply(b1)
solve(a1 == RHS_x, u, bc)
#CPUtime_solve = time.time() - start_time3
#print('Time in seconds:', CPUtime_solve)
error = abs(norm(u) - norm(u_x_n))
print(error)

u_x_n.assign(u)

CPUtime = time.time() - start_time​