A simple time-dependent 2d code


230
views
-1
6 months ago by
Hello dear FEniCS users, I am extremely new to FEniCS and I have a problem (that is probably so easy for you)

I wanna solve the PDE ut+ux+uy=0 and I would like to ask that the following code is correct or not (actually I am asking a= u*v*dx + dt*v*(u.dx(0) + u.dx(1))*dx is correct or not). Could you please help ? Thank you so much in advance:


from fenics import *

import numpy as np

T = 2.0 # final time

num_steps = 10 # number of time steps

dt = T / num_steps # time step size

# Create mesh and define function space

nx = ny = 10

mesh = UnitSquareMesh(nx, ny)

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

# Define boundary condition

u_D = Expression('exp(x[0]-x[1])',

degree=1, t=0)

def boundary(x, on_boundary):

return on_boundary

bc = DirichletBC(V, u_D, boundary)

# Define initial value

#u_n = interpolate(u_D, V)

u_n = project(u_D, V)

# Define variational problem

u = TrialFunction(V)

v = TestFunction(V)

a= u*v*dx + dt*v*(u.dx(0) + u.dx(1))*dx

L= (u_n)*v*dx

# Time-stepping

u = Function(V)

t=0

for n in range(num_steps):

# Update current time

t = t+dt

u_D.t = t

# Compute solution

solve(a == L, u, bc)

# Plot solution

plot(u)

# Compute error at vertices

u_e = interpolate(u_D, V)

error = np.abs(u_e.vector().array() - u.vector().array()).max()

print('t = %.2f: error = %.3g' % (t, error))

# Update previous solution

u_n.assign(u)

vtkfile = File('solution.pvd')

vtkfile << (u,t)
Community: FEniCS Project
Hi, if you edit your post, you will find that there is a button to format your code in the proper way. That would certainly help in having your question answered. Plus, you might want to add details about what you are trying to do and what you have done.
written 6 months ago by François Lapointe  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »