The is already implemented in the GEOS library (GEOSNearestPoints_r) and is now wrapped in the rgeos package (>= 0.3-21), see gNearestPoints:

R> gDistance(rect, poly.lake)
[1] 33806
R> gNearestPoints(rect, poly.lake)
class       : SpatialPoints 
features    : 2 
extent      : 1036613, 1070419, 825354, 825354  (xmin, xmax, ymin, ymax)
coord. ref. : NA 
R> coordinates(gNearestPoints(rect, poly.lake))
        x      y
1 1036613 825354
2 1070419 825354

The is already implemented in the GEOS library (GEOSNearestPoints_r) and is now wrapped in the rgeos package (>= 0.3-21), see gNearestPoints:

gDistance(rect, poly.lake)
# [1] 33806
gNearestPoints(rect, poly.lake)
# class       : SpatialPoints 
# features    : 2 
# extent      : 1036613, 1070419, 825354, 825354  (xmin, xmax, ymin, ymax)
# coord. ref. : NA 
coordinates(gNearestPoints(rect, poly.lake))
#         x      y
# 1 1036613 825354
# 2 1070419 825354

plot(gUnion(rect, poly.lake))
plot(gNearestPoints(rect, poly.lake), add=TRUE, col="red", pch=4)