Reproduce divide by 0 bug in FiPy

fipy_issue.py

import fipy
import numpy as np

# Problem definition
num_cell = 2
num_face = num_cell + 1
init_u = np.zeros(num_cell)
dx = 1.0
dt = 1.0
tc_cell = np.ones(num_cell)
d_face = np.zeros(num_face)
v_face = np.array([1.0, 1.0, 0.0])
left_face_constraint = 1.0
right_face_constraint = 1.0

# Solve the problem with FiPy
mesh = fipy.Grid1D(
    nx=num_cell,
    dx=dx,
)
u = fipy.CellVariable(
    mesh=mesh,
    value=init_u,
)
u.constrain(left_face_constraint, where=mesh.facesLeft)
u.constrain(right_face_constraint, where=mesh.facesRight)
transient_coeff = fipy.CellVariable(mesh=mesh, value=tc_cell)
transient = fipy.TransientTerm(coeff=transient_coeff, var=u)
diffusion_coeff = fipy.FaceVariable(mesh=mesh, value=d_face)
diffusion = fipy.DiffusionTerm(coeff=diffusion_coeff, var=u)
convection_coeff = fipy.FaceVariable(mesh=mesh, value=np.expand_dims(v_face, 0))
convection = fipy.ConvectionTerm(coeff=convection_coeff, var=u)
eq = transient + convection == diffusion
eq.solve(dt=dt)
print(u)