C++ SMESH_subMesh::GetComputeError方法代码示例

bool NETGENPlugin_NETGEN_3D::Evaluate(SMESH_Mesh& aMesh,
                                      const TopoDS_Shape& aShape,
                                      MapShapeNbElems& aResMap)
{
    int nbtri = 0, nbqua = 0;
    double fullArea = 0.0;
    for (TopExp_Explorer exp(aShape, TopAbs_FACE); exp.More(); exp.Next()) {
        TopoDS_Face F = TopoDS::Face( exp.Current() );
        SMESH_subMesh *sm = aMesh.GetSubMesh(F);
        MapShapeNbElemsItr anIt = aResMap.find(sm);
        if( anIt==aResMap.end() ) {
            SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
            smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
            return false;
        }
        std::vector<int> aVec = (*anIt).second;
        nbtri += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
        nbqua += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
        GProp_GProps G;
        BRepGProp::SurfaceProperties(F,G);
        double anArea = G.Mass();
        fullArea += anArea;
    }

    // collect info from edges
    int nb0d_e = 0, nb1d_e = 0;
    bool IsQuadratic = false;
    bool IsFirst = true;
    TopTools_MapOfShape tmpMap;
    for (TopExp_Explorer exp(aShape, TopAbs_EDGE); exp.More(); exp.Next()) {
        TopoDS_Edge E = TopoDS::Edge(exp.Current());
        if( tmpMap.Contains(E) )
            continue;
        tmpMap.Add(E);
        SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
        MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
        if( anIt==aResMap.end() ) {
            SMESH_ComputeErrorPtr& smError = aSubMesh->GetComputeError();
            smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,
                                                  "Submesh can not be evaluated",this));
            return false;
        }
        std::vector<int> aVec = (*anIt).second;
        nb0d_e += aVec[SMDSEntity_Node];
        nb1d_e += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
        if(IsFirst) {
            IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
            IsFirst = false;
        }
    }
    tmpMap.Clear();

    double ELen_face = sqrt(2.* ( fullArea/(nbtri+nbqua*2) ) / sqrt(3.0) );
    double ELen_vol = pow( 72, 1/6. ) * pow( _maxElementVolume, 1/3. );
    double ELen = Min(ELen_vol,ELen_face*2);

    GProp_GProps G;
    BRepGProp::VolumeProperties(aShape,G);
    double aVolume = G.Mass();
    double tetrVol = 0.1179*ELen*ELen*ELen;
    double CoeffQuality = 0.9;
    int nbVols = (int)aVolume/tetrVol/CoeffQuality;
    int nb1d_f = (nbtri*3 + nbqua*4 - nb1d_e) / 2;
    int nb1d_in = (int) ( nbVols*6 - nb1d_e - nb1d_f ) / 5;
    std::vector<int> aVec(SMDSEntity_Last);
    for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
    if( IsQuadratic ) {
        aVec[SMDSEntity_Node] = nb1d_in/6 + 1 + nb1d_in;
        aVec[SMDSEntity_Quad_Tetra] = nbVols - nbqua*2;
        aVec[SMDSEntity_Quad_Pyramid] = nbqua;
    }
    else {
        aVec[SMDSEntity_Node] = nb1d_in/6 + 1;
        aVec[SMDSEntity_Tetra] = nbVols - nbqua*2;
        aVec[SMDSEntity_Pyramid] = nbqua;
    }
    SMESH_subMesh *sm = aMesh.GetSubMesh(aShape);
    aResMap.insert(std::make_pair(sm,aVec));

    return true;
}

