Matrix-free `getDiagonal` returns 1 on EquationBC rows

Đang mở
  • Có gì không ổn?
#5,641 0 bình luận 0 reaction 0 người được giao Xem trên GitHub

Maintainer thường phản hồi trong vòng 1 ngày

Chưa có ai nhận issue này.

Đánh giá

Độ khó
3/5
Thời gian dự kiến
Nửa ngày
Mức phù hợp với người mới
62/100
Loại issue
Lỗi
Độ rõ ràng
Đặc tả rõ ràng
Mức độ hoạt động
Sôi nổi
Công nghệ
python
Lĩnh vực
backend

Hướng nghiên cứu

The bug is in ImplicitMatrixContext.getDiagonal, which calls bc.set(diagonal, 1) for every bc, including EquationBCSplit. Run the reproduction script in the issue and compare the matfree and aij diagonals on ebc.nodes. Then read how the aij path handles EquationBC rows and check tests/firedrake/multigrid/test_equation_bc.py. Done when the matfree diagonal matches aij on EquationBC rows and the nested test takes 8 and 9 iterations with matfree.

Do mô hình lập chỉ mục viết ra từ nội dung của issue.

Mô tả

For an operator with an EquationBC, getDiagonal of the matrix-free operator returns 1 on every EquationBC row that is not a nested DirichletBC row. The aij matrix has the diagonal of the EquationBC form on these rows. All other rows agree, also the rows of a DirichletBC nested in the EquationBC.

from firedrake import *

mesh = UnitSquareMesh(8, 8)
V = FunctionSpace(mesh, "CG", 2)
u, v = TrialFunction(V), TestFunction(V)
a = inner(grad(u), grad(v)) * dx
e2 = as_vector([0, 1])
a1 = (-inner(dot(grad(u), e2), dot(e2, grad(v))) + 4 * pi**2 * inner(u, v)) * ds(1)
corners = DirichletBC(V, 0, ((1, 3), (1, 4)))

for nested in ([], [corners]):
    ebc = EquationBC(a1 == inner(Constant(0), v) * ds(1), Function(V), 1, bcs=nested).extract_form("J")
    A = assemble(a, bcs=[ebc], mat_type="aij").petscmat
    M = assemble(a, bcs=[ebc], mat_type="matfree").petscmat
    dA, dM = A.getDiagonal(), M.getDiagonal()
    print("aij", dA.array[ebc.nodes[:4]], "matfree", dM.array[ebc.nodes[:4]])

prints

aij [-40.03477216 -36.01738608 -18.00869304 -40.03477216] matfree [1. 1. 1. 1.]
aij [-40.03477216 -36.01738608   1.         -40.03477216] matfree [1. 1. 1. 1.]

ImplicitMatrixContext.getDiagonal assembles the diagonal of the form without the bcs and then calls bc.set(diagonal, 1) for every bc, also for an EquationBCSplit.

This changes the convergence of Jacobi with a matrix-free preconditioning matrix. The converged solution is the same. For example, test_nested_equation_bc in tests/firedrake/multigrid/test_equation_bc.py with mat_type: matfree and the coarse level through AssembledPC takes 22 Krylov iterations instead of 8 in 2D and 15 instead of 9 in 3D. With pmat_type: aij it takes 8 and 9 again.

I see this on release and on main. getDiagonal has been like this since it was added in 2019. An LLM (Claude Code) found this while I was digging into #5609 and the redesign discussion on #5620. G-ADOPT does not use this path but I am flagging it anyway.

Ngôn ngữ chính
Python
Star
684
Fork
204
Merge trung bình
4 ngày 14 giờ
Pull request đã merge (30 ngày)
46

Chuẩn bị môi trường

  • Không có Dockerfile hay tệp Docker Compose
  • Có mẫu pull request
  • Không có hướng dẫn đóng góp

Bắt đầu từ đâu

  1. Đọc hết issue, rồi đọc hướng dẫn đóng góp của dự án.
  2. Bình luận trên issue rằng bạn sẽ nhận — tránh hai người làm cùng một việc.
  3. Fork repository và làm thay đổi trên một nhánh.
  4. Mở pull request có tham chiếu số hiệu của issue.

Issue khác của firedrakeproject/firedrake

Tất cả issue của firedrakeproject/firedrake

Issue tương tự

Thêm issue về Python

Nhận issue mới trong hộp thư của bạn

Bản tóm tắt ngắn những issue GitHub phù hợp với người mới.