R Under development (unstable) (2023-07-09 r84667 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2023 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> suppressPackageStartupMessages(library(sf))
> suppressPackageStartupMessages(library(units))
> 
> x = st_sfc(
+ st_point(c(0,0)),
+ st_point(c(1,0)),
+ st_point(c(2,0)),
+ st_point(c(3,0)),
+ crs = 4326
+ )
> 
> y = st_sfc(
+ st_point(c(0,10)),
+ st_point(c(1,0)),
+ st_point(c(2,0)),
+ st_point(c(3,0)),
+ st_point(c(4,0)),
+ crs = 4326
+ )
> 
> if (suppressPackageStartupMessages(require(sp, quietly = TRUE))) {
+  d.sf = st_distance(x, y)
+  d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
+  units(d.sp) = as_units("km")
+  print(round(d.sf - d.sp, 7))
+ 
+ #summary(unclass(d.sf) - d.sp)
+ 
+  st_crs(x) = st_crs(y) = NA
+  d.sf = st_distance(x, y)
+  d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
+  print(round(d.sf - d.sp, 7))
+ }
Units: [m]
         [,1]      [,2]      [,3]      [,4]      [,5]
[1,] 6107.765 -124.3896 -248.7792 -373.1688 -497.5585
[2,] 6065.138    0.0000 -124.3896 -248.7792 -373.1688
[3,] 5940.569 -124.3896    0.0000 -124.3896 -248.7792
[4,] 5743.252 -248.7792 -124.3896    0.0000 -124.3896
     [,1] [,2] [,3] [,4] [,5]
[1,]    0    0    0    0    0
[2,]    0    0    0    0    0
[3,]    0    0    0    0    0
[4,]    0    0    0    0    0
> 
> # st_length:
> st_crs(y) = 4326
> (z = st_sfc(st_linestring(rbind(c(0,10), c(1,0), c(2,0), c(3,0), c(4,0))), crs = 4326))
Geometry set for 1 feature 
Geometry type: LINESTRING
Dimension:     XY
Bounding box:  xmin: 0 ymin: 0 xmax: 4 ymax: 10
Geodetic CRS:  WGS 84
LINESTRING (0 10, 1 0, 2 0, 3 0, 4 0)
> d = st_distance(y, y)
> round(d, 7)
Units: [m]
        [,1]      [,2]      [,3]      [,4]      [,5]
[1,]       0 1117440.6 1133750.1 1160423.1 1196767.0
[2,] 1117441       0.0  111195.1  222390.2  333585.3
[3,] 1133750  111195.1       0.0  111195.1  222390.2
[4,] 1160423  222390.2  111195.1       0.0  111195.1
[5,] 1196767  333585.3  222390.2  111195.1       0.0
> st_length(z)
1451026 [m]
> st_length(z) - sum(d[1,2], d[2,3], d[3,4], d[4,5])
0 [m]
> 
> # st_line_sample:
> ls = st_sfc(st_linestring(rbind(c(0,0),c(0,1))),
+ 	st_linestring(rbind(c(0,0),c(10,0))))
> # set.seed(135)
> st_line_sample(ls, density = 1)
Geometry set for 2 features 
Geometry type: MULTIPOINT
Dimension:     XY
Bounding box:  xmin: 0 ymin: 0 xmax: 9.5 ymax: 0.5
CRS:           NA
MULTIPOINT ((0 0.5))
MULTIPOINT ((0.5 0), (1.5 0), (2.5 0), (3.5 0),...
> 
> ls = st_sfc(st_linestring(rbind(c(0,0),c(0,1))),
+ 	st_linestring(rbind(c(0,0),c(.1,0))), crs = 4326)
> 
> st_length(ls)
Units: [m]
[1] 111195.10  11119.51
> try(st_line_sample(ls, density = 1/1000))
Error in st_line_sample(ls, density = 1/1000) : 
  st_line_sample for longitude/latitude not supported; use st_segmentize?
> x = st_line_sample(st_transform(ls, 3857), density = 1/1000) # one per km
> 
> proc.time()
   user  system elapsed 
   1.15    0.09    1.25