This method is based on the original article. First we need to convert list of spring connectivity (spring) into a graph:
g = Graph[Apply[UndirectedEdge, spring, 1], VertexCoordinates -> mass];
Then Kirchhoff (Laplacian/admittance/discrete Laplacian) matrix should be found using KirchhoffMatrix:
m = KirchhoffMatrix[g] // Normal;
We also need to incorporate pinned coordinates into m by replacing the corresponding rows to zero except the diagonal element:
f[i_] := Array[If[#1 == i, 1, 0] &, Last@Dimensions[m]];
m = m /. (m[[#]] -> f[#] & /@ pinned)
Finally the problem is reduced to solve a matrix equation ($m\mathbf{r}=\mathbf{b}$) in $x$ and $y$ direction:
bx = If[MemberQ[pinned, #], mass[[#, 1]], 0] & /@Range[First@Dimensions[m]];
by = If[MemberQ[pinned, #], mass[[#, 2]], 0] & /@Range[First@Dimensions[m]];
xrelaxed = LinearSolve[m, bx];
yrelaxed = LinearSolve[m, by];
And relaxed positions are simply:
relaxedposition = Thread[{xrelaxed, yrelaxed}];