how to give Jocabian matrix?

[tu1620808124]

When solving a nonlinear system such as F(X) == 0, I want to give a Jacobian explicitly, but I don't know how to express it in the code here. if there have a simple example?

[wence-]

You can provide an explicit Jacobian by passing J=your_bilinear_form to the solve call (or, if you're building a NonlinearVariationalProblem to that).

See https://firedrakeproject.org/solving-interface.html#nonlinear-variational-problems, so something like:

solve(F == 0, X, J=some_bilinear_form)

[tu1620808124]

For example, a simple Poisson problem is solved by Newton's method. Is this J = derivative(F, u, du) an analytic jacobian matrix?

import os
import sys
import datetime
Nx =int(sys.argv[1])
mesh = UnitSquareMesh(Nx,Nx) 
V = FunctionSpace(mesh, "CG", 1)

bc = DirichletBC(V, Constant(0), [1,2,3,4])
x,y= mesh.coordinates
f = interpolate(8*pow(pi,2)*sin(2*pi*x)*sin(2*pi*y), V)

u = Function(V)
v = TestFunction(V)
a = inner(grad(u), grad(v))*dx
l = f*v*dx 
F = a - l


du = TrialFunction(V)
J = derivative(F, u, du)
problem = NonlinearVariationalProblem(F, u, bc, J)
solver = NonlinearVariationalSolver(problem, solver_parameters = parameters)

starttime = datetime.datetime.now() 
solver.solve()
endtime = datetime.datetime.now()
print(endtime -starttime)

truth = sin(2*pi*x)*sin(2*pi*y)
#File("cg_tet_result/u_poisson.pvd").write(u)
print(errornorm(truth, u, norm_type="L1")/norm(truth, norm_type = "L1"))

[wence-]

Yes, derivative(F, u, du) provides the analytic Gateaux derivative.