Re: [firedrake] Linear problems with zero right hand side
Hi Daniel, Nice to see you on here. The "F == 0" syntax is usually reserved for nonlinear problems, though it can also be used for linear problems (which I think is what you have?) However, when using this syntax, F should only contain Functions and TestFunctions, not TrialFunctions. The Jacobian (which corresponds to the LHS if F was linear) is then generated automatically using automatic differentiation of F with respect to an appropriate variable. E.g., in the Burgers example, we have F = (inner((u - u_)/timestep, v) + inner(dot(u,nabla_grad(u)), v) + nu*inner(grad(u), grad(v)))*dx Note that u was declared as a Function. Then a few lines later, while (t <= end): solve(F == 0, u) The second argument to the solve call, u, is the variable to solve for (and where the result gets put), and the AD of F is done with respect to this variable. However, when using the "a == L" syntax for linear problems, a must contain TrialFunctions and TestFunctions (and, optionally, Functions), while L must contain TestFunctions (and usually Functions), and so your other variant works. Best, Andrew On 25 April 2016 at 08:59, Daniel Ruprecht <D.Ruprecht@leeds.ac.uk> wrote:
Dear Firedrake team,
I just began to use firedrake and tried to implement a simple heat equation to get started. I encountered some strange (strange for me that is) behaviour when using zero as a right hand side in a linear problem.
Following the nonlinear Burgers' equation example on the website, I tried to write the diffusion equation u_t = nu*Delta(u) with backward Euler as
a = (inner((u - u_)/timestep, v) + nu*inner(grad(u), grad(v))) * dx
(removing the nonlinear terms essentially) but when solving
solve(a == 0, temp)
this generates a ValueError: Provided residual is not a linear form.
If instead I put the u_ term on the right hand side
a = (inner(u, v) + timestep * nu * inner(grad(u), grad(v))) * dx L = inner(u_, v) * dx solve(a == L, temp)
everything works fine. Is this working as intended? I could not find a comment about this on the website, I was confused because the nonlinear example seemed to suggest the first variant should be the way to do it.
I attached the script I used to provide a minimum example. This is with firedrake/2016-03, i.e. the version from March this year and Python 2.7.9.
Cheers, Daniel
-- Dr. Daniel Ruprecht University Academic Fellow School of Mechanical Engineering, Room 4.50 University of Leeds Leeds LS2 9JT, UK Email: d.ruprecht@leeds.ac.uk Phone: +44 (0)113 343 22 01
Hi Andrew, Thanks for clarifying this - I suspected there might be a deeper reason, but could not figure it out. I am quite happy using the second variant, I just wanted to understand why the first one was not working. Best, Daniel On 25/04/16 09:14, Andrew McRae wrote:
Hi Daniel,
Nice to see you on here.
The "F == 0" syntax is usually reserved for nonlinear problems, though it can also be used for linear problems (which I think is what you have?)
However, when using this syntax, F should only contain Functions and TestFunctions, not TrialFunctions. The Jacobian (which corresponds to the LHS if F was linear) is then generated automatically using automatic differentiation of F with respect to an appropriate variable. E.g., in the Burgers example, we have
F = (inner((u - u_)/timestep, v) + inner(dot(u,nabla_grad(u)), v) + nu*inner(grad(u), grad(v)))*dx
Note that u was declared as a Function. Then a few lines later,
while (t <= end): solve(F == 0, u)
The second argument to the solve call, u, is the variable to solve for (and where the result gets put), and the AD of F is done with respect to this variable.
However, when using the "a == L" syntax for linear problems, a must contain TrialFunctions and TestFunctions (and, optionally, Functions), while L must contain TestFunctions (and usually Functions), and so your other variant works.
Best, Andrew
On 25 April 2016 at 08:59, Daniel Ruprecht <D.Ruprecht@leeds.ac.uk <mailto:D.Ruprecht@leeds.ac.uk>> wrote:
Dear Firedrake team,
I just began to use firedrake and tried to implement a simple heat equation to get started. I encountered some strange (strange for me that is) behaviour when using zero as a right hand side in a linear problem.
Following the nonlinear Burgers' equation example on the website, I tried to write the diffusion equation u_t = nu*Delta(u) with backward Euler as
a = (inner((u - u_)/timestep, v) + nu*inner(grad(u), grad(v))) * dx
(removing the nonlinear terms essentially) but when solving
solve(a == 0, temp)
this generates a ValueError: Provided residual is not a linear form.
If instead I put the u_ term on the right hand side
a = (inner(u, v) + timestep * nu * inner(grad(u), grad(v))) * dx L = inner(u_, v) * dx solve(a == L, temp)
everything works fine. Is this working as intended? I could not find a comment about this on the website, I was confused because the nonlinear example seemed to suggest the first variant should be the way to do it.
I attached the script I used to provide a minimum example. This is with firedrake/2016-03, i.e. the version from March this year and Python 2.7.9.
Cheers, Daniel
-- Dr. Daniel Ruprecht University Academic Fellow School of Mechanical Engineering, Room 4.50 University of Leeds Leeds LS2 9JT, UK Email: d.ruprecht@leeds.ac.uk <mailto:d.ruprecht@leeds.ac.uk> Phone: +44 (0)113 343 22 01 <tel:%2B44%20%280%29113%20343%2022%2001>
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
-- Dr. Daniel Ruprecht University Academic Fellow School of Mechanical Engineering, Room 4.50 University of Leeds Leeds LS2 9JT, UK Email: d.ruprecht@leeds.ac.uk Phone: +44 (0)113 343 22 01
On 25/04/16 09:25, Daniel Ruprecht wrote:
Hi Andrew,
Thanks for clarifying this - I suspected there might be a deeper reason, but could not figure it out. I am quite happy using the second variant, I just wanted to understand why the first one was not working.
For what it's worth, I would recommend always formulating your problem in residual form (even if it is linear). The linear variational solvers just wrap up the nonlinear variational ones anyway (with appropriate options). For a linear problem, newton will converge in one step and you might (depending on what you're monitoring with petsc) pay at most one extra residual evaluation. Writing your formulation in this way makes it much easier to deal with nonlinearities once you have them: no need to go back and turn all TrialFunctions into Functions and so on. To be fair, none of our demos do this! Cheers, Lawrence
participants (3)
-
Andrew McRae
-
Daniel Ruprecht
-
Lawrence Mitchell