Yes, the whole thing is meant to insert the necessary singleton dimensions so that broadcasting works. A use case could be the following:
# I have ten displacement vectors a 3 elements:
displacements = np.random.rand(10, 3)
# I have fifteen locations, also 3 element vectors
locations = np.random.rand(15, 3)
# now, the displacements are actually supposed to be applied with four different scaling factors
scalars = np.array([0.2, 0.5, 1, 1.5])
# so what is the array holding a vector for each location displaced by each scaled displacement?
# this doesn't work because the dimensions (10, 3), (15, 3) and (4,) don't match easily:
resulting_locations_incorrect = locations + scalars * displacements
# so instead I introduce singleton dimensions for each array
resulting_locations_correct = (
locations[None, None, :, :] + scalars[:, None, None, None] * displacements[None, :, None, :])
# the resulting shape is (4, 10, 15, 3)
# 4 scalings, 10 displacements, 15 locations, 3 vector elements