//.........这里部分代码省略.........

    case TopAbs_FACE: { // FACE
      // ----------------------
      const TopoDS_Face& geomFace  = TopoDS::Face( sm->GetSubShape() );
      helper.SetSubShape( geomFace );

      // Find solids the geomFace bounds
      int solidID1 = 0, solidID2 = 0;
      const TopTools_ListOfShape& ancestors = _mesh->GetAncestors( geomFace );
      TopTools_ListIteratorOfListOfShape ancestorIt ( ancestors );
      for ( ; ancestorIt.More(); ancestorIt.Next() )
      {
        const TopoDS_Shape & solid = ancestorIt.Value();
        if ( solid.ShapeType() == TopAbs_SOLID  ) {
          int id = occgeom.somap.FindIndex ( solid );
          if ( solidID1 && id != solidID1 ) solidID2 = id;
          else                              solidID1 = id;
        }
      }
      faceID++;
      _faceDescriptors[ faceID ].first  = solidID1;
      _faceDescriptors[ faceID ].second = solidID2;

      // Orient the face correctly in solidID1 (issue 0020206)
      bool reverse = false;
      if ( solidID1 ) {
        TopoDS_Shape solid = occgeom.somap( solidID1 );
        for ( TopExp_Explorer f( solid, TopAbs_FACE ); f.More(); f.Next() ) {
          if ( geomFace.IsSame( f.Current() )) {
            reverse = SMESH_Algo::IsReversedSubMesh( TopoDS::Face( f.Current()), helper.GetMeshDS() );
            break;
          }
        }
      }

      // Add surface elements
      SMDS_ElemIteratorPtr faces = smDS->GetElements();
      while ( faces->more() ) {

        const SMDS_MeshElement* f = faces->next();
        if ( f->NbNodes() % 3 != 0 ) { // not triangle
          for ( ancestorIt.Initialize(ancestors); ancestorIt.More(); ancestorIt.Next() )
            if ( ancestorIt.Value().ShapeType() == TopAbs_SOLID  ) {
              sm = _mesh->GetSubMesh( ancestorIt.Value() );
              break;
            }
          SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
          smError.reset( new SMESH_ComputeError(COMPERR_BAD_INPUT_MESH,"Not triangle submesh"));
          smError->myBadElements.push_back( f );
          return false;
        }

        netgen::Element2d tri(3);
        tri.SetIndex ( faceID );

        for ( int i = 0; i < 3; ++i ) {
          const SMDS_MeshNode* node = f->GetNode( i ), * inFaceNode=0;
          if ( helper.IsSeamShape( node->GetPosition()->GetShapeId() ))
            if ( helper.IsSeamShape( f->GetNodeWrap( i+1 )->GetPosition()->GetShapeId() ))
              inFaceNode = f->GetNodeWrap( i-1 );
            else 
              inFaceNode = f->GetNodeWrap( i+1 );

          gp_XY uv = helper.GetNodeUV( geomFace, node, inFaceNode );
          if ( reverse ) {
            tri.GeomInfoPi(3-i).u = uv.X();
            tri.GeomInfoPi(3-i).v = uv.Y();
            tri.PNum      (3-i) = ngNodeId( node, ngMesh, nodeNgIdMap );
          } else {
            tri.GeomInfoPi(i+1).u = uv.X();
            tri.GeomInfoPi(i+1).v = uv.Y();
            tri.PNum      (i+1) = ngNodeId( node, ngMesh, nodeNgIdMap );
          }
        }

        ngMesh.AddSurfaceElement (tri);

      }
      break;
    } //

    case TopAbs_VERTEX: { // VERTEX
      // --------------------------
      SMDS_NodeIteratorPtr nodeIt = smDS->GetNodes();
      if ( nodeIt->more() )
        ngNodeId( nodeIt->next(), ngMesh, nodeNgIdMap );
      break;
    }
    default:;
    } // switch
  } // loop on submeshes

  // fill nodeVec
  nodeVec.resize( ngMesh.GetNP() + 1 );
  TNode2IdMap::iterator node_NgId, nodeNgIdEnd = nodeNgIdMap.end();
  for ( node_NgId = nodeNgIdMap.begin(); node_NgId != nodeNgIdEnd; ++node_NgId)
    nodeVec[ node_NgId->second ] = (SMDS_MeshNode*) node_NgId->first;

  return true;
}

