Skip to content

Bug: s_axis QBMM buffer reflection uses wrong azimuthal index for z_cc(l) < pi #1357

@sbryngelson

Description

@sbryngelson

Description

In s_axis() (src/common/m_boundary_common.fpp), the primitive variable reflection correctly branches on z_cc(l) < pi to choose the azimuthal source index:

if (z_cc(l) < pi) then
    q_prim_vf(i)%sf(k, -j, l) = q_prim_vf(i)%sf(k, j - 1, l + ((p + 1)/2))
else
    q_prim_vf(i)%sf(k, -j, l) = q_prim_vf(i)%sf(k, j - 1, l - ((p + 1)/2))
end if

However, the QBMM buffer reflection block below it always uses l - ((p + 1)/2) regardless of z_cc(l):

if (qbmm .and. .not. polytropic .and. present(pb_in) .and. present(mv_in)) then
    do i = 1, nb
        do q = 1, nnode
            do j = 1, buff_size
                pb_in(k, -j, l, q, i) = pb_in(k, j - 1, l - ((p + 1)/2), q, i)  ! always l - half
                mv_in(k, -j, l, q, i) = mv_in(k, j - 1, l - ((p + 1)/2), q, i)
            end do
        end do
    end do
end if

For z_cc(l) < pi, this reads from l - ((p+1)/2) which can be an invalid/negative index and copies from the wrong azimuthal location. The QBMM block should mirror the same z_cc(l) < pi branching as the primitive variables.

Suggested fix

if (qbmm .and. .not. polytropic .and. present(pb_in) .and. present(mv_in)) then
    do i = 1, nb
        do q = 1, nnode
            do j = 1, buff_size
                if (z_cc(l) < pi) then
                    pb_in(k, -j, l, q, i) = pb_in(k, j - 1, l + ((p + 1)/2), q, i)
                    mv_in(k, -j, l, q, i) = mv_in(k, j - 1, l + ((p + 1)/2), q, i)
                else
                    pb_in(k, -j, l, q, i) = pb_in(k, j - 1, l - ((p + 1)/2), q, i)
                    mv_in(k, -j, l, q, i) = mv_in(k, j - 1, l - ((p + 1)/2), q, i)
                end if
            end do
        end do
    end do
end if

Impact

Affects 3D cylindrical QBMM simulations with axis boundary conditions. Pre-existing bug on master, not introduced by any recent PR. Surfaced during code review of #1331 and #1356.

Metadata

Metadata

Assignees

Labels

bugSomething isn't working or doesn't seem right

Type

No type
No fields configured for issues without a type.

Projects

No projects

Milestone

No milestone

Relationships

None yet

Development

No branches or pull requests

Issue actions