C++ SpMat::setFromTriplets方法代码示例

void Assembler::operKernel(EOExpr<A> oper,const MeshHandler<ORDER>& mesh,
	                     FiniteElement<Integrator, ORDER>& fe, SpMat& OpMat)
{
	std::vector<coeff> triplets;


  	for(auto t=0; t<mesh.num_triangles(); t++)
  	{
		fe.updateElement(mesh.getTriangle(t));

		// Vector of vertices indices (link local to global indexing system)
		std::vector<UInt> identifiers;
		identifiers.resize(ORDER*3);
		for( auto q=0; q<ORDER*3; q++)
			identifiers[q]=mesh.getTriangle(t)[q].id();


		//localM=localMassMatrix(currentelem);
		for(int i = 0; i < 3*ORDER; i++)
		{
			for(int j = 0; j < 3*ORDER; j++)
			{
				Real s=0;

				for(int l = 0;l < Integrator::NNODES; l++)
				{
					s += oper(fe,i,j,l) * fe.getDet() * fe.getAreaReference()* Integrator::WEIGHTS[l];//(*)
				}
			  triplets.push_back(coeff(identifiers[i],identifiers[j],s));
			}
		}

	}

  	UInt nnodes = mesh.num_nodes();
  	OpMat.resize(nnodes, nnodes);
	OpMat.setFromTriplets(triplets.begin(),triplets.end());
	//cout<<"done!"<<endl;;
}

// The Real Core Function doing the actual mesh processing.
bool FilterHarmonicPlugin::applyFilter(QAction * action, MeshDocument & md, RichParameterSet & par, vcg::CallBackPos * cb)
{
	switch(ID(action))
	{
	case FP_SCALAR_HARMONIC_FIELD :
	{
		typedef vcg::GridStaticPtr<CMeshO::VertexType, CMeshO::ScalarType> VertexGrid;

		typedef double                           CoeffScalar; // TODO, when moving the code to a class make it a template (CoeffScalar = double)

		typedef CMeshO::ScalarType               ScalarType;
		typedef CMeshO::CoordType                CoordType;
		typedef CMeshO::VertexType               VertexType;
		typedef CMeshO::FaceType                 FaceType;

		typedef Eigen::Triplet<CoeffScalar>      T;
		typedef Eigen::SparseMatrix<CoeffScalar> SpMat; //sparse matrix type of double


		CMeshO & m = md.mm()->cm;

		vcg::tri::Allocator<CMeshO>::CompactFaceVector(m);
		vcg::tri::Allocator<CMeshO>::CompactVertexVector(m);

		md.mm()->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_VERTMARK);
		vcg::tri::UpdateBounding<CMeshO>::Box(m);
		vcg::tri::UpdateTopology<CMeshO>::FaceFace(m);

		int n  = m.VN();
		int fn = m.FN();

		std::vector<T>             coeffs; // coefficients of the system
		std::map<size_t,CoeffScalar> sums; // row sum of the coefficient

		SpMat laplaceMat; // the system to be solved
		laplaceMat.resize(n, n);

		Log("Generating coefficients.`");
		cb(0, "Generating coefficients...");
		vcg::tri::UpdateFlags<CMeshO>::FaceClearV(m);
		// Iterate over the faces
		for (size_t i = 0; i < m.face.size(); ++i)
		{
			CMeshO::FaceType & f = m.face[i];

			if (f.IsD())
			{
				assert(int(i) == fn);
				break; // TODO FIX the indexing of vertices
			}

			assert(!f.IsV());
			f.SetV();

			// Generate coefficients for each edge
			for (int idx = 0; idx < 3; ++idx)
			{
				CoeffScalar weight;
				WeightInfo res = ComputeWeight<FaceType, CoeffScalar>(f, idx, weight);

				switch (res)
				{
				case EdgeAlreadyVisited : continue;
				case Success            : break;
				case BorderEdge         :
					this->errorMessage = "Mesh not closed, cannot compute harmonic field on mesh containing holes or borders";
					return false;
				default: assert(0);
				}

//				if (weight < 0) weight = 0; // TODO check if negative weight may be an issue

				// Add the weight to the coefficients vector for both the vertices of the considered edge
				size_t v0_idx = vcg::tri::Index(m, f.V0(idx));
				size_t v1_idx = vcg::tri::Index(m, f.V1(idx));

				coeffs.push_back(T(v0_idx, v1_idx, -weight));
				coeffs.push_back(T(v1_idx, v0_idx, -weight));

				// Add the weight to the row sum
				sums[v0_idx] += weight;
				sums[v1_idx] += weight;
			}

			f.SetV();
		}

		// Fill the system matrix
		Log("Fill the system matrix");
		cb(10, "Filling the system matrix...");
		laplaceMat.reserve(coeffs.size());
		for (std::map<size_t,CoeffScalar>::const_iterator it = sums.begin(); it != sums.end(); ++it)
		{
			coeffs.push_back(T(it->first, it->first, it->second));
		}
		laplaceMat.setFromTriplets(coeffs.begin(), coeffs.end());

		// Get the two vertices with value set
		VertexGrid vg;
//.........这里部分代码省略.........

void ADMMCut::CalculateQ(const VX _D, SpMat &Q)
{
	// Construct Hessian Matrix for D-Qp problem
	Q.resize(6 * Ns_, Nd_);
	vector<Eigen::Triplet<double>> Q_list;

	for (int i = 0; i < Ns_; i++)
	{
		int u = ptr_dualgraph_->v_orig_id(i);
		WF_edge *edge = ptr_frame_->GetNeighborEdge(u);

		while (edge != NULL)
		{
			int e_id = ptr_dualgraph_->e_dual_id(edge->ID());
			if (e_id != -1)
			{
				int v = edge->pvert_->ID();
				int j = ptr_dualgraph_->v_dual_id(v);

				MX eKuu = ptr_stiff_->eKv(edge->ID());
				MX eKeu = ptr_stiff_->eKe(edge->ID());
				VX Fe = ptr_stiff_->Fe(edge->ID());
				VX Di(6);
				VX Dj(6);

				if (i < Ns_ && j < Ns_)
				{
					for (int k = 0; k < 6; k++)
					{
						Di[k] = _D[6 * i + k];
						Dj[k] = _D[6 * j + k];
					}
				}
				else
				{
					if (i < Ns_)
					{
						for (int k = 0; k < 6; k++)
						{
							Di[k] = _D[6 * i + k];
							Dj[k] = 0;
						}
					}

					if (j < Ns_)
					{
						for (int k = 0; k < 6; k++)
						{
							Di[k] = 0;
							Dj[k] = _D[6 * j + k];
						}
					}
				}
				VX Gamma = eKuu * Di + eKeu * Dj - Fe;

				for (int k = 0; k < 6; k++)
				{
					Q_list.push_back(Triplet<double>(6 * i + k, e_id, Gamma[k]));
				}
			}

			edge = edge->pnext_;
		}
	}

	Q.setFromTriplets(Q_list.begin(), Q_list.end());
}