//.........这里部分代码省略.........
      // ------------------------------------------------
      // sort list of sub-meshes according to mesh order
      // ------------------------------------------------
      smVec.assign( smWithAlgoSupportingSubmeshes[ aShapeDim ].begin(),
                    smWithAlgoSupportingSubmeshes[ aShapeDim ].end() );
      aMesh.SortByMeshOrder( smVec );

      // ------------------------------------------------------------
      // compute sub-meshes with local uni-dimensional algos under
      // sub-meshes with all-dimensional algos
      // ------------------------------------------------------------
      // start from lower shapes
      for ( size_t i = 0; i < smVec.size(); ++i )
      {
        sm = smVec[i];

        // get a shape the algo is assigned to
        if ( !GetAlgo( sm, & algoShape ))
          continue; // strange...

        // look for more local algos
        smIt = sm->getDependsOnIterator(!includeSelf, !complexShapeFirst);
        while ( smIt->more() )
        {
          SMESH_subMesh* smToCompute = smIt->next();

          const TopoDS_Shape& aSubShape = smToCompute->GetSubShape();
          const int aShapeDim = GetShapeDim( aSubShape );
          //if ( aSubShape.ShapeType() == TopAbs_VERTEX ) continue;
          if ( aShapeDim < 1 ) continue;

          // check for preview dimension limitations
          if ( aShapesId && GetShapeDim( aSubShape.ShapeType() ) > (int)aDim )
            continue;

          SMESH_HypoFilter filter( SMESH_HypoFilter::IsAlgo() );
          filter
            .And( SMESH_HypoFilter::IsApplicableTo( aSubShape ))
            .And( SMESH_HypoFilter::IsMoreLocalThan( algoShape, aMesh ));

          if ( SMESH_Algo* subAlgo = (SMESH_Algo*) aMesh.GetHypothesis( smToCompute, filter, true))
          {
            if ( ! subAlgo->NeedDiscreteBoundary() ) continue;
            SMESH_Hypothesis::Hypothesis_Status status;
            if ( subAlgo->CheckHypothesis( aMesh, aSubShape, status ))
              // mesh a lower smToCompute starting from vertices
              Compute( aMesh, aSubShape, aShapeOnly, /*anUpward=*/true, aDim, aShapesId );
          }
        }
      }
      // --------------------------------
      // apply the all-dimensional algos
      // --------------------------------
      for ( size_t i = 0; i < smVec.size(); ++i )
      {
        sm = smVec[i];
        if ( sm->GetComputeState() == SMESH_subMesh::READY_TO_COMPUTE)
        {
          const TopAbs_ShapeEnum shapeType = sm->GetSubShape().ShapeType();
          // check for preview dimension limitations
          if ( aShapesId && GetShapeDim( shapeType ) > (int)aDim )
            continue;

          if (_compute_canceled)
            return false;
          setCurrentSubMesh( sm );
          sm->ComputeStateEngine( computeEvent );
          setCurrentSubMesh( NULL );
          if ( aShapesId )
            aShapesId->insert( sm->GetId() );
        }
      }
    } // loop on shape dimensions

    // -----------------------------------------------
    // mesh the rest sub-shapes starting from vertices
    // -----------------------------------------------
    ret = Compute( aMesh, aShape, aShapeOnly, /*anUpward=*/true, aDim, aShapesId );
  }

  MESSAGE( "VSR - SMESH_Gen::Compute() finished, OK = " << ret);
  MEMOSTAT;

  SMESHDS_Mesh *myMesh = aMesh.GetMeshDS();
  MESSAGE("*** compactMesh after compute");
  myMesh->compactMesh();

  // fix quadratic mesh by bending iternal links near concave boundary
  if ( aShape.IsSame( aMesh.GetShapeToMesh() ) &&
       !aShapesId && // not preview
       ret ) // everything is OK
  {
    SMESH_MesherHelper aHelper( aMesh );
    if ( aHelper.IsQuadraticMesh() != SMESH_MesherHelper::LINEAR )
    {
      aHelper.FixQuadraticElements( sm->GetComputeError() );
    }
  }
  return ret;
}

