(* Mathematica code to draw a region, and compute the curvatures, Oliver Knill, July 2, 2009        *)
(* See the article: Oliver Knill: A discrete Gauss-Bonnet type theorem, Submitted February 11, 2010 *)

GG={{0,1},{0,0},{1,0},{1,1},{1,2},{0,3},{0,4},{2,2}};                          (*  the set of interior points          *)
M=10;                                                                          (*  size of region to compute           *)
R=5; Rx=2; Ry=2;                                                               (*  center (Rx,Ry) and size of set      *) 

v={1,0}; w = {1/2,Sqrt[3]/2};  r[{i_,j_}]:=v*i+w*j;
G0=Partition[Flatten[Table[{i,j},{i,-1,1},{j,-1,1}]],2];G1=Table[G0[[k]],{k,2,8}]; G0=G1;
Cluster[A_]:=Union[Partition[Flatten[Table[G0[[j]]+A[[i]],{j,Length[G0]},{i,Length[A]}]],2]];

DR[G_]:=Module[{},
InRegion[x_]:=MemberQ[Map[r,G],x];
P={}; Do[ x=r[{i,j}]; If[InRegion[x],P=Append[P,Point[x]]],{i,-M,M},{j,-M,M}]; SP = {RGBColor[1,0.5,0],P};
L={};Do[x=r[{i,j}];y=r[{i+s,j+t}];
   If[({s,t}=={1,0}||{s,t}=={0,1}||{s,t}=={-1,1}||{s,t}=={1,-1}||{s,t}=={-1,0}||{s,t}=={0,-1})
   &&InRegion[x] && InRegion[y], L=Append[L,Line[{x,y}]]],{i,-M,M},{j,-M,M},{s,-1,1},{t,-1,1}];
SL = {RGBColor[0.2,0.2,0.2],Thickness[0.001],L};
myline[o_]:=Table[If[(o[[Mod[i-1,Length[o]]+1]]-o[[Mod[i,Length[o]]+1]]).(o[[Mod[i-1,Length[o]]+1]]-o[[Mod[i,Length[o]]+1]])<1.1,
     Line[{o[[Mod[i-1,Length[o]]+1]],o[[Mod[i,Length[o]]+1]]}],{}],{i,0,Length[o]}];
Sphere0[x_]:=If[InRegion[x],{RGBColor[1,1,0],PointSize[0.02],Point[x]},{}];
Sphere1[x_]:=Module[{},U={x+v,x+w,x-v+w,x-v,x-w,x+v-w};
   V={}; Do[If[InRegion[x] && InRegion[U[[i]]],V=Append[V,U[[i]]]],{i,Length[U]}];
   If[Length[V]>0,If[(Last[V]-First[V]).(Last[V]-First[V])<1.1,V=Append[V,First[V]]]];
   If[V=={},{},{{RGBColor[1,0,0],Thickness[0.01],myline[V]},{RGBColor[0,1,0],PointSize[0.005],Map[Point,V]}}]];
Sphere2[x_]:=Module[{},U={x+2 v,x+v+w,x+2w,x+2w-v,x+2w-2v,x-2v+w,x-2v,x-v-w,x-2w,x-2w+v,x-2w+2v,x+2 v-w};
   V={}; Do[If[InRegion[x] && InRegion[U[[i]]],V=Append[V,U[[i]]]],{i,Length[U]}];
   If[Length[V]>0,If[(Last[V]-First[V]).(Last[V]-First[V])<1.1,V=Append[V,First[V]]]];
   If[V=={},{},{{RGBColor[0.5,1,0],Thickness[0.01],myline[V]},{RGBColor[0,0,1],PointSize[0.005],Map[Point,V]}}]];
Neighborhood[x_]:={Sphere0[x],Sphere1[x],Sphere2[x]};
Curvature[x_]:=Module[{},c1=0; c2=0;
  U1={x+v,x+w,x-v+w,x-v,x-w,x+v-w};   (* cc1 = number of neighbors, c1 = length of |S_1| *)
  V1={}; Do[If[InRegion[x] && InRegion[U1[[i]]],V1=Append[V1,U1[[i]]]],{i,Length[U1]}]; cc1=Length[V1];
  If[cc1==2,c1=1,If[cc1>0,If[0<((Last[V1]-First[V1]).(Last[V1]-First[V1]))<1.1,V1=Append[V1,First[V1]]]];
    Do[vv=V1[[i+1]]-V1[[i]]; If[vv.vv<1.1,c1++],{i,Length[V1]-1}]; ];
  U2={x+2v,x+v+w,x+2w,x+2w-v,x+2w-2v,x-2v+w,x-2v,x-v-w,x-2w,x-2w+v,x-2w+2v,x+2v-w};
  V2={}; Do[If[InRegion[x] && InRegion[U2[[i]]],V2=Append[V2,U2[[i]]]],{i,Length[U2]}];
  If[Length[V2]==2,c2=1,
    If[Length[V2]>0,If[0<((Last[V2]-First[V2]).(Last[V2]-First[V2]))<1.1,V2=Append[V2,First[V2]]]];
    Do[vv=V2[[i+1]]-V2[[i]]; If[vv.vv<1.1,c2++],{i,Length[V2]-1}];
  ]; {cc1,c1,c2,2*c1-c2}];   (*  { neighbors, |S1|, |S2|, K } *)
Dimension[x_]:=Module[{}, {cc1,c1,c2,c}=Curvature[x]; dim=0; If[cc1==6,dim=2]; If[6>cc1>0,dim=1]; dim];
BoundaryPoint[x_]:=Module[{}, bnd=False; If[Dimension[x]==1 && (
      Curvature[x+v][[2]]==6 || Curvature[x+w][[2]]==6 || Curvature[x-v+w][[2]]==6 ||
      Curvature[x-v][[2]]==6 || Curvature[x-w][[2]]==6 || Curvature[x+v-w][[2]]==6), bnd=True]; bnd];
TK = 0; K=0; T={}; PL={};  (* K boundary curvature, T text, PL  Point labels *)
Do[x=r[{i,j}]; If[InRegion[x],  dim=Dimension[x]; bnd=BoundaryPoint[x];
   If[bnd,{red,green,blue}={0,1,0},{red,green,blue}={1,1,0}]; If[dim==2,{red,green,blue}={1,0,0}];
   TK=TK+Curvature[x][[4]];
   If[dim==1,K=K+Curvature[x][[4]]];
   PL=Append[PL,{RGBColor[red,green,blue],Thickness[0.01],Circle[x,(dim+1)*0.1]}];
   T=Append[T,{FontSize->20,Text[Curvature[x][[4]],x+{0.2+0.1*dim,0.1+0.1*dim}]}]],{i,-M,M},{j,-M,M}]; (* Print[K]; *)
KK1= {FontSize->50,Text[K,{R-0.2,R-0.2}]}; curvature=K; KK2= {FontSize->20,Text[TK,{R-0.4,-R+0.4}]};
S01= Neighborhood[r[{3,1}]];
S=Show[Graphics[{SP,SL,T,PL,KK1}],PlotRange->{{Rx-R,Rx+R},{Ry-R,Ry+R}}]
]

DR[Cluster[GG]]
