MeshCheck.cpp

Go to the documentation of this file.

/*
 * Copyright 2011-2016 Universiteit Antwerpen
 *
 * Licensed under the EUPL, Version 1.1 or  as soon they will be approved by
 * the European Commission - subsequent versions of the EUPL (the "Licence");
 * You may not use this work except in compliance with the Licence.
 * You may obtain a copy of the Licence at: http://ec.europa.eu/idabc/eupl5
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the Licence is distributed on an "AS IS" basis,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the Licence for the specific language governing
 * permissions and limitations under the Licence.
 */
#include "MeshCheck.h"

#include "Cell.h"
#include "Edge.h"
#include "Mesh.h"
#include "NeighborNodes.h"
#include "Node.h"

#include "util/container/CircularIterator.h"
#include "util/container/ConstCircularIterator.h"
#include "util/misc/log_debug.h"

#include <boost/property_tree/ptree.hpp>
#include <algorithm>
#include <array>
#include <cmath>
#include <functional>
#include <iomanip>
#include <ios>
#include <iostream>
#include <iterator>
#include <list>
#include <memory>
#include <set>
#include <sstream>
#include <string>
#include <tuple>
#include <vector>

namespace SimPT_Sim { class Wall; }

using namespace boost::property_tree;
using namespace std;
using namespace SimPT_Sim;
using namespace SimPT_Sim::Container;