bool StdMeshers_RadialPrism_3D::Evaluate(SMESH_Mesh& aMesh,
                                         const TopoDS_Shape& aShape,
                                         MapShapeNbElems& aResMap)
{
  // get 2 shells
  TopoDS_Solid solid = TopoDS::Solid( aShape );
  TopoDS_Shell outerShell = BRepClass3d::OuterShell( solid );
  TopoDS_Shape innerShell;
  int nbShells = 0;
  for ( TopoDS_Iterator It (solid); It.More(); It.Next(), ++nbShells )
    if ( !outerShell.IsSame( It.Value() ))
      innerShell = It.Value();
  if ( nbShells != 2 ) {
    std::vector<int> aResVec(SMDSEntity_Last);
    for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
    SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
    aResMap.insert(std::make_pair(sm,aResVec));
    SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
    smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
    return false;
  }

  // Associate sub-shapes of the shells
  ProjectionUtils::TShapeShapeMap shape2ShapeMap;
  if ( !ProjectionUtils::FindSubShapeAssociation( outerShell, &aMesh,
                                                  innerShell, &aMesh,
                                                  shape2ShapeMap) ) {
    std::vector<int> aResVec(SMDSEntity_Last);
    for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
    SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
    aResMap.insert(std::make_pair(sm,aResVec));
    SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
    smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
    return false;
  }

  // get info for outer shell
  int nb0d_Out=0, nb2d_3_Out=0, nb2d_4_Out=0;
  //TopTools_SequenceOfShape FacesOut;
  for (TopExp_Explorer exp(outerShell, TopAbs_FACE); exp.More(); exp.Next()) {
    //FacesOut.Append(exp.Current());
    SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
    MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
    std::vector<int> aVec = (*anIt).second;
    nb0d_Out += aVec[SMDSEntity_Node];
    nb2d_3_Out += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
    nb2d_4_Out += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
  }
  int nb1d_Out = 0;
  TopTools_MapOfShape tmpMap;
  for (TopExp_Explorer exp(outerShell, TopAbs_EDGE); exp.More(); exp.Next()) {
    if( tmpMap.Contains( exp.Current() ) )
      continue;
    tmpMap.Add( exp.Current() );
    SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
    MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
    std::vector<int> aVec = (*anIt).second;
    nb0d_Out += aVec[SMDSEntity_Node];
    nb1d_Out += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
  }
  tmpMap.Clear();
  for (TopExp_Explorer exp(outerShell, TopAbs_VERTEX); exp.More(); exp.Next()) {
    if( tmpMap.Contains( exp.Current() ) )
      continue;
    tmpMap.Add( exp.Current() );
    nb0d_Out++;
  }

  // get info for inner shell
  int nb0d_In=0, nb2d_3_In=0, nb2d_4_In=0;
  //TopTools_SequenceOfShape FacesIn;
  for (TopExp_Explorer exp(innerShell, TopAbs_FACE); exp.More(); exp.Next()) {
    //FacesIn.Append(exp.Current());
    SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
    MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
    std::vector<int> aVec = (*anIt).second;
    nb0d_In += aVec[SMDSEntity_Node];
    nb2d_3_In += Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
    nb2d_4_In += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
  }
  int nb1d_In = 0;
  tmpMap.Clear();
  bool IsQuadratic = false;
  bool IsFirst = true;
  for (TopExp_Explorer exp(innerShell, TopAbs_EDGE); exp.More(); exp.Next()) {
    if( tmpMap.Contains( exp.Current() ) )
      continue;
    tmpMap.Add( exp.Current() );
    SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
    MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
    std::vector<int> aVec = (*anIt).second;
    nb0d_In += aVec[SMDSEntity_Node];
    nb1d_In += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
    if(IsFirst) {
      IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
      IsFirst = false;
    }
  }
  tmpMap.Clear();
  for (TopExp_Explorer exp(innerShell, TopAbs_VERTEX); exp.More(); exp.Next()) {
//.........这里部分代码省略.........

//.........这里部分代码省略.........
    LinEdge2 = E3;
    double fp = f1;
    double lp = l1;
    if( aCirc.IsNull() ) {
      aCirc = Handle(Geom_Circle)::DownCast(C2);
      CircEdge = E2;
      LinEdge1 = E3;
      LinEdge2 = E1;
      fp = f2;
      lp = l2;
      aLine1 = Handle(Geom_Line)::DownCast(C3);
      aLine2 = Handle(Geom_Line)::DownCast(C1);
      if( aCirc.IsNull() ) {
        aCirc = Handle(Geom_Circle)::DownCast(C3);
        CircEdge = E3;
        LinEdge1 = E1;
        LinEdge2 = E2;
        fp = f3;
        lp = l3;
        aLine1 = Handle(Geom_Line)::DownCast(C1);
        aLine2 = Handle(Geom_Line)::DownCast(C2);
      }
    }
    bool ok = !aCirc.IsNull() && !aLine1.IsNull() && !aLine1.IsNull();
    SMESH_subMesh* sm = aMesh.GetSubMesh(LinEdge1);
    MapShapeNbElemsItr anIt = aResMap.find(sm);
    if( anIt!=aResMap.end() ) {
      ok = false;
    }
    sm = aMesh.GetSubMesh(LinEdge2);
    anIt = aResMap.find(sm);
    if( anIt!=aResMap.end() ) {
      ok = false;
    }
    if(ok) {
      ok = _gen->Evaluate( aMesh, CircEdge, aResMap );
    }
    if(ok) {
      SMESH_subMesh * sm = aMesh.GetSubMesh(CircEdge);
      MapShapeNbElemsItr anIt = aResMap.find(sm);
      vector<int> aVec = (*anIt).second;
      isQuadratic = aVec[SMDSEntity_Quad_Edge]>aVec[SMDSEntity_Edge];
      if(isQuadratic) {
        // main nodes
        nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
        // radial medium nodes
        nb0d += aVec[SMDSEntity_Node] * (myLayerPositions.size()+1);
        // other medium nodes
        nb0d += (aVec[SMDSEntity_Node]+1) * myLayerPositions.size();
      }
      else {
        nb0d = aVec[SMDSEntity_Node] * myLayerPositions.size();
      }
      nb2d_tria = aVec[SMDSEntity_Node] + 1;
      nb2d_quad = nb2d_tria * myLayerPositions.size();
      // add evaluation for edges
      vector<int> aResVec(SMDSEntity_Last);
      for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
      if(isQuadratic) {
        aResVec[SMDSEntity_Node] = 2*myLayerPositions.size() + 1;
        aResVec[SMDSEntity_Quad_Edge] = myLayerPositions.size() + 1;
      }
      else {
        aResVec[SMDSEntity_Node] = myLayerPositions.size();
        aResVec[SMDSEntity_Edge] = myLayerPositions.size() + 1;
      }
      sm = aMesh.GetSubMesh(LinEdge1);
      aResMap.insert(make_pair(sm,aResVec));
      sm = aMesh.GetSubMesh(LinEdge2);
      aResMap.insert(make_pair(sm,aResVec));
    }
  }

  vector<int> aResVec(SMDSEntity_Last);
  for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
  SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);

  //cout<<"nb0d = "<<nb0d<<"   nb2d_tria = "<<nb2d_tria<<"   nb2d_quad = "<<nb2d_quad<<endl;
  if(nb0d>0) {
    aResVec[0] = nb0d;
    if(isQuadratic) {
      aResVec[SMDSEntity_Quad_Triangle] = nb2d_tria;
      aResVec[SMDSEntity_Quad_Quadrangle] = nb2d_quad;
    }
    else {
      aResVec[SMDSEntity_Triangle] = nb2d_tria;
      aResVec[SMDSEntity_Quadrangle] = nb2d_quad;
    }
    aResMap.insert(make_pair(sm,aResVec));
    return true;
  }

  // invalid case
  aResMap.insert(make_pair(sm,aResVec));
  SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
  smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,
                                        "Submesh can not be evaluated",this));
  return false;

}

