Skip to content

LoopyIndexError: Real coefficient inside jump() on interior facet integral #4991

@sghelichkhani

Description

@sghelichkhani

After commit 1bd729a ("Real x extruded fixes #4970"), assembling an interior facet integral where a Real-space coefficient appears inside a jump() (or any facet restriction) fails with a loopy out-of-bounds error.

The commit removed the TSFC-side special casing that prevented Real coefficients from being facet-doubled, but firedrake/assemble.py was not updated to match. Specifically, _as_global_kernel_arg_coefficient (line ~1863) still returns op2.GlobalKernelArg((V.value_size,)) with shape (1,), and _as_parloop_arg_coefficient (line ~2213) still returns op2.GlobalParloopArg(coeff.dat) with undoubled data. The kernel now expects shape (2,) since TSFC generates facet-side indexing for the Real coefficient (w[0] for +, w[1] for -).

Minimal reproducer

from firedrake import *

mesh = UnitSquareMesh(2, 2)
K = FunctionSpace(mesh, "DG", 1)
R = FunctionSpace(mesh, "Real", 0)

q = Function(K)
v = TestFunction(K)
c = Function(R)
c.assign(1.0)

F = jump(c * q) * jump(v) * dS
assemble(F)  # LoopyIndexError

Note that c * jump(q) * jump(v) * dS (Real outside the jump) works fine, because c is never restricted.

Error

loopy.diagnostic.LoopyIndexError: 'w_1[1]' in instruction 'insn_4' accesses
out-of-bounds array element (could not establish '{ [1] }' is a subset of
'{ [i0] : i0 = 0 }').

Impact

This breaks Irksome's DIRK time steppers for DG advection. Irksome constructs stage values like q_old + dt * a11 * q where dt and a11 are Real-space Functions, and the DG upwind form restricts the whole expression with (+)/(-), putting the Real coefficients inside the restriction. This is currently failing in G-ADOPT's CI for all multi-material / level-set tests.

Metadata

Metadata

Assignees

Labels

Type

Projects

No projects

Milestone

No milestone

Relationships

None yet

Development

No branches or pull requests

Issue actions