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).

Thank you for your help.

Community: FEniCS Project

1 Answer


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')

a1 = 1.0/dt*u*v*dx + dot(nabla_grad(u), nabla_grad(v))*dx

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​
Please login to add an answer/comment or follow this question.

Similar posts:
Search »