bool NETGENPlugin_NETGEN_2D_ONLY::Evaluate(SMESH_Mesh& aMesh,
        const TopoDS_Shape& aShape,
        MapShapeNbElems& aResMap)
{
    TopoDS_Face F = TopoDS::Face(aShape);
    if(F.IsNull())
        return false;

    // collect info from edges
    int nb0d = 0, nb1d = 0;
    bool IsQuadratic = false;
    bool IsFirst = true;
    double fullLen = 0.0;
    TopTools_MapOfShape tmpMap;
    for (TopExp_Explorer exp(F, TopAbs_EDGE); exp.More(); exp.Next()) {
        TopoDS_Edge E = TopoDS::Edge(exp.Current());
        if( tmpMap.Contains(E) )
            continue;
        tmpMap.Add(E);
        SMESH_subMesh *aSubMesh = aMesh.GetSubMesh(exp.Current());
        MapShapeNbElemsItr anIt = aResMap.find(aSubMesh);
        if( anIt==aResMap.end() ) {
            SMESH_subMesh *sm = aMesh.GetSubMesh(F);
            SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
            smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
            return false;
        }
        std::vector<int> aVec = (*anIt).second;
        nb0d += aVec[SMDSEntity_Node];
        nb1d += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
        double aLen = SMESH_Algo::EdgeLength(E);
        fullLen += aLen;
        if(IsFirst) {
            IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
            IsFirst = false;
        }
    }
    tmpMap.Clear();

    // compute edge length
    double ELen = 0;
    if (_hypLengthFromEdges || (!_hypLengthFromEdges && !_hypMaxElementArea)) {
        if ( nb1d > 0 )
            ELen = fullLen / nb1d;
    }
    if ( _hypMaxElementArea ) {
        double maxArea = _hypMaxElementArea->GetMaxArea();
        ELen = sqrt(2. * maxArea/sqrt(3.0));
    }
    GProp_GProps G;
    BRepGProp::SurfaceProperties(F,G);
    double anArea = G.Mass();

    const int hugeNb = numeric_limits<int>::max()/10;
    if ( anArea / hugeNb > ELen*ELen )
    {
        SMESH_subMesh *sm = aMesh.GetSubMesh(F);
        SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
        smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated.\nToo small element length",this));
        return false;
    }
    int nbFaces = (int) ( anArea / ( ELen*ELen*sqrt(3.) / 4 ) );
    int nbNodes = (int) ( ( nbFaces*3 - (nb1d-1)*2 ) / 6 + 1 );
    std::vector<int> aVec(SMDSEntity_Last);
    for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aVec[i]=0;
    if( IsQuadratic ) {
        aVec[SMDSEntity_Node] = nbNodes;
        aVec[SMDSEntity_Quad_Triangle] = nbFaces;
    }
    else {
        aVec[SMDSEntity_Node] = nbNodes;
        aVec[SMDSEntity_Triangle] = nbFaces;
    }
    SMESH_subMesh *sm = aMesh.GetSubMesh(F);
    aResMap.insert(std::make_pair(sm,aVec));

    return true;
}

