This is not dumb, in fact it is rather interesting. For these three approaches should in theory give equivalent results:
library(sf)
library(spData)
# dramatis personae: shapes & centroids
data("nz")
cent_nz <- st_centroid(st_geometry(nz))
# base solution
plot(nz["Population"])
plot(cent_nz, pch = 4, col = "red", add = T)
# tidyverse solution
library(ggplot2)
ggplot() +
geom_sf(data = nz, aes(fill = Population)) +
geom_sf(data = cent_nz, pch = 4, color = "red")
# tmap solution
library(tmap)
tm_shape(nz) + tm_fill(col = "Population") +
tm_shape(cent_nz) + tm_dots(col = "red", shape = 4)
In practice the result is not equivalent - the base plot "steals" the centroids.
I suspect it has something to do with the projection (st_crs(nz)) and unless you are absolutely positively required to use base plot I would choose either {ggplot2} or {tmap} for your final vizualization.
Most people do, which is why many edge cases in the base plot scenario are left unresolved.