help with a problem #2867
Unanswered
Andre17420
asked this question in
Firedrake support
Replies: 2 comments 2 replies
-
Hello Andre,
What is your question?
all the best
cjc
On 9 Apr 2023, at 11:22, Andre17420 ***@***.***> wrote:
Hello to everybody,
I am a student from Milan and I have just started facing numerical methods. I gotta solve this parabolic problem, but I'm stuck in checking the stability and the convergency. I attached two pictures of the problem. Thank you all.
Here is what I have done so far:
Define the mesh
n = 20
mesh = UnitSquareMesh(n, n, diagonal='crossed')
Define the function space
V = FunctionSpace(mesh, 'Lagrange', 1)
Boundary condition
g = Constant(1.0)
bc = DirichletBC(V, g, "on_boundary")
Problem
x = SpatialCoordinate(mesh)
f_fun = lambda t: (pi/2cos(pit) + pi**2sin(pit)) * sin(pix[0]) * sin(pix[1])
T = 1.
dt = 0.1
theta = 0.5
Exact solution
u_ex = lambda t: 1 + 0.5 * sin(pit) * sin(pix[0]) * sin(pi*x[1])
Variational formulation
u = TrialFunction(V)
v = TestFunction(V)
t = 0
u_old = interpolate(u_ex(t), V)
a = 1/dt * u * v * dx + theta * dot(grad(u), grad(v)) * dx
Problem definition
u_h = Function(V, name='Temperature')
u_h.assign(u_old)
start_time = perf_counter()
A = assemble(a, bcs=bc)
solver = LinearSolver(A,
solver_parameters={"ksp_type": "preonly",
"pc_type": "lu"})
Time stepping
errL2 = 0
errH1 = 0
while t+dt/2<T:
t_old = t
t += dt
print('t = ', t)
L = ( 1/dt * u_old * v * dx - (1-theta) * (dot(grad(u_old), grad(v)) * dx)
+ theta * f_fun(t) * v * dx + (1-theta) * (f_fun(t_old) * v * dx) )
b = assemble(L, bcs=bc)
solver.solve(u_h, b)
u_old = u_h
#outfile.write(u_h, time=t)
# Compute the error
errL2 = max(errL2, errornorm(u_ex(t), u_h, 'L2'))
errH1 = max(errH1, errornorm(u_ex(t), u_h, 'H1'))
exec_time = perf_counter() - start_time
[g1]<https://user-images.githubusercontent.com/130237366/230764888-bc73995f-3cac-4745-99f9-a10f70958752.png>
[image]<https://user-images.githubusercontent.com/130237366/230764921-40ebe7f8-75a2-4046-9bde-45452e9cf8a3.png>
—
Reply to this email directly, view it on GitHub<#2867>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/ABOSV4QOAQ6ZJS5HO7IX7F3XAJ5T5ANCNFSM6AAAAAAWX75RMU>.
You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>
|
Beta Was this translation helpful? Give feedback.
1 reply
-
It sounds like this is a conceptual problem rather than a Firedrake one. I think that the best approach to learning this would be to look at the definitions of stability and convergence that you have been provided with and think about what you would need to measure in order to check if they are satisfied.
On 9 Apr 2023, at 10:56, Andre17420 ***@***.***> wrote:
I'm not able to check the stability and the convergence, I'm asking wheter anybody could help me out
—
Reply to this email directly, view it on GitHub<#2867 (reply in thread)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/ABOSV4XGDZ2FCQS5KHLNKYTXAKBVVANCNFSM6AAAAAAWX75RMU>.
You are receiving this because you commented.Message ID: ***@***.***>
|
Beta Was this translation helpful? Give feedback.
1 reply
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Uh oh!
There was an error while loading. Please reload this page.
-
Hello to everybody,
I am a student from Milan and I have just started facing numerical methods. I gotta solve this parabolic problem, but I'm stuck in checking the stability and the convergency. I attached two pictures of the problem. Thank you all.
Here is what I have done so far:
Define the mesh
n = 20
mesh = UnitSquareMesh(n, n, diagonal='crossed')
Define the function space
V = FunctionSpace(mesh, 'Lagrange', 1)
Boundary condition
g = Constant(1.0)
bc = DirichletBC(V, g, "on_boundary")
Problem
x = SpatialCoordinate(mesh)
f_fun = lambda t: (pi/2cos(pit) + pi**2sin(pit)) * sin(pix[0]) * sin(pix[1])
T = 1.
dt = 0.1
theta = 0.5
Exact solution
u_ex = lambda t: 1 + 0.5 * sin(pit) * sin(pix[0]) * sin(pi*x[1])
Variational formulation
u = TrialFunction(V)
v = TestFunction(V)
t = 0
u_old = interpolate(u_ex(t), V)
a = 1/dt * u * v * dx + theta * dot(grad(u), grad(v)) * dx
Problem definition
u_h = Function(V, name='Temperature')
u_h.assign(u_old)
start_time = perf_counter()
A = assemble(a, bcs=bc)
solver = LinearSolver(A,
solver_parameters={"ksp_type": "preonly",
"pc_type": "lu"})
Time stepping
errL2 = 0
errH1 = 0
while t+dt/2<T:
t_old = t
t += dt
print('t = ', t)
exec_time = perf_counter() - start_time
Beta Was this translation helpful? Give feedback.
All reactions