bool StdMeshers_Hexa_3D::Evaluate(SMESH_Mesh & aMesh,
                                  const TopoDS_Shape & aShape,
                                  MapShapeNbElems& aResMap)
{
  vector < SMESH_subMesh * >meshFaces;
  TopTools_SequenceOfShape aFaces;
  for (TopExp_Explorer exp(aShape, TopAbs_FACE); exp.More(); exp.Next()) {
    aFaces.Append(exp.Current());
    SMESH_subMesh *aSubMesh = aMesh.GetSubMeshContaining(exp.Current());
    ASSERT(aSubMesh);
    meshFaces.push_back(aSubMesh);
  }
  if (meshFaces.size() != 6) {
    //return error(COMPERR_BAD_SHAPE, TComm(meshFaces.size())<<" instead of 6 faces in a block");
    static StdMeshers_CompositeHexa_3D compositeHexa(-10, 0, aMesh.GetGen());
    return compositeHexa.Evaluate(aMesh, aShape, aResMap);
  }
  
  int i = 0;
  for(; i<6; i++) {
    //TopoDS_Shape aFace = meshFaces[i]->GetSubShape();
    TopoDS_Shape aFace = aFaces.Value(i+1);
    SMESH_Algo *algo = _gen->GetAlgo(aMesh, aFace);
    if( !algo ) {
      std::vector<int> aResVec(SMDSEntity_Last);
      for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
      SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
      aResMap.insert(std::make_pair(sm,aResVec));
      SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
      smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Submesh can not be evaluated",this));
      return false;
    }
    string algoName = algo->GetName();
    bool isAllQuad = false;
    if (algoName == "Quadrangle_2D") {
      MapShapeNbElemsItr anIt = aResMap.find(meshFaces[i]);
      if( anIt == aResMap.end() ) continue;
      std::vector<int> aVec = (*anIt).second;
      int nbtri = Max(aVec[SMDSEntity_Triangle],aVec[SMDSEntity_Quad_Triangle]);
      if( nbtri == 0 )
        isAllQuad = true;
    }
    if ( ! isAllQuad ) {
      return EvaluatePentahedralMesh(aMesh, aShape, aResMap);
    }
  }
  
  // find number of 1d elems for 1 face
  int nb1d = 0;
  TopTools_MapOfShape Edges1;
  bool IsQuadratic = false;
  bool IsFirst = true;
  for (TopExp_Explorer exp(aFaces.Value(1), TopAbs_EDGE); exp.More(); exp.Next()) {
    Edges1.Add(exp.Current());
    SMESH_subMesh *sm = aMesh.GetSubMesh(exp.Current());
    if( sm ) {
      MapShapeNbElemsItr anIt = aResMap.find(sm);
      if( anIt == aResMap.end() ) continue;
      std::vector<int> aVec = (*anIt).second;
      nb1d += Max(aVec[SMDSEntity_Edge],aVec[SMDSEntity_Quad_Edge]);
      if(IsFirst) {
        IsQuadratic = (aVec[SMDSEntity_Quad_Edge] > aVec[SMDSEntity_Edge]);
        IsFirst = false;
      }
    }
  }
  // find face opposite to 1 face
  int OppNum = 0;
  for(i=2; i<=6; i++) {
    bool IsOpposite = true;
    for(TopExp_Explorer exp(aFaces.Value(i), TopAbs_EDGE); exp.More(); exp.Next()) {
      if( Edges1.Contains(exp.Current()) ) {
        IsOpposite = false;
        break;
      }
    }
    if(IsOpposite) {
      OppNum = i;
      break;
    }
  }
  // find number of 2d elems on side faces
  int nb2d = 0;
  for(i=2; i<=6; i++) {
    if( i == OppNum ) continue;
    MapShapeNbElemsItr anIt = aResMap.find( meshFaces[i-1] );
    if( anIt == aResMap.end() ) continue;
    std::vector<int> aVec = (*anIt).second;
    nb2d += Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
  }
  
  MapShapeNbElemsItr anIt = aResMap.find( meshFaces[0] );
  std::vector<int> aVec = (*anIt).second;
  int nb2d_face0 = Max(aVec[SMDSEntity_Quadrangle],aVec[SMDSEntity_Quad_Quadrangle]);
  int nb0d_face0 = aVec[SMDSEntity_Node];

  std::vector<int> aResVec(SMDSEntity_Last);
  for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
  if(IsQuadratic) {
    aResVec[SMDSEntity_Quad_Hexa] = nb2d_face0 * ( nb2d/nb1d );
//.........这里部分代码省略.........