+
+ CONVINTERSECTOR_TEMPLATE
+ double CONVEX_INTERSECTOR_::intersectGeometryGeneral(const std::vector<double>& targetCoords,
+ const std::vector<double>& sourceCoords)
+ {
+ double result = 0;
+ int nbOfNodesS=sourceCoords.size()/SPACEDIM;
+ int nbOfNodesT=targetCoords.size()/SPACEDIM;
+ /*** Compute the intersection area ***/
+ INTERP_KERNEL::PolygonAlgorithms<SPACEDIM> P(_epsilon, PlanarIntersector<MyMeshType,MyMatrix>::_precision);
+ std::deque<double> inter = P.intersectConvexPolygons(&targetCoords[0], &sourceCoords[0],
+ nbOfNodesT, nbOfNodesS);
+ double area[SPACEDIM];
+ int nb_inter =((int)inter.size())/SPACEDIM;
+ for(int i = 1; i<nb_inter-1; i++)
+ {
+ INTERP_KERNEL::crossprod<SPACEDIM>(&inter[0],&inter[SPACEDIM*i],&inter[SPACEDIM*(i+1)],area);
+ result +=0.5*norm<SPACEDIM>(area);
+ }
+ return result;
+ }
+
+ //================================================================================
+ /*!
+ * \brief Intersect a triangle and a polygon for P1P0 barycentric algorithm
+ * \param targetCell - list of coordinates of target polygon in full interlace
+ * \param targetCellQuadratic - specifies if target polygon is quadratic or not
+ * \param sourceTria - list of coordinates of source triangle
+ * \param res - coefficients a,b and c associated to nodes of sourceTria
+ */
+ //================================================================================
+
+ CONVINTERSECTOR_TEMPLATE
+ double CONVEX_INTERSECTOR_::intersectGeoBary(const std::vector<double>& targetCell,
+ bool targetCellQuadratic,
+ const double * sourceTria,
+ std::vector<double>& res)
+ {
+ double area = 0;
+ double barycenter[SPACEDIM] = {0., 0.};
+ int nbOfNodesT=targetCell.size()/SPACEDIM;
+
+ /*** Compute the intersection area ***/
+ INTERP_KERNEL::PolygonAlgorithms<SPACEDIM> P(_epsilon, PlanarIntersector<MyMeshType,MyMatrix>::_precision);
+ std::deque<double> inter = P.intersectConvexPolygons(sourceTria, &targetCell[0], 3, nbOfNodesT);
+ double cross[SPACEDIM];
+ int nb_inter =((int)inter.size())/SPACEDIM;
+ for(int i = 1; i<nb_inter-1; i++)
+ {
+ INTERP_KERNEL::crossprod<SPACEDIM>(&inter[0],&inter[SPACEDIM*i],&inter[SPACEDIM*(i+1)],cross);
+ area += 0.5*norm<SPACEDIM>(cross);
+ barycenter[0] += inter[SPACEDIM*i];
+ barycenter[1] += inter[SPACEDIM*i+1];
+ }
+ if ( area > std::numeric_limits<double>::min() )
+ {
+ barycenter[0] = ( barycenter[0] + inter[0] + inter[SPACEDIM*(nb_inter-1)] ) / nb_inter;
+ barycenter[1] = ( barycenter[1] + inter[1] + inter[SPACEDIM*(nb_inter-1)+1]) / nb_inter;
+ res.resize(3);
+ barycentric_coords<2>( sourceTria, &barycenter[0], &res[0]);
+ res[0] *= area;
+ res[1] *= area;
+ res[2] *= area;
+ }
+ else
+ {
+ area = 0;
+ }
+ return area;
+ }