Error while saving results as xdmf when solving Stokes Equation

Errors
1

I’m trying to solve Lid driven cavity using Stokes equation. I am unable to export the results as xdmf file.

from mpi4py import MPI
import numpy as np
from dolfinx import mesh
from dolfinx import fem
import ufl
from ufl import inner,grad,div,dx
from petsc4py import PETSc
from dolfinx.io import XDMFFile,VTKFile

msh = mesh.create_unit_square(MPI.COMM_WORLD, 40, 40, mesh.CellType.triangle)

# Function to mark x = 0, x = 1 and y = 0
def noslip_boundary(x):
    return np.logical_or(np.logical_or(np.isclose(x[0], 0.0), np.isclose(x[0], 1.0)), np.isclose(x[1], 0.0))

# Function to mark the lid (y = 1)
def lid(x):
    return np.isclose(x[1], 1.0)

# Lid velocity
def lid_velocity_expression(x):
    return np.stack((np.ones(x.shape[1]), np.zeros(x.shape[1])))

P2 = ufl.VectorElement("Lagrange", msh.ufl_cell(), 2)
P1 = ufl.FiniteElement("Lagrange", msh.ufl_cell(), 1)
V, Q = fem.FunctionSpace(msh, P2), fem.FunctionSpace(msh, P1)

# Create the function space
TH = P2 * P1
W = fem.FunctionSpace(msh, TH)
W0, _ = W.sub(0).collapse()

# No slip boundary condition
noslip = fem.Function(V)
facets = mesh.locate_entities_boundary(msh, 1, noslip_boundary)
dofs = fem.locate_dofs_topological((W.sub(0), V), 1, facets)
bc0 = fem.dirichletbc(noslip, dofs, W.sub(0))

# Driving velocity condition u = (1, 0) on top boundary (y = 1)
lid_velocity = fem.Function(W0)
lid_velocity.interpolate(lid_velocity_expression)
facets = mesh.locate_entities_boundary(msh, 1, lid)
dofs = fem.locate_dofs_topological((W.sub(0), V), 1, facets)
bc1 = fem.dirichletbc(lid_velocity, dofs, W.sub(0))

# Since for this problem the pressure is only determined up to a constant, we pin the pressure at the point (0, 0)
zero = fem.Function(Q)
zero.x.set(0.0)
dofs = fem.locate_dofs_geometrical((W.sub(1), Q), lambda x: np.isclose(x.T, [0, 0, 0]).all(axis=1))
bc2 = fem.dirichletbc(zero, dofs, W.sub(1))

# Collect Dirichlet boundary conditions
bcs = [bc0, bc1, bc2]

# Define variational problem
(u, p) = ufl.TrialFunctions(W)
(v, q) = ufl.TestFunctions(W)
f = fem.Function(W0)
a = inner(grad(u), grad(v)) * dx - p * div(v) * dx + div(u) * q * dx
L = inner(f, v) * dx

# Compute solution
solver = fem.petsc.LinearProblem(
    a, L, bcs, petsc_options={"ksp_type": "preonly", "pc_type": "lu",
                              "pc_factor_mat_solver_type": "mumps"})
solver.solve()

# Split the mixed solution and collapse
u = W.sub(0).collapse()
p = W.sub(1).collapse()

# Write the solution to file
with XDMFFile(MPI.COMM_WORLD, "stokesx/pressure.xdmf", "w") as pfile_xdmf:
    p.x.scatter_forward()
    pfile_xdmf.write_mesh(msh)
    pfile_xdmf.write_function(p)

with XDMFFile(MPI.COMM_WORLD, "stokesx/velocity.xdmf", "w") as ufile_xdmf:
    u.x.scatter_forward()
    ufile_xdmf.write_mesh(msh)
    ufile_xdmf.write_function(u)

print('Execution Complete')

This gives an error

p.x.scatter_forward()
AttributeError: 'tuple' object has no attribute 'x'

It would be very helpful if anyone could give me the code for exporting the results.

2

W is not your solution (it should be returned by solver.solve(), it is your function space. Thus the error message.

1 Like