As I commented, your problem might be solved by buffering -which means to extend the polygon homogeneously-, dissolve and then use a negative buffer to approximately restore it to its original size; your polygons had remarkably large gaps of around 100 m, usually 0.1 m works for cleaning these kind of polygons.
At the end, change in total area is 0.0007%, which may be negligible.
cities <- read_municipality(code_muni = "all", year = 2010)
# first we buffer
dissolved = cities |>
st_transform(3857) |> st_buffer(100)
# then dissolve
dissolved = dissolved |> group_by(code_state) |> summarise()
# buffer back
dissolved = dissolved |> st_buffer(-100) |>
st_area(dissolved) |> sum()
st_make_valid(cities) |> st_area() |>