[Code scan] Preserve finite s-orbital cube values at the nucleus

[njzjz]

This issue comes from a Codex global scan of deepmodeling/tbplas at commit 4d3652b.

Severity: Low

set_cube() forces every grid point with r < 1.0e-7 to zero. That avoids division by zero for angular functions with l > 0, but it also zeros l == 0 s orbitals at the nucleus, where the value is finite. The resulting cube output contains an artificial node at the nucleus for s states.

Code reference:

tbplas/tbplas/cython/atom.pyx

Lines 193 to 209 in 4d3652b

	xmin, ymin, zmin = cube_origin[0], cube_origin[1], cube_origin[2]
	z, n, l, m = qn[0], qn[1], qn[2], qn[3]
	factor = _rnl_factor(z, n, l) * _ylm_factor(l, m)
	for i in range(num_grid[0]):
	dx = xmin + res * i - c0[0]
	dx2 = dx**2
	for j in range(num_grid[1]):
	dy = ymin + res * j - c0[1]
	dx2y2 = dx2 + dy**2
	for k in range(num_grid[2]):
	dz = zmin + res * k - c0[2]
	r = sqrt(dx2y2 + dz**2)
	if r < 1.0e-7:
	gridval = 0.0
	else:
	gridval = _radial(z, n, l, r) * _ylm(l, m, dx, dy, dz, r) * factor

Suggested fix: handle l == 0 separately at small r, and only force zero for angular components that are genuinely undefined or vanish at the origin.