C++ NodeMap::SetNewIndex方法代码示例

void MutableMesh<ELEMENT_DIM, SPACE_DIM>::ReMesh(NodeMap& map)
{
    // Make sure that we are in the correct dimension - this code will be eliminated at compile time
    #define COVERAGE_IGNORE
    assert( ELEMENT_DIM == SPACE_DIM );
    #undef COVERAGE_IGNORE

    // Avoid some triangle/tetgen errors: need at least four
    // nodes for tetgen, and at least three for triangle
    if (GetNumNodes() <= SPACE_DIM)
    {
        EXCEPTION("The number of nodes must exceed the spatial dimension.");
    }

    // Make sure the map is big enough
    map.Resize(this->GetNumAllNodes());
    if (mAddedNodes || !mDeletedNodeIndices.empty())
    {
        // Size of mesh is about to change
        if (this->mpDistributedVectorFactory)
        {
            delete this->mpDistributedVectorFactory;
            this->mpDistributedVectorFactory = new DistributedVectorFactory(this->GetNumNodes());
        }
    }
    if (SPACE_DIM==1)
    {
        // Store the node locations
        std::vector<c_vector<double, SPACE_DIM> > old_node_locations;
        unsigned new_index = 0;
        for (unsigned i=0; i<this->GetNumAllNodes(); i++)
        {
            if (this->mNodes[i]->IsDeleted())
            {
                map.SetDeleted(i);
            }
            else
            {
                map.SetNewIndex(i, new_index);
                old_node_locations.push_back(this->mNodes[i]->rGetLocation());
                new_index++;
            }
        }

        // Remove current data
        Clear();

        // Construct the nodes and boundary nodes
        for (unsigned node_index=0; node_index<old_node_locations.size(); node_index++)
        {
            Node<SPACE_DIM>* p_node = new Node<SPACE_DIM>(node_index, old_node_locations[node_index], false);
            this->mNodes.push_back(p_node);

            // As we're in 1D, the boundary nodes are simply at either end of the mesh
            if ( node_index==0 || node_index==old_node_locations.size()-1 )
            {
                this->mBoundaryNodes.push_back(p_node);
            }
        }

        // Create a map between node indices and node locations
        std::map<double, unsigned> location_index_map;
        for (unsigned i=0; i<this->mNodes.size(); i++)
        {
            location_index_map[this->mNodes[i]->rGetLocation()[0]] = this->mNodes[i]->GetIndex();
        }

        // Use this map to generate a vector of node indices that are ordered spatially
        std::vector<unsigned> node_indices_ordered_spatially;
        for (std::map<double, unsigned>::iterator iter = location_index_map.begin();
             iter != location_index_map.end();
             ++iter)
        {
            node_indices_ordered_spatially.push_back(iter->second);
        }

        // Construct the elements
        this->mElements.reserve(old_node_locations.size()-1);
        for (unsigned element_index=0; element_index<old_node_locations.size()-1; element_index++)
        {
            std::vector<Node<SPACE_DIM>*> nodes;
            for (unsigned j=0; j<2; j++)
            {
                unsigned global_node_index = node_indices_ordered_spatially[element_index + j];
                assert(global_node_index < this->mNodes.size());
                nodes.push_back(this->mNodes[global_node_index]);
            }
            this->mElements.push_back(new Element<ELEMENT_DIM, SPACE_DIM>(element_index, nodes));
        }

        // Construct the two boundary elements - as we're in 1D, these are simply at either end of the mesh
        std::vector<Node<SPACE_DIM>*> nodes;
        nodes.push_back(this->mNodes[0]);
        this->mBoundaryElements.push_back(new BoundaryElement<ELEMENT_DIM-1, SPACE_DIM>(0, nodes));

        nodes.clear();
        nodes.push_back(this->mNodes[old_node_locations.size()-1]);
        this->mBoundaryElements.push_back(new BoundaryElement<ELEMENT_DIM-1, SPACE_DIM>(1, nodes));

        this->RefreshJacobianCachedData();
//.........这里部分代码省略.........

//.........这里部分代码省略.........
    assert(mRightOriginals.size() == mRightImages.size());

    for (unsigned i=0; i<mLeftOriginals.size(); i++)
    {
        mLeftOriginals[i] = big_map.GetNewIndex(mLeftOriginals[i]);
        mLeftImages[i] = big_map.GetNewIndex(mLeftImages[i]);
        mImageToLeftOriginalNodeMap[mLeftImages[i]] = mLeftOriginals[i];
    }

    for (unsigned i=0; i<mRightOriginals.size(); i++)
    {
        mRightOriginals[i] = big_map.GetNewIndex(mRightOriginals[i]);
        mRightImages[i] = big_map.GetNewIndex(mRightImages[i]);
        mImageToRightOriginalNodeMap[mRightImages[i]] = mRightOriginals[i];
    }

    for (unsigned i=0; i<mTopHaloNodes.size(); i++)
    {
        mTopHaloNodes[i] = big_map.GetNewIndex(mTopHaloNodes[i]);
        mBottomHaloNodes[i] = big_map.GetNewIndex(mBottomHaloNodes[i]);
    }

    /*
     * Check that elements crossing the periodic boundary have been meshed
     * in the same way at each side.
     */
    CorrectNonPeriodicMesh();

    /*
     * Take the double-sized mesh, with its new boundary elements, and
     * remove the relevant nodes, elements and boundary elements to leave
     * a proper periodic mesh.
     */
    ReconstructCylindricalMesh();

    DeleteHaloNodes();

    /*
     * Create a random boundary element between two nodes of the first
     * element if it is not deleted. This is a temporary measure to get
     * around re-index crashing when there are no boundary elements.
     */
    unsigned num_elements = GetNumAllElements();
    bool boundary_element_made = false;
    unsigned elem_index = 0;

    while (elem_index<num_elements && !boundary_element_made)
    {
        Element<2,2>* p_element = GetElement(elem_index);
        if (!p_element->IsDeleted())
        {
            boundary_element_made = true;
            std::vector<Node<2>*> nodes;
            nodes.push_back(p_element->GetNode(0));
            nodes.push_back(p_element->GetNode(1));
            BoundaryElement<1,2>* p_boundary_element = new BoundaryElement<1,2>(0, nodes);
            p_boundary_element->RegisterWithNodes();
            mBoundaryElements.push_back(p_boundary_element);
            this->mBoundaryElementWeightedDirections.push_back(zero_vector<double>(2));
            this->mBoundaryElementJacobianDeterminants.push_back(0.0);
        }
        elem_index++;
    }

    // Now call ReIndex() to remove the temporary nodes which are marked as deleted
    NodeMap reindex_map(GetNumAllNodes());
    ReIndex(reindex_map);
    assert(!reindex_map.IsIdentityMap());  // maybe don't need this

    /*
     * Go through the reindex map and use it to populate the original NodeMap
     * (the one that is returned to the user).
     */
    for (unsigned i=0; i<rMap.GetSize(); i++) // only going up to be size of map, not size of reindex_map
    {
        if (reindex_map.IsDeleted(i))
        {
            /*
             * i < num_original_nodes and node is deleted, this should correspond to
             * a node that was labelled as before the remeshing, so should have already
             * been set as deleted in the map above.
             */
            assert(rMap.IsDeleted(i));
        }
        else
        {
            rMap.SetNewIndex(i, reindex_map.GetNewIndex(i));
        }
    }

    // We can now clear the index vectors and maps; they are only used for remeshing
    mLeftOriginals.clear();
    mLeftImages.clear();
    mImageToLeftOriginalNodeMap.clear();
    mRightOriginals.clear();
    mRightImages.clear();
    mImageToRightOriginalNodeMap.clear();
    mLeftPeriodicBoundaryElementIndices.clear();
    mRightPeriodicBoundaryElementIndices.clear();
}

void MutableMesh<ELEMENT_DIM, SPACE_DIM>::ReIndex(NodeMap& map)
{
    assert(!mAddedNodes);
    map.Resize(this->GetNumAllNodes());

    std::vector<Element<ELEMENT_DIM, SPACE_DIM> *> live_elements;

    for (unsigned i=0; i<this->mElements.size(); i++)
    {
        assert(i==this->mElements[i]->GetIndex()); // We need this to be true to be able to reindex the Jacobian cache
        if (this->mElements[i]->IsDeleted())
        {
            delete this->mElements[i];
        }
        else
        {
            live_elements.push_back(this->mElements[i]);

            unsigned this_element_index = this->mElements[i]->GetIndex();
            if (SPACE_DIM == ELEMENT_DIM)
            {
                this->mElementJacobians[live_elements.size()-1] = this->mElementJacobians[this_element_index];
                this->mElementInverseJacobians[live_elements.size()-1] = this->mElementInverseJacobians[this_element_index];
            }
            else
            {
                this->mElementWeightedDirections[live_elements.size()-1] = this->mElementWeightedDirections[this_element_index];
            }
            this->mElementJacobianDeterminants[live_elements.size()-1] = this->mElementJacobianDeterminants[this_element_index];
        }
    }

    assert(mDeletedElementIndices.size() == this->mElements.size()-live_elements.size());
    mDeletedElementIndices.clear();
    this->mElements = live_elements;
    unsigned num_elements = this->mElements.size();

    if (SPACE_DIM == ELEMENT_DIM)
    {
        this->mElementJacobians.resize(num_elements);
        this->mElementInverseJacobians.resize(num_elements);
    }
    else
    {
        this->mElementWeightedDirections.resize(num_elements);
    }
    this->mElementJacobianDeterminants.resize(num_elements);

    std::vector<Node<SPACE_DIM> *> live_nodes;
    for (unsigned i=0; i<this->mNodes.size(); i++)
    {
        if (this->mNodes[i]->IsDeleted())
        {
            delete this->mNodes[i];
            map.SetDeleted(i);
        }
        else
        {
            live_nodes.push_back(this->mNodes[i]);
            // the nodes will have their index set to be the index into the live_nodes
            // vector further down
            map.SetNewIndex(i, (unsigned)(live_nodes.size()-1));
        }
    }

    assert(mDeletedNodeIndices.size() == this->mNodes.size()-live_nodes.size());
    this->mNodes = live_nodes;
    mDeletedNodeIndices.clear();

    std::vector<BoundaryElement<ELEMENT_DIM-1, SPACE_DIM> *> live_boundary_elements;
    for (unsigned i=0; i<this->mBoundaryElements.size(); i++)
    {
        if (this->mBoundaryElements[i]->IsDeleted())
        {
            delete this->mBoundaryElements[i];
        }
        else
        {
            live_boundary_elements.push_back(this->mBoundaryElements[i]);

            this->mBoundaryElementWeightedDirections[live_boundary_elements.size()-1] = this->mBoundaryElementWeightedDirections[this->mBoundaryElements[i]->GetIndex()];
            this->mBoundaryElementJacobianDeterminants[live_boundary_elements.size()-1] = this->mBoundaryElementJacobianDeterminants[this->mBoundaryElements[i]->GetIndex()];
        }
    }

    assert(mDeletedBoundaryElementIndices.size() == this->mBoundaryElements.size()-live_boundary_elements.size());
    this->mBoundaryElements = live_boundary_elements;
    mDeletedBoundaryElementIndices.clear();

    unsigned num_boundary_elements = this->mBoundaryElements.size();

    this->mBoundaryElementWeightedDirections.resize(num_boundary_elements);
    this->mBoundaryElementJacobianDeterminants.resize(num_boundary_elements);

    for (unsigned i=0; i<this->mNodes.size(); i++)
    {
        this->mNodes[i]->SetIndex(i);
    }

    for (unsigned i=0; i<this->mElements.size(); i++)
//.........这里部分代码省略.........