namespace SimPT_Sim {

MeshCheck::MeshCheck(const Mesh& mesh)
        :  m_mesh(mesh)
{
}

bool MeshCheck::CheckAll() const
{
        auto overall = true;
        overall = overall && CheckAreas();
        overall = overall && CheckAtBoundaryNodes();
        overall = overall && CheckCellBoundaryWallNodes();
        overall = overall && CheckCellBoundaryWalls();
        overall = overall && CheckCellIdsSequence();
        overall = overall && CheckCellIdsUnique();
        overall = overall && CheckEdgeOwners();
        overall = overall && CheckMutuallyNeighbors();
        overall = overall && CheckNodeIdsSequence();
        overall = overall && CheckNodeIdsUnique();
        overall = overall && CheckNodeOwningNeighbors();
        overall = overall && CheckNodeOwningWalls();
        overall = overall && CheckWallIdsSequence();
        overall = overall && CheckWallIdsUnique();
        overall = overall && CheckWallNeighborsList();
        return overall;
}

bool MeshCheck::CheckAreas() const
{
        double area = 0.0;
        for (const auto& c: m_mesh.GetCells()) {
                area += c->GetArea();
        }
        const double delta = abs(area - m_mesh.GetBoundaryPolygon()->GetArea());

        const auto status = (delta < 1.0e-7);
        assert( status && "Sum of cell areas disagrees with boundary polygon area!");

        return status;
}

bool MeshCheck::CheckAtBoundaryNodes() const
{
        auto status = true;
        for (const auto& node : m_mesh.GetNodes()) {
                if (!(node->IsAtBoundary() == m_mesh.IsInBoundaryPolygon(node))) {
                        status = false;
                        break;
                }
        }
        assert( status && "Nodes and boundary_polgon disagree on at_boundary status!");
        return status;
}

bool MeshCheck::CheckCellBoundaryWallNodes() const
{
        auto status = true;
        list<Cell*> bw_cells;

        // Select cells with a boundary wall.
        for (auto cell : m_mesh.GetCells()) {
                if (cell->IsWallNeighbor(m_mesh.GetBoundaryPolygon())) {
                        bw_cells.push_back(cell);
                }
        }

        // For each cell, isolate the boundary wall and check that
        // the intervening nodes are really on the boundary polygon.
        for (auto cell : bw_cells) {
                auto is_bw = [this](Wall* w) {
                        return this->m_mesh.IsAtBoundary(w);
                };
                const auto& walls = cell->GetWalls();

                //TODO: What if there are multiple boundary walls?
                const auto it = find_if(walls.begin(), walls.end(), is_bw);
                if (it == walls.end()) {
                        status = false;
                        break;
                }
                const auto bw = *it;
                const Node* start = bw->GetN1();
                const Node* stop  = bw->GetN2();
                if (cell == bw->GetC2()) {
                        swap(start, stop);
                }
                auto& nodes = cell->GetNodes();
                const auto start_it = find(nodes.begin(), nodes.end(), start);
                auto c_it = make_circular(nodes.begin(), nodes.end(), start_it);
                do {
                        if (!m_mesh.IsInBoundaryPolygon(*c_it)) {
                                status = false;
                                break;
                        };
                } while(*(++c_it) != stop);
                if (!status) {
                        break;
                }
        }

        assert( status && "Boundary walls have nodes not on boundary!");
        return status;
}

bool MeshCheck::CheckCellBoundaryWalls() const
{
        auto status = true;
        if (m_mesh.GetCells().size() > 1U) {
                for (auto cell : m_mesh.GetCells()) {
                        if (cell->IsWallNeighbor(m_mesh.GetBoundaryPolygon()) != cell->HasBoundaryWall()) {
                                status = false;
                                break;
                        }
                }
        }
        assert( status && "Node and wall criteria disagree in at_boundary status!");
        return status;
}

bool MeshCheck::CheckCellIdsSequence() const
{
        bool status = true;

        auto cells = m_mesh.GetCells();
        for(unsigned int i = 0; i < cells.size(); ++i) {
                if (static_cast<unsigned int>(cells[i]->GetIndex()) != i) {
                        status = false;
                        break;
                }
        }

        assert( status && "Cell Ids not in sequence.");
        return status;
}

bool MeshCheck::CheckCellIdsUnique() const
{
        bool status = true;

        for (const auto& c1: m_mesh.GetCells()) {
                for (const auto& c2: m_mesh.GetCells()) {
                        if (c1->GetIndex() == c2->GetIndex() && c1 != c2) {
                                status = false;
                                break;
                        }
                }
                if (status == false) {
                        break;
                }
        };

        assert( status && "Cell Ids not unique.");
        return status;
}
bool MeshCheck::CheckEdgeOwners() const
{
        auto status = true;
        for (auto cell : m_mesh.GetCells()) {
                const auto& nodes = cell->GetNodes();
                const auto& walls = cell->GetWalls();
                for (const auto& wall : walls) {
                        const auto cit1      = find(nodes.begin(), nodes.end(), wall->GetN1());
                        const auto cit2      = find(nodes.begin(), nodes.end(), wall->GetN2());
                        const auto cit_start = (wall->GetC1() == cell) ? cit1 : cit2;
                        const auto cit_stop  = (wall->GetC1() == cell) ? cit2 : cit1;
                        const auto nghb_cell = (wall->GetC1() == cell) ? wall->GetC2() : wall->GetC1();

                        auto cit = make_const_circular(nodes, cit_start);
                        do {
                                const Edge e(*cit, *next(cit));
                                const auto edge_owners = m_mesh.GetEdgeOwners(e);
                                status = cell->HasEdge(e) && nghb_cell->HasEdge(e);
                                status = status && (edge_owners.size() == 4 );
                                for (const auto& owner : edge_owners) {
                                        status = status &&
                                                (owner.GetCell() == cell || owner.GetCell() == nghb_cell);
                                }
                                if (!status) {
                                        break;
                                }
                        } while(++cit != cit_stop);

                        if (!status) {
                                break;
                        }
                }
                if (!status) {
                        break;
                }
        }

        assert( status && "EdgeOwners relation not OK!");
        return status;
}

bool MeshCheck::CheckMutuallyNeighbors() const
{
        auto status = true;
        for (auto cell1 : m_mesh.GetCells()) {
                for (auto cell2 : m_mesh.GetCells()) {
                        if (cell1->IsWallNeighbor(cell2) != cell2->IsWallNeighbor(cell1)) {
                                status = false;
                                break;
                        }
                }
        }

        assert( status && "Neighbor relation is not symmetric!");
        return status;
}

bool MeshCheck::CheckNodeIdsSequence() const
{
        bool status = true;

        const auto& nodes = m_mesh.GetNodes();
        for(unsigned int i = 0; i < nodes.size(); ++i) {
                if (static_cast<unsigned int>(nodes[i]->GetIndex()) != i) {
                        status = false;
                        break;
                }
        }

        assert( status && "Cell Ids not in sequence.");
        return status;
}

bool MeshCheck::CheckNodeIdsUnique() const
{
        bool status = true;
        for (const auto& n1: m_mesh.GetNodes()) {
                for (const auto& n2: m_mesh.GetNodes()) {
                        if (n1->GetIndex() == n2->GetIndex() && n1 != n2) {
                                status = false;
                                break;
                        }
                }
                if (status == false) {
                        break;
                }
        };

        assert( status && "Node Ids not unique.");
        return status;
}

bool MeshCheck::CheckNodeOwningNeighbors() const
{
        auto status = true;
        for (auto cell : m_mesh.GetCells()) {
                const auto& nodes = cell->GetNodes();
                auto cit_node  = make_const_circular(nodes, nodes.begin());
                do {
                        NeighborNodes nn(cell, *prev(cit_node), *next(cit_node));
                        const auto&  owner = m_mesh.GetNodeOwningNeighbors(*cit_node);
                        const auto owner_it = find(owner.begin(), owner.end(), nn);
                        const auto owner_it2 = find(next(owner_it), owner.end(), nn);
                        status = status && (owner_it != owner.end()) && (owner_it2 == owner.end());
                } while(++cit_node != nodes.begin());
        }

        assert( status && "NodeOwningNeighbors relation inconsistent!");
        return status;
}

bool MeshCheck::CheckNodeOwningWalls() const
{
        auto status = true;

        // Make a copy of the node owning walls so that we can remove every wall
        // and check afterwards whether the map is empty
        auto n_o_w_copy = m_mesh.m_node_owning_walls;
        for (const auto& wall : m_mesh.GetWalls()) {
                const auto& nodes = wall->GetC1()->GetNodes();
                auto from = make_const_circular(nodes, find(nodes.begin(), nodes.end(), wall->GetN1()));
                auto to   = make_const_circular(nodes, find(nodes.begin(), nodes.end(), wall->GetN2()));

                // Check whether the endpoints of the wall are contained in C1
                assert(from != nodes.end() && "Couldn't find N1 in C1, wall inconsistent");
                assert(to   != nodes.end() && "Couldn't find N2 in C1, wall inconsistent");
                if (from == nodes.end() || to == nodes.end()) {
                        status = false;
                        continue;
                }

                // Check whether this wall is in the node owning walls of all his nodes
                do {
                        auto& owning_walls = n_o_w_copy[*from];
                        const auto it = std::find(owning_walls.begin(), owning_walls.end(), wall);
                        if (it == owning_walls.end()) {
                                status = false;
                        } else {
                                owning_walls.erase(it);
                        }
                } while (from++ != to);
        }

        assert(status && "Not every wall is associated with all its nodes");

        // Now check whether there are no more entries in the node owning walls than the ones we just found
        for (const auto& pair : n_o_w_copy) {
                if (!pair.second.empty()) {
                        status = false;
                }
        }

        assert(status && "Some node has an owning walls of which it is not part");
        return status;
}

bool MeshCheck::CheckWallIdsSequence() const
{
        bool status = true;

        const auto& walls = m_mesh.GetWalls();
        for(unsigned int i = 0; i < walls.size(); ++i) {
                if (static_cast<unsigned int>(walls[i]->GetIndex()) != i) {
                        status = false;
                        break;
                }
        }

        assert( status && "Wall Ids not in sequence.");
        return status;
}

bool MeshCheck::CheckWallIdsUnique() const
{
        bool status = true;
        for (const auto& w1: m_mesh.GetNodes()) {
                for (const auto& w2: m_mesh.GetNodes()) {
                        if (w1->GetIndex() == w2->GetIndex() && w1 != w2) {
                                status = false;
                                break;
                        }
                }
                if (status == false) {
                        break;
                }
        };

        assert( status && "Wall Ids not unique.");
        return status;
}

bool MeshCheck::CheckWallNeighborsList() const
{
        auto status = true;
        for (auto cell1 : m_mesh.GetCells()) {
                for (auto cell2 : m_mesh.GetCells()) {
                        if (cell1->IsWallNeighbor(cell2)) {
                                if (m_mesh.m_wall_neighbors.count(cell1) != 1) {
                                        status = false;
                                        break;
                                }
                                const auto& nb_list = m_mesh.m_wall_neighbors.at(cell1);
                                const auto it = find(nb_list.begin(), nb_list.end(), cell2);
                                if (it == nb_list.end()) {
                                        status = false;
                                        break;
                                }
                        }
                }
        }

        assert( status && "Wall_neighbor list in mesh has inconsistent info!");
        return status;
}

}