public static TripletMatrix CreateCotanLaplacian(TriangleMesh mesh, double lambda = 1.0, double eye = 0.0, bool normalized = false)
        {
            var n = mesh.Vertices.Count;
            int nz = mesh.Vertices.Aggregate(0, (c, x) => c += x.VertexCount());
            var L = new TripletMatrix(n, n, nz + n);

            /// TODO_A2 Task 3:
            /// 3.	Implement cotangents Laplacian matrix (medium)
            double sumF = 0;

            for (int i = 0; i < n; i++)
            {
                foreach(var f in mesh.Vertices.ElementAt(i).Faces)
                {
                    foreach(var he in f.Halfedges)
                    {
                        if(he.FromVertex == mesh.Vertices.ElementAt(i))
                        {
                            sumF += he.Traits.VoronoiRegionArea;
                        }
                    }

                }
                var Half = mesh.Vertices.ElementAt(i).Halfedges;
                for (int j = 0; j < n; j++)
                {
                    if (i == j)
                    {
                        L[i, i] = sumF;
                    }
                    else if (mesh.Vertices.ElementAt(i).Vertices.Contains(mesh.Vertices.ElementAt(j)))
                    {
                        double sumhef = 0;
                        foreach(var edge in mesh.Vertices.ElementAt(i).Edges)
                        {
                            if(edge.Vertex0 == mesh.Vertices.ElementAt(i) && edge.Vertex1 == mesh.Vertices.ElementAt(j) || edge.Vertex0 == mesh.Vertices.ElementAt(j) && edge.Vertex1 == mesh.Vertices.ElementAt(i))
                            {
                                var hef1 = edge.Face0.Halfedges;
                                var hef2 = edge.Face1.Halfedges;
                                foreach(var hef in hef1)
                                {
                                    if(hef.FromVertex == mesh.Vertices.ElementAt(i))
                                    {
                                        sumhef += hef.Traits.VoronoiRegionArea;
                                    }
                                }
                                foreach (var hef in hef2)
                                {
                                    if (hef.FromVertex == mesh.Vertices.ElementAt(i))
                                    {
                                        sumhef += hef.Traits.VoronoiRegionArea;
                                    }
                                }
                            }
                        }
                        if (normalized)
                            L[i, j] = sumhef / 2 / sumF;
                        else
                            L[i, j] = sumhef / 2;
                    }
                    else
                        L[i, j] = 0.0;
                }
            }


            return L;
        }