Xyst test code coverage report
Current view: top level - Inciter/AMR - mesh_adapter.cpp (source / functions) Hit Total Coverage
Commit: 5689ba12dc66a776d3d75f1ee48cc7d78eaa18dc Lines: 417 634 65.8 %
Date: 2024-11-22 19:17:03 Functions: 17 24 70.8 %
Legend: Lines: hit not hit | Branches: + taken - not taken # not executed Branches: 322 763 42.2 %

           Branch data     Line data    Source code
       1                 :            : #include "mesh_adapter.hpp"
       2                 :            : 
       3                 :            : #include <assert.h>                        // for assert
       4                 :            : #include <cstddef>                         // for size_t
       5                 :            : #include <iostream>                        // for operator<<, endl, basic_os...
       6                 :            : #include <set>                             // for set
       7                 :            : #include <utility>                         // for pair
       8                 :            : #include "AMR/AMR_types.hpp"                 // for Edge_Refinement, edge_list_t
       9                 :            : #include "AMR/Loggers.hpp"                   // for trace_out
      10                 :            : #include "AMR/Refinement_State.hpp"          // for Refinement_Case, Refinemen...
      11                 :            : #include "AMR/edge.hpp"                      // for operator<<, edge_t
      12                 :            : #include "AMR/edge_store.hpp"                // for edge_store_t
      13                 :            : #include "AMR/marked_refinements_store.hpp"  // for marked_refinements_store_t
      14                 :            : #include "AMR/node_connectivity.hpp"         // for node_connectivity_t
      15                 :            : #include "AMR/refinement.hpp"                // for refinement_t
      16                 :            : #include "AMR/tet_store.hpp"                 // for tet_store_t
      17                 :            : 
      18                 :            : #if defined(__clang__)
      19                 :            :   #pragma clang diagnostic push
      20                 :            :   #pragma clang diagnostic ignored "-Wunreachable-code"
      21                 :            :   #pragma clang diagnostic ignored "-Wdocumentation"
      22                 :            : #endif
      23                 :            : 
      24                 :            : namespace AMR {
      25                 :            : 
      26                 :            : 
      27                 :            : #ifdef ENABLE_NODE_STORE
      28                 :            :     /**
      29                 :            :      * @brief This accepts external coord arrays and allows the node_store to
      30                 :            :      * track the new node positions as they are added
      31                 :            :      *
      32                 :            :      * @param m_x X coodinates
      33                 :            :      * @param m_y Y coodinates
      34                 :            :      * @param m_z Z coodinates
      35                 :            :      * @param graph_size Total number of nodes
      36                 :            :      */
      37                 :            :     // TODO: remove graph size and use m.size()
      38                 :            :     // TODO: remove these pointers
      39                 :            :     //void mesh_adapter_t::init_node_store(coord_type* m_x, coord_type* m_y, coord_type* m_z)
      40                 :            :     //{
      41                 :            :     //    assert( m_x->size() == m_y->size() );
      42                 :            :     //    assert( m_x->size() == m_z->size() );
      43                 :            : 
      44                 :            :     //    node_store.set_x(*m_x);
      45                 :            :     //    node_store.set_y(*m_y);
      46                 :            :     //    node_store.set_z(*m_z);
      47                 :            :     //}
      48                 :            : #endif
      49                 :            : 
      50                 :         22 :     std::pair< bool, std::size_t > mesh_adapter_t::check_same_face(
      51                 :            :       std::size_t tet_id,
      52                 :            :       const std::unordered_set<std::size_t>& inactive_nodes)
      53                 :            :     {
      54         [ +  - ]:         22 :        edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
      55                 :            : 
      56 [ +  + ][ -  + ]:         22 :        Assert(inactive_nodes.size()==3 || inactive_nodes.size()==2,
         [ -  - ][ -  - ]
                 [ -  - ]
      57                 :            :          "Incorrectly sized inactive nodes set");
      58                 :            : 
      59                 :            :        // for a tet ABCD, the keys (edges) are ordered
      60                 :            :        // A-B, A-C, A-D, B-C, B-D, C-D
      61                 :            :        // 0-1, 0-2, 0-3, 1-2, 1-3, 2-3
      62                 :            : 
      63                 :            :        std::array< std::array< std::size_t, 3 >, 4 >
      64                 :            :          edges_on_face;
      65                 :            : 
      66                 :            :        // A-B-C
      67                 :         22 :        edges_on_face[0][0] =
      68         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[0].get_data());
      69                 :         22 :        edges_on_face[0][1] =
      70         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[1].get_data());
      71                 :         22 :        edges_on_face[0][2] =
      72         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[3].get_data());
      73                 :            : 
      74                 :            :        // A-B-D
      75                 :         22 :        edges_on_face[1][0] =
      76         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[0].get_data());
      77                 :         22 :        edges_on_face[1][1] =
      78         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[2].get_data());
      79                 :         22 :        edges_on_face[1][2] =
      80         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[4].get_data());
      81                 :            : 
      82                 :            :        // B-C-D
      83                 :         22 :        edges_on_face[2][0] =
      84         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[3].get_data());
      85                 :         22 :        edges_on_face[2][1] =
      86         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[4].get_data());
      87                 :         22 :        edges_on_face[2][2] =
      88         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[5].get_data());
      89                 :            : 
      90                 :            :        // A-C-D
      91                 :         22 :        edges_on_face[3][0] =
      92         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[1].get_data());
      93                 :         22 :        edges_on_face[3][1] =
      94         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[2].get_data());
      95                 :         22 :        edges_on_face[3][2] =
      96         [ +  - ]:         22 :          tk::cref_find(node_connectivity.data(),edge_list[5].get_data());
      97                 :            : 
      98                 :            :        //Iterate over edges to determine if inactive_nodes are all part of a face
      99                 :         22 :        bool same_face(false);
     100                 :         22 :        [[maybe_unused]] bool tnode_set(false);
     101                 :         22 :        std::size_t third_node = 0;
     102         [ +  + ]:        110 :        for(const auto& face : edges_on_face)
     103                 :            :        {
     104                 :         88 :          std::size_t icount = 0;
     105         [ +  + ]:        352 :          for (const auto& np_node : face) {
     106 [ +  - ][ +  + ]:        264 :            if (inactive_nodes.count(np_node)) ++icount;
     107                 :            :          }
     108         [ +  + ]:         88 :          if (inactive_nodes.size() == icount) {
     109                 :         16 :            same_face = true;
     110                 :            :            // if the two inactive_nodes being checked are on the same parent
     111                 :            :            // face, determine the third node on that face
     112         [ +  + ]:         16 :            if (inactive_nodes.size() == 2) {
     113         [ +  - ]:          2 :              for (auto fn:face) {
     114 [ +  - ][ +  - ]:          2 :                if (inactive_nodes.count(fn) == 0) {
     115                 :          2 :                  third_node = fn;
     116                 :          2 :                  tnode_set = true;
     117                 :          2 :                  break;
     118                 :            :                }
     119                 :            :              }
     120                 :            :            }
     121                 :            :          }
     122                 :            :        }
     123                 :            : 
     124 [ +  + ][ +  + ]:         22 :        if (same_face && inactive_nodes.size() == 2)
                 [ +  + ]
     125 [ -  + ][ -  - ]:          2 :          Assert(tnode_set, "Third node on face not set in derefine");
         [ -  - ][ -  - ]
     126                 :            : 
     127                 :         22 :        return {same_face, third_node};
     128                 :            :     }
     129                 :            : 
     130                 :            :     /** @brief Consume an existing mesh, and turn it into the AMRs
     131                 :            :      * representations of tets and nodes
     132                 :            :      *
     133                 :            :      * @param tetinpoel Vector of nodes grouped together in blocks of 4 to
     134                 :            :      * represent tets
     135                 :            :      */
     136                 :       2468 :     void mesh_adapter_t::consume_tets(const std::vector< std::size_t >& tetinpoel )
     137                 :            :     {
     138         [ +  + ]:     927918 :         for (size_t i = 0; i < tetinpoel.size(); i+=4)
     139                 :            :         {
     140                 :            :             tet_t t = {
     141                 :            :                 {
     142                 :     925450 :                     tetinpoel[i],
     143                 :     925450 :                     tetinpoel[i+1],
     144                 :     925450 :                     tetinpoel[i+2],
     145                 :     925450 :                     tetinpoel[i+3]
     146                 :            :                 }
     147                 :     925450 :             };
     148                 :            : 
     149                 :     925450 :             trace_out << "Consume tet " << i << std::endl;
     150         [ +  - ]:     925450 :             tet_store.add(t, AMR::Refinement_Case::initial_grid);
     151                 :            :         }
     152                 :       2468 :     }
     153                 :            : 
     154                 :            :     /**
     155                 :            :      * @brief Place holder function to evaluate error estimate at
     156                 :            :      * nodes, and therefor mark things as needing to be refined
     157                 :            :      */
     158                 :            :     //void mesh_adapter_t::evaluate_error_estimate() {
     159                 :            :     //    for (auto& kv : tet_store.edge_store.edges)
     160                 :            :     //    {
     161                 :            :     //        // Mark them as needing refinement
     162                 :            :     //        if (kv.second.refinement_criteria > refinement_cut_off)
     163                 :            :     //        {
     164                 :            :     //            kv.second.needs_refining = 1;
     165                 :            :     //        }
     166                 :            :     //        else
     167                 :            :     //        {
     168                 :            :     //            // TODO: Check this won't be overwriting valuable
     169                 :            :     //            // information from last iteration
     170                 :            :     //            kv.second.needs_refining = 0;
     171                 :            :     //        }
     172                 :            :     //    }
     173                 :            :     //}
     174                 :            : 
     175                 :            :     /**
     176                 :            :      * @brief Helper function to apply uniform refinement to all tets
     177                 :            :      */
     178                 :         31 :     void mesh_adapter_t::mark_uniform_refinement()
     179                 :            :     {
     180         [ +  + ]:      33650 :         for (auto& kv : tet_store.edge_store.edges) {
     181                 :      33619 :            auto& local = kv.second;
     182         [ +  - ]:      33619 :            if (local.lock_case == Edge_Lock_Case::unlocked)
     183                 :      33619 :              local.needs_refining = 1;
     184                 :            :         }
     185                 :         31 :         mark_refinement();
     186                 :         31 :     }
     187                 :            : 
     188                 :            :     /**
     189                 :            :      * @brief Helper function to apply uniform derefinement to all tets
     190                 :            :      */
     191                 :         24 :     void mesh_adapter_t::mark_uniform_derefinement()
     192                 :            :     {
     193                 :         24 :       const auto& inp = tet_store.get_active_inpoel();
     194                 :         24 :       auto& edge_store = tet_store.edge_store;
     195         [ +  + ]:      49409 :       for (std::size_t t=0; t<inp.size()/4; ++t) {
     196                 :            :         const auto edges =
     197         [ +  - ]:      49385 :           edge_store.generate_keys(
     198                 :      49385 :             {inp[t*4+0], inp[t*4+1], inp[t*4+2], inp[t*4+3]});
     199         [ +  + ]:     345695 :         for (const auto& tetedge : edges) {
     200         [ +  - ]:     296310 :           auto e = edge_store.edges.find(tetedge);
     201         [ +  - ]:     296310 :           if (e != end(edge_store.edges)) {
     202                 :     296310 :             auto& local = e->second;
     203                 :     296310 :             local.needs_derefining = 1;
     204                 :            :           }
     205                 :            :         }
     206                 :            :       }
     207                 :         24 :       mark_derefinement();
     208                 :         24 :     }
     209                 :            : 
     210                 :            :     /**
     211                 :            :      * @brief For a given set of edges, set their refinement criteria for
     212                 :            :      * refinement
     213                 :            :      *
     214                 :            :      * @param remote Vector of edges and edge tags
     215                 :            :      */
     216                 :         16 :     void mesh_adapter_t::mark_error_refinement(
     217                 :            :             const std::vector< std::pair< edge_t, edge_tag > >& remote )
     218                 :            :     {
     219         [ +  + ]:       5360 :        for (const auto& r : remote) {
     220         [ +  - ]:       5344 :          auto& local = tet_store.edge_store.get( r.first );
     221         [ +  + ]:       5344 :          if (r.second == edge_tag::REFINE) {
     222         [ -  + ]:       3612 :            if (local.lock_case > Edge_Lock_Case::unlocked) {
     223                 :          0 :              local.needs_refining = 0;
     224                 :            :            } else {
     225                 :       3612 :              local.needs_refining = 1;
     226                 :            :              // an edge deemed to be 'needing refinement', cannot be derefined
     227                 :       3612 :              local.needs_derefining = 0;
     228                 :            :            }
     229         [ +  - ]:       1732 :          } else if (r.second == edge_tag::DEREFINE) {
     230                 :       1732 :            local.needs_derefining = 1;
     231                 :            :          }
     232                 :            :        }
     233                 :            : 
     234                 :         16 :        mark_refinement();
     235                 :         16 :        mark_derefinement();
     236                 :         16 :     }
     237                 :            : 
     238                 :          0 :    void mesh_adapter_t::mark_error_refinement_corr( const EdgeData& edges )
     239                 :            :     {
     240         [ -  - ]:          0 :        for (const auto& r : edges)
     241                 :            :        {
     242         [ -  - ]:          0 :            auto& edgeref = tet_store.edge_store.get( edge_t(r.first) );
     243                 :          0 :            edgeref.needs_refining = std::get<0>(r.second);
     244                 :          0 :            edgeref.needs_derefining = std::get<1>(r.second);
     245 [ -  - ][ -  - ]:          0 :            Assert(edgeref.lock_case == Edge_Lock_Case::unlocked ?
         [ -  - ][ -  - ]
         [ -  - ][ -  - ]
         [ -  - ][ -  - ]
         [ -  - ][ -  - ]
         [ -  - ][ -  - ]
                 [ -  - ]
     246                 :            :              edgeref.lock_case <= std::get<2>(r.second) : true,
     247                 :            :              "Edge " + std::to_string(r.first[0]) +
     248                 :            :              "-" + std::to_string(r.first[1]) +
     249                 :            :              " : current edge-lock " + std::to_string(edgeref.lock_case) +
     250                 :            :              " stricter than received edge-lock " +
     251                 :            :              std::to_string(std::get<2>(r.second)));
     252                 :          0 :            edgeref.lock_case = std::get<2>(r.second);
     253                 :            :        }
     254                 :          0 :        mark_refinement();
     255                 :          0 :        mark_derefinement();
     256                 :          0 :     }
     257                 :            : 
     258                 :            :     /**
     259                 :            :      * @brief Function to detect the compatibility class (1,
     260                 :            :      * 2, or 3) based on the number of locked edges and the existence
     261                 :            :      * of intermediate edges
     262                 :            :      *
     263                 :            :      * @param num_locked_edges The number of locked edges
     264                 :            :      * @param num_intermediate_edges The number of intermediate edges
     265                 :            :      * @param refinement_case The refinement case of the tet
     266                 :            :      * @param normal TODO: Document this!
     267                 :            :      *
     268                 :            :      * @return The compatibili4y class of the current scenario
     269                 :            :      */
     270                 :      50634 :     int mesh_adapter_t::detect_compatibility(
     271                 :            :             int num_locked_edges,
     272                 :            :             int num_intermediate_edges,
     273                 :            :             AMR::Refinement_Case refinement_case,
     274                 :            :             int normal
     275                 :            :     )
     276                 :            :     {
     277                 :      50634 :         int compatibility = 0;
     278                 :      50634 :         num_locked_edges += num_intermediate_edges;
     279                 :            : 
     280                 :            :         /*
     281                 :            :         // Split this into three categories
     282                 :            :         // 1. Normal elements without locked edges. => 1
     283                 :            :         //if (normal) {
     284                 :            : 
     285                 :            :             // 3. Intermediate elements with at least one edge marked for refinement => 3
     286                 :            :             if (num_intermediate_edges > 0)
     287                 :            :             {
     288                 :            :                 compatibility = 3;
     289                 :            :             }
     290                 :            :             else if (num_locked_edges == 0) {
     291                 :            :                 compatibility = 1;
     292                 :            :             }
     293                 :            :         // 2. Normal elements with locked edges. => 2
     294                 :            :             else {
     295                 :            :                 compatibility = 2;
     296                 :            :             }
     297                 :            :         //}
     298                 :            :         */
     299                 :            : 
     300                 :            :         //else {
     301                 :            :             //if (num_intermediate_edges > 0) { compatibility = 3; }
     302                 :            :         //}
     303                 :            : 
     304                 :            : 
     305                 :            : 
     306                 :            :         // Only 1:2 and 1:4 are intermediates and eligible for class3 // NOT TRUE!
     307                 :            :         /*
     308                 :            :         if (num_intermediate_edges > 0)
     309                 :            :         {
     310                 :            :             if (!normal) {
     311                 :            :                 trace_out << " not normal 3 " << std::endl;
     312                 :            :                 compatibility = 3;
     313                 :            :             }
     314                 :            :             else { // Attempt to allow for "normal" 1:4 and 1:8
     315                 :            :                 compatibility = 2;
     316                 :            :                 trace_out << " normal 3 " << std::endl;
     317                 :            :             }
     318                 :            : 
     319                 :            :         }
     320                 :            :         else {
     321                 :            :             if (num_locked_edges == 0) {
     322                 :            :                 trace_out << " no lock 1 " << std::endl;
     323                 :            :                 compatibility = 1;
     324                 :            :             }
     325                 :            :             else {
     326                 :            :                 trace_out << " lock 2 " << std::endl;
     327                 :            :                 compatibility = 2;
     328                 :            :             }
     329                 :            :         }
     330                 :            :         */
     331                 :            : 
     332                 :            : 
     333                 :            :         // Old implementation
     334                 :            :         // Only 1:2 and 1:4 are intermediates and eligible for class3 // NOT TRUE!
     335         [ +  - ]:      50634 :         if (
     336         [ -  + ]:      50634 :                 (refinement_case == AMR::Refinement_Case::one_to_two) or
     337                 :            :                 (refinement_case == AMR::Refinement_Case::one_to_four)
     338                 :            :            )
     339                 :            :         {
     340         [ -  - ]:          0 :             if (!normal) {
     341                 :          0 :                 trace_out << " not normal 3 " << std::endl;
     342                 :          0 :                 compatibility = 3;
     343                 :            :             }
     344                 :            :             else { // Attempt to allow for "normal" 1:4 and 1:8
     345                 :          0 :                 compatibility = 2;
     346                 :          0 :                 trace_out << " normal 3 " << std::endl;
     347                 :            :             }
     348                 :            : 
     349                 :            :         }
     350                 :            :         else {
     351         [ +  - ]:      50634 :             if (num_locked_edges == 0) {
     352                 :      50634 :                 trace_out << " no lock 1 " << std::endl;
     353                 :      50634 :                 compatibility = 1;
     354                 :            :             }
     355                 :            :             else {
     356                 :          0 :                 trace_out << " lock 2 " << std::endl;
     357                 :          0 :                 compatibility = 2;
     358                 :            :             }
     359                 :            :         }
     360                 :            : 
     361         [ -  + ]:      50634 :         assert(compatibility > 0);
     362         [ -  + ]:      50634 :         assert(compatibility < 4);
     363                 :      50634 :         return compatibility;
     364                 :            :     }
     365                 :            : 
     366                 :            :     /**
     367                 :            :      * @brief Function which implements the main refinement algorithm from
     368                 :            :      * the paper Iterating over the cells, deciding which refinement and
     369                 :            :      * compatibility types are appropriate etc
     370                 :            :      */
     371                 :         47 :     void mesh_adapter_t::mark_refinement() {
     372                 :            : 
     373                 :            : #ifndef AMR_MAX_ROUNDS
     374                 :            :         // Paper says the average actual num rounds will be 5-15
     375                 :            : #define AMR_MAX_ROUNDS 50
     376                 :            : #endif
     377                 :         47 :         const size_t max_num_rounds = AMR_MAX_ROUNDS;
     378                 :            : 
     379                 :            :         // Mark refinements
     380                 :            :         size_t iter;
     381                 :            :         //Iterate until convergence
     382         [ +  - ]:        107 :         for (iter = 0; iter < max_num_rounds; iter++)
     383                 :            :         {
     384                 :            : 
     385                 :        107 :             tet_store.marked_refinements.get_state_changed() = false;
     386                 :            : 
     387                 :            :             // Loop over Tets.
     388         [ +  + ]:      56949 :             for (const auto& kv : tet_store.tets)
     389                 :            :             {
     390                 :      56842 :                 size_t tet_id = kv.first;
     391                 :            : 
     392                 :      56842 :                 trace_out << "Process tet " << tet_id << std::endl;
     393                 :            : 
     394                 :            :                 // Only apply checks to tets on the active list
     395 [ +  - ][ +  + ]:      56842 :                 if (tet_store.is_active(tet_id)) {
     396                 :      51555 :                     int num_locked_edges = 0;
     397                 :      51555 :                     int num_intermediate_edges = 0;
     398                 :            : 
     399                 :            :                     // Loop over nodes and count the number which need refining
     400                 :      51555 :                     int num_to_refine = 0;
     401                 :            : 
     402                 :            :                     // This is useful for later inspection
     403         [ +  - ]:      51555 :                     edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
     404                 :            : 
     405                 :            :                     //Iterate over edges
     406         [ +  + ]:     360885 :                     for(auto & key : edge_list)
     407                 :            :                     {
     408                 :            : 
     409                 :     309330 :                         trace_out << "Edge " << key << std::endl;
     410                 :            : 
     411                 :            :                         //Count locked edges and edges in need of
     412                 :            :                         // refinement
     413                 :            :                         // Count Locked Edges
     414 [ +  - ][ -  + ]:     309330 :                         if(tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::locked)
     415                 :            :                         {
     416                 :          0 :                             trace_out << "Found locked edge " << key << std::endl;
     417                 :          0 :                             trace_out << "Locked :" << tet_store.edge_store.get(key).lock_case << std::endl;
     418                 :          0 :                             num_locked_edges++;
     419                 :            :                         }
     420         [ +  - ]:     309330 :                         else if(tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::intermediate
     421 [ +  - ][ +  - ]:     309330 :                           || tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::temporary)
         [ -  + ][ -  + ]
     422                 :            :                         {
     423                 :          0 :                             trace_out << "Found intermediate edge " << key << std::endl;
     424                 :          0 :                             num_intermediate_edges++;
     425                 :            :                         }
     426                 :            :                         else
     427                 :            :                         {
     428                 :            :                             // Count edges which need refining
     429                 :            :                             //  We check in here as we won't refine a
     430                 :            :                             //  locked edge and will thus ignore it
     431 [ +  - ][ +  + ]:     309330 :                             if (tet_store.edge_store.get(key).needs_refining == 1)
     432                 :            :                             {
     433                 :     298129 :                                 num_to_refine++;
     434                 :     298129 :                                 trace_out << "key needs ref " << key << std::endl;
     435                 :            :                             }
     436                 :            :                         }
     437                 :            :                     }
     438                 :            : 
     439                 :            :                     // TODO: Should this be a reference?
     440         [ +  - ]:      51555 :                     AMR::Refinement_Case refinement_case = tet_store.get_refinement_case(tet_id);
     441         [ +  - ]:      51555 :                     int normal = tet_store.is_normal(tet_id);
     442                 :            : 
     443                 :      51555 :                     trace_out << "Checking " << tet_id <<
     444                 :            :                         " ref case " << refinement_case <<
     445                 :            :                         " num ref " << num_to_refine <<
     446                 :            :                         " normal " << normal <<
     447                 :            :                         std::endl;
     448                 :            : 
     449                 :            : 
     450                 :            : 
     451                 :            :                     //If we have any tets to refine
     452         [ +  + ]:      51555 :                     if (num_to_refine > 0)
     453                 :            :                     {
     454                 :            :                         //Determine compatibility case
     455                 :            : 
     456         [ +  - ]:      50634 :                         int compatibility = detect_compatibility(num_locked_edges,
     457                 :            :                                 num_intermediate_edges, refinement_case, normal);
     458                 :            : 
     459                 :      50634 :                         trace_out << "Compat " << compatibility << std::endl;
     460                 :            : 
     461                 :            :                         // Now check num_to_refine against situations
     462         [ +  - ]:      50634 :                         if (compatibility == 1)
     463                 :            :                         {
     464         [ +  - ]:      50634 :                             refinement_class_one(num_to_refine, tet_id);
     465                 :            :                         }
     466         [ -  - ]:          0 :                         else if (compatibility == 2)
     467                 :            :                         {
     468         [ -  - ]:          0 :                             refinement_class_two(edge_list, tet_id);
     469                 :            :                         }
     470         [ -  - ]:          0 :                         else if (compatibility == 3)
     471                 :            :                         {
     472         [ -  - ]:          0 :                             refinement_class_three(tet_id);
     473                 :            :                         }
     474                 :            : 
     475                 :            :                         /*
     476                 :            :                         // Write temp mesh out
     477                 :            :                         std::string temp_file =  "temp." +
     478                 :            :                         std::to_string(iter) + "." +
     479                 :            :                         std::to_string(tet_id) + ".exo";
     480                 :            : 
     481                 :            :                         std::cout << "Writing " << temp_file << std::endl;
     482                 :            :                         Adaptive_UnsMesh outmesh(
     483                 :            :                         get_active_inpoel(), x(), y(), z()
     484                 :            :                         );
     485                 :            :                         tk::ExodusIIMeshWriter( temp_file, tk::ExoWriter::CREATE ).
     486                 :            :                         writeMesh(outmesh);
     487                 :            :                         */
     488                 :            : 
     489                 :            :                     } // if num_to_refine
     490                 :            :                     else {
     491                 :            :                         // If we got here, we don't want to refine this guy
     492         [ +  - ]:        921 :                         tet_store.marked_refinements.add(tet_id, AMR::Refinement_Case::none);
     493                 :            :                     }
     494                 :            :                 } // if active
     495                 :            :                 else {
     496                 :       5287 :                     trace_out << "Inactive" << std::endl;
     497                 :            :                 }
     498                 :            :             } // For
     499                 :            : 
     500                 :            :             // If nothing changed during that round, break
     501         [ +  + ]:        107 :             if (!tet_store.marked_refinements.get_state_changed())
     502                 :            :             {
     503                 :         47 :                 trace_out << "Terminating loop at iter " << iter << std::endl;
     504                 :         47 :                 break;
     505                 :            :             }
     506                 :         60 :             trace_out << "End iter " << iter << std::endl;
     507                 :            :         }
     508                 :         47 :         trace_out << "Loop took " << iter << " rounds." << std::endl;
     509                 :            : 
     510                 :            :         //std::cout << "Print Tets" << std::endl;
     511                 :            :         //print_tets();
     512                 :         47 :     }
     513                 :            : 
     514                 :            :     /**
     515                 :            :      * @brief Helper function to print tet information when needed
     516                 :            :      */
     517                 :          0 :     void mesh_adapter_t::print_tets() {
     518                 :          0 :         tet_store.print_tets();
     519                 :          0 :     }
     520                 :            : 
     521                 :            :     /**
     522                 :            :      * @brief Function to call refinement after each tet has had it's
     523                 :            :      * refinement case marked and calculated
     524                 :            :      */
     525                 :         47 :     void mesh_adapter_t::perform_refinement()
     526                 :            :     {
     527                 :            :         // Track tets which need to be deleted this iteration
     528                 :         47 :         std::set<size_t> round_two;
     529                 :            : 
     530                 :         47 :         trace_out << "Perform ref" << std::endl;
     531                 :            : 
     532                 :            :         // Do refinements
     533         [ +  + ]:     137891 :         for (const auto& kv : tet_store.tets)
     534                 :            :         {
     535                 :     137844 :             size_t tet_id = kv.first;
     536                 :            : 
     537                 :     137844 :             trace_out << "Do refine of " << tet_id << std::endl;
     538 [ +  - ][ +  + ]:     137844 :             if (tet_store.has_refinement_decision(tet_id))
     539                 :            :             {
     540 [ +  - ][ +  + ]:      23419 :                 switch(tet_store.marked_refinements.get(tet_id))
         [ +  - ][ -  + ]
                    [ - ]
     541                 :            :                 {
     542                 :        130 :                     case AMR::Refinement_Case::one_to_two:
     543         [ +  - ]:        130 :                         refiner.refine_one_to_two(tet_store,node_connectivity,tet_id);
     544                 :        130 :                         break;
     545                 :        176 :                     case AMR::Refinement_Case::one_to_four:
     546         [ +  - ]:        176 :                         refiner.refine_one_to_four(tet_store,node_connectivity,tet_id);
     547                 :        176 :                         break;
     548                 :      22907 :                     case AMR::Refinement_Case::one_to_eight:
     549         [ +  - ]:      22907 :                         refiner.refine_one_to_eight(tet_store,node_connectivity,tet_id);
     550                 :      22907 :                         break;
     551                 :          0 :                     case AMR::Refinement_Case::two_to_eight:
     552 [ -  - ][ -  - ]:          0 :                         round_two.insert( tet_store.get_parent_id(tet_id) );
     553                 :            :                         //std::cout << "2->8\n";
     554                 :          0 :                         break;
     555                 :          0 :                     case AMR::Refinement_Case::four_to_eight:
     556 [ -  - ][ -  - ]:          0 :                         round_two.insert( tet_store.get_parent_id(tet_id));
     557                 :            :                         //std::cout << "4->8\n";
     558                 :          0 :                         break;
     559                 :        206 :                     case AMR::Refinement_Case::initial_grid:
     560                 :            :                         // Do nothing
     561                 :            :                     case AMR::Refinement_Case::none:
     562                 :            :                         // Do nothing
     563                 :        206 :                         break;
     564                 :            :                         // No need for default as enum is explicitly covered
     565                 :            :                 }
     566                 :            :                 // Mark tet as not needing refinement
     567         [ +  - ]:      23419 :                 tet_store.marked_refinements.erase(tet_id);
     568                 :            :             }
     569                 :            :         }
     570                 :            : 
     571                 :         47 :         trace_out << "round_two size " << round_two.size() << std::endl;
     572         [ -  + ]:         47 :         for (const auto i : round_two)
     573                 :            :         {
     574                 :          0 :             trace_out << "round two i " << i << std::endl;
     575                 :            : 
     576                 :            :             // Cache children as we're about to change this data
     577 [ -  - ][ -  - ]:          0 :             auto former_children = tet_store.data(i).children;
     578                 :            : 
     579         [ -  - ]:          0 :             AMR::Refinement_State& element = tet_store.data(i);
     580                 :            : 
     581         [ -  - ]:          0 :             if (element.children.size() == 2)
     582                 :            :             {
     583                 :          0 :                 trace_out << "perform 2:8" << std::endl;
     584         [ -  - ]:          0 :                 refiner.derefine_two_to_one(tet_store,node_connectivity,i);
     585                 :            :             }
     586         [ -  - ]:          0 :             else if (element.children.size() == 4)
     587                 :            :             {
     588                 :          0 :                 trace_out << "perform 4:8" << std::endl;
     589         [ -  - ]:          0 :                 refiner.derefine_four_to_one(tet_store,node_connectivity,i);
     590                 :            :             }
     591                 :            :             else {
     592 [ -  - ][ -  - ]:          0 :                 std::cout << "num children " << element.children.size() << std::endl;
                 [ -  - ]
     593                 :          0 :                 assert(0);
     594                 :            :             }
     595                 :            :             // remove tets and edges marked for deletion above
     596         [ -  - ]:          0 :             refiner.delete_intermediates_of_children(tet_store);
     597         [ -  - ]:          0 :             tet_store.process_delete_list();
     598                 :            : 
     599         [ -  - ]:          0 :             refiner.refine_one_to_eight(tet_store,node_connectivity,i);
     600                 :            : 
     601                 :            :             // Grab children after it has been updated
     602 [ -  - ][ -  - ]:          0 :             auto current_children = tet_store.data(i).children;
     603                 :            : 
     604                 :            :             // I want to set the children stored in *my* own children, to be
     605                 :            :             // the value of my new children....
     606                 :            :             //refiner.overwrite_children(tet_store, former_children, current_children);
     607                 :            : 
     608         [ -  - ]:          0 :             tet_store.unset_marked_children(i); // FIXME: This will not work well in parallel
     609                 :          0 :             element.refinement_case = AMR::Refinement_Case::one_to_eight;
     610                 :          0 :         }
     611                 :            : 
     612                 :            :         // Clean up dead edges
     613                 :            :         // clean_up_dead_edges(); // Nothing get's marked as "dead" atm?
     614                 :            : 
     615                 :            :         //std::cout << "Total Edges : " << tet_store.edge_store.size() << std::endl;
     616                 :            :         //std::cout << "Total Tets : " << tet_store.size() << std::endl;
     617                 :            :         //std::cout << "Total Nodes : " << m_x.size() << std::endl;
     618                 :            : 
     619                 :         47 :         trace_out << "Done ref" << std::endl;
     620                 :         47 :         node_connectivity.print();
     621                 :         47 :         node_connectivity.print();
     622         [ +  - ]:         47 :         tet_store.print_node_types();
     623         [ +  - ]:         47 :         tet_store.print_tets();
     624                 :            :         //node_connectivity.print();
     625                 :            : 
     626                 :            :         //reset_intermediate_edges();
     627         [ +  - ]:         47 :         remove_edge_locks(1);
     628         [ +  - ]:         47 :         remove_normals();
     629                 :            : 
     630         [ +  - ]:         47 :         lock_intermediates();
     631                 :            : 
     632         [ +  + ]:     194921 :         for (auto& kv : tet_store.edge_store.edges) {
     633                 :     194874 :            auto& local = kv.second;
     634         [ +  + ]:     194874 :            if (local.needs_refining == 1) local.needs_refining = 0;
     635                 :            :         }
     636                 :         47 :     }
     637                 :            : 
     638                 :         71 :     void mesh_adapter_t::lock_intermediates()
     639                 :            :     {
     640                 :            :         /*
     641                 :            :         for (auto k : tet_store.intermediate_list)
     642                 :            :         {
     643                 :            :             refiner.lock_edges_from_node(tet_store,k, Edge_Lock_Case::intermediate);
     644                 :            :         }
     645                 :            :         */
     646                 :            :         // TODO: Passing tet_store twice probably isn't the best
     647                 :         71 :         refiner.lock_intermediates(tet_store, tet_store.intermediate_list, Edge_Lock_Case::intermediate);
     648                 :         71 :     }
     649                 :            : 
     650                 :            :     /**
     651                 :            :      * @brief A method implementing "Algorithm 1" from the paper
     652                 :            :      *
     653                 :            :      * @param num_to_refine Number of edges to refine
     654                 :            :      * @param tet_id The id of the given tet
     655                 :            :      */
     656                 :      50634 :     void mesh_adapter_t::refinement_class_one(int num_to_refine, size_t tet_id)
     657                 :            :     {
     658                 :      50634 :         trace_out << "Refinement Class One" << std::endl;
     659                 :            : 
     660                 :            :         // "If nrefine = 1
     661                 :            :         // Accept as a 1:2 refinement"
     662         [ +  + ]:      50634 :         if (num_to_refine == 1)
     663                 :            :         {
     664                 :        529 :             tet_store.mark_one_to_two(tet_id);
     665                 :            :         }
     666                 :            : 
     667                 :            :         // "Else if nrefine = 2 OR nrefine = 3"
     668 [ +  - ][ +  + ]:      50105 :         else if (num_to_refine > 1 && num_to_refine < 4)
     669                 :            :         {
     670                 :            : 
     671                 :            :             // We need to detect if the edges which need to refine are
     672                 :            :             // on the same face
     673                 :            :             // and if so which face so we know how to 1:4
     674                 :            : 
     675         [ +  - ]:        686 :             face_list_t face_list = tet_store.generate_face_lists(tet_id);
     676                 :        686 :             bool edges_on_same_face = false;
     677                 :        686 :             size_t face_refine_id = 0;
     678                 :            : 
     679                 :            :             // Iterate over each face
     680         [ +  + ]:       2099 :             for (size_t face = 0; face < NUM_TET_FACES; face++)
     681                 :            :             {
     682                 :       2018 :                 int num_face_refine_edges = 0;
     683                 :       2018 :                 face_ids_t face_ids = face_list[face];
     684                 :            : 
     685                 :       2018 :                 trace_out << "Face is " <<
     686                 :            :                     face_ids[0] << ", " <<
     687                 :            :                     face_ids[1] << ", " <<
     688                 :            :                     face_ids[2] << ", " <<
     689                 :            :                     std::endl;
     690                 :            : 
     691         [ +  - ]:       2018 :                 edge_list_t face_edge_list = AMR::edge_store_t::generate_keys_from_face_ids(face_ids);
     692                 :            :                 // For this face list, see which ones need refining
     693         [ +  + ]:       8072 :                 for (size_t k = 0; k < NUM_FACE_NODES; k++)
     694                 :            :                 {
     695                 :       6054 :                     edge_t key = face_edge_list[k];
     696 [ +  - ][ +  + ]:       6054 :                     if (tet_store.edge_store.get(key).needs_refining == 1)
     697                 :            :                     {
     698                 :       3346 :                         num_face_refine_edges++;
     699                 :            :                     }
     700                 :            :                 }
     701         [ +  + ]:       2018 :                 if (num_face_refine_edges == num_to_refine)
     702                 :            :                 {
     703                 :        605 :                     edges_on_same_face = true;
     704                 :        605 :                     face_refine_id = face;
     705                 :        605 :                     trace_out << "Breaking with face value " << face << std::endl;
     706                 :        605 :                     break;
     707                 :            :                 }
     708                 :            :             }
     709                 :            : 
     710                 :            :             // "If active edges are on the same face
     711                 :            :             // Activate any inactive edges of the face
     712                 :            :             // Accept as a 1:4 // refinement"
     713         [ +  + ]:        686 :             if (edges_on_same_face)
     714                 :            :             {
     715                 :        605 :                 size_t opposite_offset = AMR::node_connectivity_t::face_list_opposite(face_list,
     716                 :            :                         face_refine_id);
     717                 :            : 
     718         [ +  - ]:        605 :                 tet_t tet = tet_store.get(tet_id);
     719                 :        605 :                 size_t opposite_id = tet[opposite_offset];
     720                 :            : 
     721                 :        605 :                 trace_out << "face_refine_id " << face_refine_id << std::endl;
     722                 :        605 :                 trace_out << "opposite_offset " << opposite_offset << std::endl;
     723                 :        605 :                 trace_out << "opposite_id " << opposite_id << std::endl;
     724                 :            : 
     725                 :            :                 // Activate edges on this face
     726         [ +  - ]:        605 :                 edge_list_t face_edge_list = AMR::edge_store_t::generate_keys_from_face_ids(face_list[face_refine_id]);
     727                 :            : 
     728         [ +  + ]:       2420 :                 for (size_t k = 0; k < NUM_FACE_NODES; k++)
     729                 :            :                 {
     730                 :       1815 :                     edge_t key = face_edge_list[k];
     731         [ +  - ]:       1815 :                     tet_store.edge_store.mark_for_refinement(key);
     732                 :            :                 }
     733                 :            : 
     734                 :            :                 //refiner.refine_one_to_four(tet_id, face_list[face_refine_id],
     735                 :            :                 //opposite_id);
     736         [ +  - ]:        605 :                 tet_store.mark_one_to_four(tet_id);
     737                 :            :             }
     738                 :            :             // "Else if active edges are not on the same face
     739                 :            :             // Activate all edges
     740                 :            :             // Accept as a 1:8 refinement"
     741                 :            :             else {
     742                 :            :                 //refiner.refine_one_to_eight(tet_id);
     743         [ +  - ]:         81 :                 tet_store.mark_edges_for_refinement(tet_id);
     744         [ +  - ]:         81 :                 tet_store.mark_one_to_eight(tet_id);
     745                 :            :             }
     746                 :            : 
     747                 :        686 :         }
     748                 :            : 
     749                 :            :         // "Else if nrefine > 3
     750                 :            :         // Activate any inactive edges
     751                 :            :         // Accept as a 1:8 refinement"
     752         [ +  - ]:      49419 :         else if (num_to_refine > 3)
     753                 :            :         {
     754                 :            :             //refiner.refine_one_to_eight(tet_id);
     755                 :      49419 :             tet_store.mark_edges_for_refinement(tet_id);
     756                 :      49419 :             tet_store.mark_one_to_eight(tet_id);
     757                 :            :         }
     758                 :      50634 :     }
     759                 :            : 
     760                 :            :     // TODO: Document this
     761                 :          0 :     void mesh_adapter_t::lock_tet_edges(size_t tet_id) {
     762         [ -  - ]:          0 :         edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
     763         [ -  - ]:          0 :         for (size_t k = 0; k < NUM_TET_EDGES; k++)
     764                 :            :         {
     765                 :          0 :             edge_t key = edge_list[k];
     766 [ -  - ][ -  - ]:          0 :             if (tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::unlocked)
     767                 :            :             {
     768                 :          0 :                 trace_out << "LOCKING! " << key << std::endl;
     769         [ -  - ]:          0 :                 tet_store.edge_store.get(key).lock_case = AMR::Edge_Lock_Case::locked;
     770                 :            :             }
     771                 :            :         }
     772                 :          0 :     }
     773                 :            : 
     774                 :            :     // TODO: Document this
     775                 :            :     // TODO: This has too similar a name to deactivate_tet
     776                 :          0 :     void mesh_adapter_t::deactivate_tet_edges(size_t tet_id) {
     777         [ -  - ]:          0 :         edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
     778         [ -  - ]:          0 :         for (size_t k = 0; k < NUM_TET_EDGES; k++)
     779                 :            :         {
     780                 :          0 :             edge_t key = edge_list[k];
     781                 :            : 
     782         [ -  - ]:          0 :             tet_store.edge_store.unmark_for_refinement(key);
     783                 :          0 :             trace_out << "Deactivating " << key << std::endl;
     784         [ -  - ]:          0 :             tet_store.edge_store.get(key).needs_derefining = false;
     785                 :            :         }
     786                 :          0 :     }
     787                 :            : 
     788                 :            :     /**
     789                 :            :      * @brief Unmarks edges of given tet for derefinement only
     790                 :            :      *
     791                 :            :      * @param tet_id The id of the given tet
     792                 :            :      */
     793                 :      25898 :     void mesh_adapter_t::deactivate_deref_tet_edges(size_t tet_id) {
     794         [ +  - ]:      25898 :         edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
     795         [ +  + ]:     181286 :         for (size_t k = 0; k < NUM_TET_EDGES; k++)
     796                 :            :         {
     797                 :     155388 :             edge_t key = edge_list[k];
     798                 :            : 
     799                 :     155388 :             trace_out << "Deactivating " << key << std::endl;
     800         [ +  - ]:     155388 :             tet_store.edge_store.get(key).needs_derefining = false;
     801                 :            :         }
     802                 :      25898 :     }
     803                 :            : 
     804                 :            :     /**
     805                 :            :      * @brief An implementation of "Algorithm 2" from the paper
     806                 :            :      *
     807                 :            :      * @param edge_list The list of edges for the given tet
     808                 :            :      * @param tet_id The id of the given tet
     809                 :            :      */
     810                 :          0 :     void mesh_adapter_t::refinement_class_two(edge_list_t edge_list, size_t tet_id)
     811                 :            :     {
     812                 :          0 :         trace_out << "Refinement Class Two" << std::endl;
     813                 :            : 
     814                 :            : 
     815                 :            :         // "Deactivate all locked edges"
     816                 :            : 
     817                 :            :         // count number of active edges
     818                 :          0 :         int num_active_edges = 0;
     819         [ -  - ]:          0 :         for (size_t k = 0; k < NUM_TET_EDGES; k++)
     820                 :            :         {
     821                 :          0 :             edge_t key = edge_list[k];
     822 [ -  - ][ -  - ]:          0 :             if (tet_store.edge_store.get(key).lock_case != AMR::Edge_Lock_Case::unlocked)
     823                 :            :             {
     824         [ -  - ]:          0 :                 tet_store.edge_store.unmark_for_refinement(key);
     825                 :            :             }
     826                 :            :             // "Count number of active edges"
     827 [ -  - ][ -  - ]:          0 :             if (tet_store.edge_store.get(key).needs_refining == 1) {
     828                 :          0 :                 num_active_edges++;
     829                 :            :             }
     830                 :            :         }
     831                 :            : 
     832                 :            :         // Find out of two active edges live on the same face
     833                 :          0 :         bool face_refine = false;
     834                 :          0 :         size_t face_refine_id = 0; // FIXME: Does this need a better default
     835         [ -  - ]:          0 :         face_list_t face_list = tet_store.generate_face_lists(tet_id);
     836                 :            : 
     837                 :            :         // Iterate over each face
     838         [ -  - ]:          0 :         for (size_t face = 0; face < NUM_TET_FACES; face++)
     839                 :            :         {
     840                 :          0 :             trace_out << "face " << face << std::endl;
     841                 :          0 :             int num_face_refine_edges = 0;
     842                 :          0 :             int num_face_locked_edges = 0;
     843                 :            : 
     844                 :          0 :             face_ids_t face_ids = face_list[face];
     845         [ -  - ]:          0 :             edge_list_t face_edge_list = AMR::edge_store_t::generate_keys_from_face_ids(face_ids);
     846                 :            :             // For this face list, see which ones need refining
     847         [ -  - ]:          0 :             for (size_t k = 0; k < NUM_FACE_NODES; k++)
     848                 :            :             {
     849                 :          0 :                 edge_t key = face_edge_list[k];
     850                 :          0 :                 trace_out << "Checking " << key << std::endl;
     851 [ -  - ][ -  - ]:          0 :                 if (tet_store.edge_store.get(key).needs_refining == 1)
     852                 :            :                 {
     853                 :          0 :                     num_face_refine_edges++;
     854                 :          0 :                     trace_out << "ref! " << key << std::endl;
     855                 :            :                 }
     856                 :            : 
     857                 :            :                 // Check for locked edges
     858                 :            :                 // This case only cares about faces with no locks
     859 [ -  - ][ -  - ]:          0 :                 if (tet_store.edge_store.get(key).lock_case != AMR::Edge_Lock_Case::unlocked)
     860                 :            :                 {
     861                 :          0 :                     num_face_locked_edges++;
     862                 :          0 :                     trace_out << "locked! " << key << std::endl;
     863                 :            :                 }
     864                 :            :             }
     865                 :            : 
     866                 :            : 
     867                 :            :             // Decide if we want to process this face
     868 [ -  - ][ -  - ]:          0 :             if (num_face_refine_edges >= 2 && num_face_locked_edges == 0)
     869                 :            :             {
     870                 :            :                 // We can refine this face
     871                 :          0 :                 face_refine = true;
     872                 :          0 :                 face_refine_id = face;
     873                 :          0 :                 break;
     874                 :            :             }
     875                 :            :         }
     876                 :            : 
     877                 :            :         // "If nrefine = 1
     878                 :            :         // Accept as 1:2 refinement"
     879                 :            :         // TODO: can we hoist this higher
     880         [ -  - ]:          0 :         if (num_active_edges == 1)
     881                 :            :         {
     882                 :            :             //node_pair_t nodes = find_single_refinement_nodes(edge_list);
     883                 :            :             //refine_one_to_two( tet_id, nodes[0], nodes[1]);
     884         [ -  - ]:          0 :             tet_store.mark_one_to_two(tet_id);
     885                 :            :         }
     886                 :            :         // "Else if any face has nrefine >= 2 AND no locked edges
     887                 :            :         // Active any inactive edges of the face
     888                 :            :         // Accept as a 1:4 refinement"
     889         [ -  - ]:          0 :         else if (face_refine)
     890                 :            :         {
     891                 :          0 :             size_t opposite_offset = AMR::node_connectivity_t::face_list_opposite(face_list, face_refine_id);
     892                 :            : 
     893         [ -  - ]:          0 :             tet_t tet = tet_store.get(tet_id);
     894                 :          0 :             size_t opposite_id = tet[opposite_offset];
     895                 :            : 
     896                 :          0 :             trace_out << "Tet ID " << tet_id << std::endl;
     897                 :          0 :             trace_out << "Opposite offset " << opposite_offset << std::endl;
     898                 :          0 :             trace_out << "Opposite id " << opposite_id << std::endl;
     899                 :          0 :             trace_out << "Face refine id " << face_refine_id << std::endl;
     900                 :            : 
     901                 :            :             edge_list_t face_edge_list =
     902         [ -  - ]:          0 :                 AMR::edge_store_t::generate_keys_from_face_ids(face_list[face_refine_id]);
     903                 :            : 
     904         [ -  - ]:          0 :             for (size_t k = 0; k < NUM_FACE_NODES; k++)
     905                 :            :             {
     906                 :          0 :                 edge_t key = face_edge_list[k];
     907         [ -  - ]:          0 :                 tet_store.edge_store.mark_for_refinement(key);
     908                 :            :             }
     909                 :            : 
     910                 :            :             //refiner.refine_one_to_four(tet_id, face_list[face_refine_id], opposite_id);
     911         [ -  - ]:          0 :             tet_store.mark_one_to_four(tet_id);
     912                 :            :         }
     913                 :            : 
     914                 :            :         // "Else
     915                 :            :         // Deactivate all edges
     916                 :            :         // Mark all edges as locked"
     917                 :            :         else {
     918                 :          0 :             trace_out << "Class 2 causes some locking.." << std::endl;
     919         [ -  - ]:          0 :             deactivate_tet_edges(tet_id);
     920         [ -  - ]:          0 :             lock_tet_edges(tet_id);
     921                 :            :         }
     922                 :            : 
     923                 :          0 :     }
     924                 :            : 
     925                 :            :     /**
     926                 :            :      * @brief Based on a tet_id, decide if it's current state of locked
     927                 :            :      * and marked edges maps to a valid refinement case. The logic for
     928                 :            :      * this was derived from talking to JW and reading Chicoma.
     929                 :            :      *
     930                 :            :      * It basically just checks if something a 1:2 and has 3
     931                 :            :      * intermediates and 3 makred edges, or is a 1:4 and has 5/6
     932                 :            :      * intermediates
     933                 :            :      *
     934                 :            :      * @param child_id the id of the tet to check
     935                 :            :      *
     936                 :            :      * @return A bool saying if the tet is in a valid state to be refined
     937                 :            :      */
     938                 :          0 :     bool mesh_adapter_t::check_valid_refinement_case(size_t child_id) {
     939                 :            : 
     940                 :          0 :         trace_out << "check valid ref " << child_id << std::endl;
     941         [ -  - ]:          0 :         edge_list_t edge_list = tet_store.generate_edge_keys(child_id);
     942                 :            : 
     943                 :          0 :         size_t num_to_refine = 0;
     944                 :          0 :         size_t num_intermediate = 0;
     945                 :          0 :         size_t unlocked = 0;
     946                 :          0 :         size_t locked = 0;
     947                 :            : 
     948         [ -  - ]:          0 :         for (size_t k = 0; k < NUM_TET_EDGES; k++)
     949                 :            :         {
     950                 :          0 :             edge_t key = edge_list[k];
     951                 :          0 :             trace_out << "Key " << key << std::endl;
     952                 :            : 
     953                 :            :             // Count intermediate edges
     954         [ -  - ]:          0 :             if (tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::intermediate
     955 [ -  - ][ -  - ]:          0 :               || tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::temporary)
         [ -  - ][ -  - ]
     956                 :            :             {
     957                 :          0 :                 trace_out << "found intermediate" << std::endl;
     958                 :          0 :                 num_intermediate++;
     959                 :            :             }
     960                 :            : 
     961                 :            :             // Count number of marked for refinement
     962 [ -  - ][ -  - ]:          0 :             if (tet_store.edge_store.get(key).needs_refining == 1)
     963                 :            :             {
     964                 :          0 :                 trace_out << "found refine" << std::endl;
     965                 :          0 :                 num_to_refine++;
     966                 :            :             }
     967                 :            : 
     968 [ -  - ][ -  - ]:          0 :             if (tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::unlocked)
     969                 :            :             {
     970                 :          0 :                 trace_out << "found unlocked" << std::endl;
     971                 :          0 :                 unlocked++;
     972                 :            :             }
     973                 :            : 
     974 [ -  - ][ -  - ]:          0 :             if (tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::locked)
     975                 :            :             {
     976                 :          0 :                 trace_out << "found locked" << std::endl;
     977                 :          0 :                 locked++;
     978                 :            :             }
     979                 :            : 
     980                 :            :         }
     981                 :            : 
     982         [ -  - ]:          0 :         AMR::Refinement_State& element = tet_store.data(child_id);
     983                 :            : 
     984                 :          0 :         trace_out <<
     985                 :            :             "Intermediates " << num_intermediate <<
     986                 :            :             " num to refine " << num_to_refine <<
     987                 :            :             " unlocked " << unlocked <<
     988                 :            :             " locked " << locked <<
     989                 :            :             " Case " << element.refinement_case <<
     990                 :            :             std::endl;
     991                 :            : 
     992                 :            :         // check if element is 1:2
     993         [ -  - ]:          0 :         if (element.refinement_case == AMR::Refinement_Case::one_to_two)
     994                 :            :         {
     995                 :            :             // If so check it has 3 intermediates and 3 which need refining
     996 [ -  - ][ -  - ]:          0 :             if (num_intermediate != 3 || num_to_refine != 3) {
     997                 :          0 :                 return false;
     998                 :            :             }
     999                 :            :             else {
    1000                 :          0 :                 trace_out << "True " <<
    1001                 :            :                     "Intermediates " << num_intermediate <<
    1002                 :            :                     " num to refine " << num_to_refine <<
    1003                 :            :                     " Case " << element.refinement_case <<
    1004                 :            :                     " 2:1 " << AMR::Refinement_Case::one_to_two <<
    1005                 :            :                     std::endl;
    1006                 :            :             }
    1007                 :            :         }
    1008                 :            : 
    1009                 :            :         // check if element is 1:4
    1010         [ -  - ]:          0 :         else if (element.refinement_case == AMR::Refinement_Case::one_to_four)
    1011                 :            :         {
    1012                 :            :             // TODO: Check if it's a center tet for a 1:4
    1013                 :            :             // FIXME: Is this even needed? How else would you get these
    1014                 :            :             // combinations? Can't we just combine these two checks?
    1015                 :            : 
    1016         [ -  - ]:          0 :             bool is_center_tet = tet_store.is_center(child_id);
    1017                 :            : 
    1018         [ -  - ]:          0 :             if (is_center_tet)
    1019                 :            :             {
    1020 [ -  - ][ -  - ]:          0 :                 if (num_to_refine != 0 || num_intermediate != 6)
    1021                 :            :                 {
    1022                 :          0 :                     trace_out << "Fail compat 1:4 center" << std::endl;
    1023                 :          0 :                     return false;
    1024                 :            :                 }
    1025                 :            :             }
    1026                 :            :             else { // Is one of the outsides (not center)
    1027 [ -  - ][ -  - ]:          0 :                 if (num_to_refine != 1 || num_intermediate != 5)
    1028                 :            :                 {
    1029                 :          0 :                     trace_out << "Fail compat 1:4 non center" << std::endl;
    1030                 :          0 :                     return false;
    1031                 :            :                 }
    1032                 :            :             }
    1033                 :            :         }
    1034                 :            : 
    1035                 :            :         // If it makes it here, it's compatible
    1036                 :          0 :         return true;
    1037                 :            : 
    1038                 :            :     }
    1039                 :            : 
    1040                 :            :     /**
    1041                 :            :      * @brief Place holder method for the implementation of "Algorithm
    1042                 :            :      * 3" from the paper
    1043                 :            :      */
    1044                 :            :     // TODO: Does this parse a childs siblings multiple times?
    1045                 :          0 :     void mesh_adapter_t::refinement_class_three(size_t tet_id) {
    1046                 :            : 
    1047                 :          0 :         trace_out << "Refinement Class Three" << std::endl;
    1048                 :            : 
    1049                 :            :         // "Identify parent element iparent"
    1050                 :            :         // TODO: WE should either always use the id to fetch, or always do the data lookup
    1051                 :            :         //size_t parent_id = master_elements.get_parent(tet_id);
    1052         [ -  - ]:          0 :         size_t parent_id = tet_store.get_parent_id(tet_id);
    1053                 :            : 
    1054                 :          0 :         trace_out << "Parent id = " << parent_id << std::endl;
    1055                 :            : 
    1056                 :            :         // NOTE: This implies comms when we use these ids?
    1057 [ -  - ][ -  - ]:          0 :         child_id_list_t children = tet_store.data(parent_id).children;
    1058                 :            : 
    1059                 :            :         // "Do for each child element ielement
    1060                 :            :         // Activate all non-locked edges
    1061                 :            :         // Deactivate all locked edges"
    1062         [ -  - ]:          0 :         for (size_t i = 0; i < children.size(); i++)
    1063                 :            :         {
    1064                 :            :             // TODO: Is this in element or tet ids?
    1065                 :          0 :             trace_out << "Checking child " << children[i] << std::endl;
    1066         [ -  - ]:          0 :             edge_list_t edge_list = tet_store.generate_edge_keys(children[i]);
    1067         [ -  - ]:          0 :             for (size_t k = 0; k < NUM_TET_EDGES; k++)
    1068                 :            :             {
    1069                 :          0 :                 edge_t key = edge_list[k];
    1070                 :          0 :                 trace_out << "Compat 3 " << key << std::endl;
    1071 [ -  - ][ -  - ]:          0 :                 if (tet_store.edge_store.get(key).lock_case == AMR::Edge_Lock_Case::unlocked)
    1072                 :            :                 {
    1073                 :          0 :                     trace_out << "Compat 3 marking edge " << key << std::endl;
    1074         [ -  - ]:          0 :                     tet_store.edge_store.mark_for_refinement(key);
    1075                 :            :                 }
    1076                 :            :                 else {
    1077         [ -  - ]:          0 :                     tet_store.edge_store.unmark_for_refinement(key);
    1078                 :            :                 }
    1079                 :            :             }
    1080                 :            :         }
    1081                 :            : 
    1082                 :            :         // "Set compatible = TRUE
    1083                 :          0 :         bool compatible = true;
    1084                 :            : 
    1085                 :            :         // Do for each child element ielement
    1086                 :            :         // If ielement is not a valid refinement case
    1087                 :            :         // compatible = FALSE"
    1088         [ -  - ]:          0 :         for (size_t i = 0; i < children.size(); i++)
    1089                 :            :         {
    1090                 :          0 :             size_t child = children[i];
    1091 [ -  - ][ -  - ]:          0 :             if ( !check_valid_refinement_case(child) )
    1092                 :            :             {
    1093                 :          0 :                 trace_out << "Compat 3 Marking compatible false because of invalid refinement case" << std::endl;
    1094                 :            : 
    1095                 :          0 :                 compatible = false;
    1096                 :            :             }
    1097                 :            :             else {
    1098                 :          0 :                 trace_out << "Is compatible" << std::endl;
    1099                 :            :             }
    1100                 :            :         }
    1101                 :            : 
    1102                 :            :         // "If compatible = FALSE
    1103                 :            :         // Do for each child element ielement
    1104                 :            :         // Deactive all edges of ielement
    1105                 :            :         // Mark all edges of ielement as locked
    1106                 :            :         // Mark ielement as normal"
    1107         [ -  - ]:          0 :         if (compatible == false)
    1108                 :            :         {
    1109         [ -  - ]:          0 :             for (size_t i = 0; i < children.size(); i++)
    1110                 :            :             {
    1111                 :          0 :                 size_t child = children[i];
    1112         [ -  - ]:          0 :                 deactivate_tet_edges(child);
    1113         [ -  - ]:          0 :                 lock_tet_edges(child);
    1114                 :          0 :                 trace_out << "Compat 3 locking edges of " << child << std::endl;
    1115                 :            :                 // Here we interpret normal to mean "don't treat it like it has intermediates"
    1116         [ -  - ]:          0 :                 tet_store.mark_normal(child);
    1117                 :          0 :                 trace_out << "Compat 3 " << child << std::endl;
    1118                 :            :             }
    1119                 :            :         }
    1120                 :            :         else {
    1121                 :          0 :             trace_out << "TIME TO 2:8 " << tet_id << std::endl;
    1122                 :            :             // Accept as 2:8 or 4:8
    1123         [ -  - ]:          0 :             AMR::Refinement_State& element = tet_store.data(tet_id);
    1124         [ -  - ]:          0 :             if (element.refinement_case == AMR::Refinement_Case::one_to_two)
    1125                 :            :             {
    1126         [ -  - ]:          0 :                 tet_store.mark_two_to_eight(tet_id);
    1127                 :            :             }
    1128         [ -  - ]:          0 :             else if (element.refinement_case == AMR::Refinement_Case::one_to_four)
    1129                 :            :             {
    1130         [ -  - ]:          0 :                 tet_store.mark_four_to_eight(tet_id);
    1131                 :            :             }
    1132                 :            :             else {
    1133                 :          0 :                 trace_out << " I don't know what to do with this..it looks like you're trying to 2/4:8 an 8... " << std::endl;
    1134                 :            :             }
    1135                 :            : 
    1136                 :            :         }
    1137                 :          0 :     }
    1138                 :            : 
    1139                 :         47 :     void mesh_adapter_t::remove_normals()
    1140                 :            :     {
    1141         [ +  + ]:     137891 :         for (const auto& kv : tet_store.tets)
    1142                 :            :         {
    1143                 :     137844 :             size_t tet_id = kv.first;
    1144         [ +  - ]:     137844 :             tet_store.set_normal(tet_id, 0);
    1145                 :            :         }
    1146                 :         47 :     }
    1147                 :            : 
    1148                 :        118 :     void mesh_adapter_t::remove_edge_locks(int intermediate)
    1149                 :            :     {
    1150         [ +  + ]:     283420 :         for (const auto& kv : tet_store.tets)
    1151                 :            :         {
    1152                 :     283302 :             size_t tet_id = kv.first;
    1153                 :            : 
    1154                 :     283302 :             trace_out << "Process tet removelock " << tet_id << std::endl;
    1155                 :            : 
    1156                 :            :             // Only apply checks to tets on the active list
    1157 [ +  - ][ +  + ]:     283302 :             if (tet_store.is_active(tet_id)) {
    1158                 :            :                 // change it from intermediate to locked
    1159         [ +  - ]:     249425 :                 update_tet_edges_lock_type(tet_id, AMR::Edge_Lock_Case::locked, AMR::Edge_Lock_Case::unlocked);
    1160         [ +  - ]:     249425 :                 if (intermediate) {
    1161         [ +  - ]:     249425 :                     update_tet_edges_lock_type(tet_id, AMR::Edge_Lock_Case::intermediate, AMR::Edge_Lock_Case::unlocked);
    1162                 :            :                 }
    1163                 :            :             }
    1164                 :            :         }
    1165                 :            : 
    1166                 :        118 :     }
    1167                 :            :     //void mesh_adapter_t::reset_intermediate_edges()
    1168                 :            :     //{
    1169                 :            :     //    for (const auto& kv : tet_store.tets)
    1170                 :            :     //    {
    1171                 :            :     //        size_t tet_id = kv.first;
    1172                 :            : 
    1173                 :            :     //        trace_out << "Process tet reset " << tet_id << std::endl;
    1174                 :            : 
    1175                 :            :     //        // Only apply checks to tets on the active list
    1176                 :            :     //        if (tet_store.is_active(tet_id)) {
    1177                 :            :     //            // change it from intermediate to locked
    1178                 :            :     //            update_tet_edges_lock_type(tet_id, AMR::Edge_Lock_Case::intermediate, AMR::Edge_Lock_Case::locked);
    1179                 :            :     //        }
    1180                 :            :     //    }
    1181                 :            :     //}
    1182                 :            : 
    1183                 :     627017 :     void mesh_adapter_t::update_tet_edges_lock_type(size_t tet_id, AMR::Edge_Lock_Case check, AMR::Edge_Lock_Case new_case) {
    1184         [ +  - ]:     627017 :         edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
    1185         [ +  + ]:    4389119 :         for (size_t k = 0; k < NUM_TET_EDGES; k++)
    1186                 :            :         {
    1187                 :    3762102 :             edge_t key = edge_list[k];
    1188 [ +  - ][ +  + ]:    3762102 :             if (tet_store.edge_store.get(key).lock_case == check)
    1189                 :            :             {
    1190         [ +  - ]:       3796 :                 tet_store.edge_store.get(key).lock_case = new_case;
    1191                 :            :             }
    1192                 :            :         }
    1193                 :     627017 :     }
    1194                 :            : 
    1195                 :            :     /**
    1196                 :            :      * @brief This unlocks edges that were previously locked with a `temporary'
    1197                 :            :      * lock, indicating a parallel compatibility induced locking
    1198                 :            :      */
    1199                 :         71 :     void mesh_adapter_t::remove_edge_temp_locks()
    1200                 :            :     {
    1201         [ +  + ]:     145529 :       for (const auto& kv : tet_store.tets)
    1202                 :            :       {
    1203                 :     145458 :         size_t tet_id = kv.first;
    1204                 :            : 
    1205                 :     145458 :         trace_out << "Process tet remove temp lock " << tet_id << std::endl;
    1206                 :            : 
    1207                 :            :         // Only apply checks to tets on the active list
    1208 [ +  - ][ +  + ]:     145458 :         if (tet_store.is_active(tet_id)) {
    1209                 :            :           // change it from temporary to unlocked
    1210         [ +  - ]:     128167 :           update_tet_edges_lock_type(tet_id, AMR::Edge_Lock_Case::temporary,
    1211                 :            :             AMR::Edge_Lock_Case::unlocked);
    1212                 :            :         }
    1213                 :            :       }
    1214                 :         71 :     }
    1215                 :            : 
    1216                 :         40 :     void mesh_adapter_t::mark_derefinement()
    1217                 :            :     {
    1218                 :         40 :         const size_t max_num_rounds = AMR_MAX_ROUNDS;
    1219                 :            : 
    1220                 :            :         // Mark refinements
    1221                 :            :         size_t iter;
    1222                 :            :         //Iterate until convergence
    1223         [ +  - ]:         65 :         for (iter = 0; iter < max_num_rounds; iter++)
    1224                 :            :         {
    1225                 :         65 :             tet_store.marked_derefinements.get_state_changed() = false;
    1226                 :            : 
    1227                 :            :             // set of elements which have been considered for derefinement
    1228                 :         65 :             std::unordered_set< size_t > done_deref_marking;
    1229                 :            : 
    1230                 :            :             // Loop over tets
    1231         [ +  + ]:     132964 :             for (const auto& kv : tet_store.tets)
    1232                 :            :             {
    1233                 :            :                 // this loop only runs for active tets
    1234 [ +  - ][ +  + ]:     132899 :                 if (!tet_store.is_active(kv.first)) {
    1235         [ +  - ]:      16442 :                   deactivate_deref_tet_edges(kv.first);
    1236                 :     118687 :                   continue;
    1237                 :            :                 }
    1238                 :     116457 :                 size_t activetet_id = kv.first;
    1239                 :            : 
    1240                 :            :                 // check if activetet_id has a parent (assign to tet_id)
    1241                 :            :                 // if it does not, activetet_id is not a derefinement candidate
    1242                 :            :                 size_t tet_id;
    1243         [ +  - ]:     116457 :                 const auto& activetet_data = tet_store.data(activetet_id);
    1244         [ +  + ]:     116457 :                 if (!activetet_data.has_parent) {
    1245         [ +  - ]:       1551 :                   deactivate_deref_tet_edges(activetet_id);
    1246                 :       1551 :                   continue;
    1247                 :            :                 }
    1248                 :            :                 else {
    1249                 :     114906 :                   tet_id = activetet_data.parent_id;
    1250                 :            :                 }
    1251                 :            : 
    1252                 :            :                 // if already considered for deref, do not reconsider
    1253 [ +  - ][ +  + ]:     114906 :                 if (done_deref_marking.count(tet_id) > 0) {
    1254                 :      99949 :                   continue;
    1255                 :            :                 }
    1256         [ +  - ]:      14957 :                 done_deref_marking.insert(tet_id);
    1257                 :            : 
    1258 [ +  - ][ +  - ]:      14957 :                 child_id_list_t children = tet_store.data(tet_id).children;
    1259                 :            : 
    1260                 :            :                 // check if any child of tet_id (i.e. any active tet) is marked
    1261                 :            :                 // for refinement
    1262                 :      14957 :                 bool is_child_ref(false);
    1263         [ +  + ]:     130337 :                 for (size_t i=0; i<children.size(); i++) {
    1264         [ +  - ]:     115380 :                   edge_list_t chedge_list = tet_store.generate_edge_keys(children[i]);
    1265                 :            :                   // Check each edge, see if it is marked for refinement
    1266         [ +  + ]:     807660 :                   for (size_t k=0; k<NUM_TET_EDGES; k++) {
    1267                 :     692280 :                     edge_t edge = chedge_list[k];
    1268                 :            : 
    1269 [ +  - ][ +  + ]:     692280 :                     if (tet_store.edge_store.get(edge).needs_refining == 1) {
    1270                 :      30120 :                       is_child_ref = true;
    1271                 :      30120 :                       continue;
    1272                 :            :                     }
    1273                 :            :                   }
    1274                 :            :                 }
    1275                 :            :                 // deactivate from deref if marked for ref
    1276         [ +  + ]:      14957 :                 if (is_child_ref) {
    1277                 :        745 :                   trace_out << tet_id << " Looping cancelled since child marked for refinement." << std::endl;
    1278         [ +  + ]:       6705 :                   for (auto child_id : children) {
    1279         [ +  - ]:       5960 :                     deactivate_deref_tet_edges(child_id);
    1280                 :            :                   }
    1281         [ +  - ]:        745 :                   deactivate_deref_tet_edges(tet_id);
    1282         [ +  - ]:        745 :                   tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::skip);
    1283                 :        745 :                   continue;
    1284                 :        745 :                 }
    1285                 :            : 
    1286                 :            :                 // check if tet_id has been marked for deref-ref
    1287         [ +  - ]:      14212 :                 edge_list_t pedge_list = tet_store.generate_edge_keys(tet_id);
    1288                 :            :                 // Check each edge, see if it is marked for refinement
    1289         [ +  + ]:      99484 :                 for (size_t k=0; k<NUM_TET_EDGES; k++) {
    1290                 :      85272 :                   edge_t edge = pedge_list[k];
    1291                 :            : 
    1292                 :            :                   // deactivate child-edges from deref if marked '2'
    1293 [ +  - ][ +  + ]:      85272 :                   if (tet_store.edge_store.get(edge).needs_refining == 2) {
    1294                 :        236 :                     auto edge_nodes = edge.get_data();
    1295         [ +  - ]:        236 :                     auto ch_node = node_connectivity.data().at(edge_nodes);
    1296                 :        236 :                     std::array< edge_t, 2> ch_edge;
    1297         [ +  - ]:        236 :                     ch_edge[0] = {edge_nodes[0], ch_node};
    1298         [ +  - ]:        236 :                     ch_edge[1] = {edge_nodes[1], ch_node};
    1299         [ +  - ]:        236 :                     tet_store.edge_store.get(ch_edge[0]).needs_derefining = 0;
    1300         [ +  - ]:        236 :                     tet_store.edge_store.get(ch_edge[1]).needs_derefining = 0;
    1301                 :            :                   }
    1302                 :            :                 }
    1303                 :            : 
    1304                 :            :                 // This is useful for later inspection
    1305                 :            :                 //edge_list_t edge_list = tet_store.generate_edge_keys(tet_id);
    1306                 :      14212 :                 std::size_t num_to_derefine = 0; // Nodes
    1307                 :            : 
    1308         [ +  - ]:      14212 :                 AMR::Refinement_Case refinement_case = tet_store.get_refinement_case(tet_id);
    1309                 :            : 
    1310         [ +  - ]:      14212 :                 auto derefine_node_set = refiner.find_derefine_node_set(tet_store, tet_id);
    1311                 :            :                 // Find the set of nodes which are not in the parent
    1312                 :            :                 std::unordered_set<size_t> non_parent_nodes =
    1313         [ +  - ]:      14212 :                   refiner.child_exclusive_nodes(tet_store, tet_id);
    1314                 :            : 
    1315                 :            :                 //for (auto drnode: derefine_node_set)
    1316                 :            :                 //  trace_out << "derefine node: " << drnode << std::endl;
    1317                 :            : 
    1318                 :      14212 :                 num_to_derefine = derefine_node_set.size();
    1319                 :            : 
    1320         [ +  + ]:      14212 :                 if (num_to_derefine > 0) {
    1321                 :      14094 :                     trace_out << "num_to_derefine " << num_to_derefine << std::endl;
    1322                 :      14094 :                     trace_out << "ref_case " << refinement_case << std::endl;
    1323                 :      14094 :                     trace_out << "num children " << children.size() << std::endl;
    1324                 :            :                 }
    1325                 :            : 
    1326                 :            :                 //num_to_derefine = convert_derefine_edges_to_points(tet_store, tet_id, num_edges_to_derefine, refinement_case);
    1327                 :            : 
    1328                 :            :                 // "If nderefine = 1
    1329         [ +  + ]:      14212 :                 if (num_to_derefine == 1)
    1330                 :            :                 {
    1331                 :            :                     // If icase = 1:2
    1332                 :            : 
    1333                 :            :                     //if (refinement_case == AMR::Refinement_Case::one_to_two)
    1334         [ +  + ]:        390 :                     if (children.size() == 2)
    1335                 :            :                     {
    1336                 :            :                         // Accept as 2:1 derefine"
    1337                 :        378 :                         trace_out << "Accept as 2:1" << std::endl;
    1338                 :            :                         //refiner.derefine_two_to_one(tet_store, node_connectivity, tet_id);
    1339         [ +  - ]:        378 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::two_to_one);
    1340                 :            :                     }
    1341                 :            :                     // "Else
    1342                 :            :                     else {
    1343                 :            :                         // Deactivate all points"
    1344         [ +  + ]:        108 :                         for (auto child_id : children) {
    1345         [ +  - ]:         96 :                           deactivate_deref_tet_edges(child_id);
    1346                 :            :                         }
    1347         [ +  - ]:         12 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::skip);
    1348                 :         12 :                         trace_out << "giving up on deref decision. deactivate near 2:1 ntd = 1" << std::endl;
    1349                 :            :                     }
    1350                 :            :                 }
    1351                 :            : 
    1352                 :            : 
    1353                 :            :                 // "If nderefine = 2
    1354         [ +  + ]:      13822 :                 else if (num_to_derefine == 2)
    1355                 :            :                 {
    1356                 :            :                     // If icase = 1:4
    1357                 :            :                     //if (refinement_case == AMR::Refinement_Case::one_to_four)
    1358         [ -  + ]:         14 :                     if (children.size() == 4)
    1359                 :            :                     {
    1360                 :            :                         // Accept as 4:2 derefine"
    1361                 :          0 :                         trace_out << "Accept as 4:2" << std::endl;
    1362                 :            :                         //refiner.derefine_four_to_two(tet_store,  node_connectivity, tet_id);
    1363         [ -  - ]:          0 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::four_to_two);
    1364                 :            :                     }
    1365                 :            :                     // "Else
    1366                 :            :                     else {
    1367                 :            :                         // Deactivate all points"
    1368         [ +  + ]:        126 :                         for (auto child_id : children) {
    1369         [ +  - ]:        112 :                           deactivate_deref_tet_edges(child_id);
    1370                 :            :                         }
    1371         [ +  - ]:         14 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::skip);
    1372                 :         14 :                         trace_out << "giving up on deref decision. deactivate near 4:2 ntd = 2" << std::endl;
    1373                 :            :                     }
    1374                 :            :                 }
    1375                 :            : 
    1376                 :            :                 // "If nderefine = 3
    1377         [ +  + ]:      13808 :                 else if (num_to_derefine == 3)
    1378                 :            :                 {
    1379                 :            :                     // If icase = 1:4
    1380                 :            :                     //if (refinement_case == AMR::Refinement_Case::one_to_four)
    1381         [ +  + ]:        522 :                     if (children.size() == 4)
    1382                 :            :                     {
    1383                 :            :                         // Accept as 4:1 derefine"
    1384                 :        502 :                         trace_out << "Accept as 4:1" << std::endl;
    1385                 :            :                         //refiner.derefine_four_to_one(tet_store,  node_connectivity, tet_id);
    1386         [ +  - ]:        502 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::four_to_one);
    1387                 :            :                     }
    1388                 :            :                     // "Else if icase = 1:8
    1389                 :            :                     //else if (refinement_case == AMR::Refinement_Case::one_to_eight)
    1390         [ +  - ]:         20 :                     else if (children.size() == 8)
    1391                 :            :                     {
    1392                 :            :                         // we have a list of (non-parent) nodes that is marked
    1393                 :            :                         // for derefinement. First, determine the nodes that are
    1394                 :            :                         // unmarked for derefinement (or inactive_nodes). Then,
    1395                 :            :                         // determine if these are on a single face.
    1396                 :         20 :                         std::unordered_set<size_t> inactive_node_set;
    1397         [ +  + ]:        140 :                         for (auto npn : non_parent_nodes) {
    1398 [ +  - ][ +  + ]:        120 :                           if (derefine_node_set.count(npn) == 0)
    1399         [ +  - ]:         60 :                             inactive_node_set.insert(npn);
    1400                 :            :                         }
    1401 [ -  + ][ -  - ]:         20 :                         Assert(inactive_node_set.size() == 3, "Incorrectly "
         [ -  - ][ -  - ]
    1402                 :            :                           "sized inactive-node set");
    1403         [ +  - ]:         20 :                         auto same_face = check_same_face(tet_id, inactive_node_set);
    1404                 :            :                         // If inactive points lie on same face
    1405         [ +  + ]:         20 :                         if (same_face.first == true)
    1406                 :            :                         {
    1407                 :            :                             // Accept as 8:4 derefinement
    1408                 :         14 :                             trace_out << "Accept as 8:4" << std::endl;
    1409                 :            : 
    1410                 :            :                             // create a vector of node-array-pairs to mark edges
    1411                 :            :                             // for refinement 1:4
    1412                 :         14 :                             std::vector< std::array< std::size_t, 2 > > ref_edges;
    1413                 :         14 :                             trace_out << "inactive nodes on same face: ";
    1414         [ +  + ]:         56 :                             for (auto n:inactive_node_set) {
    1415                 :         42 :                               trace_out << n << ", ";
    1416 [ +  - ][ +  - ]:         42 :                               ref_edges.push_back(node_connectivity.get(n));
    1417                 :            :                             }
    1418                 :         14 :                             trace_out << std::endl;
    1419                 :            : 
    1420 [ +  - ][ +  - ]:         14 :                             tet_store.edge_store.mark_edges_for_deref_ref(ref_edges);
    1421                 :            :                             //refiner.derefine_eight_to_four(tet_store,  node_connectivity, tet_id);
    1422         [ +  - ]:         14 :                             tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::eight_to_four);
    1423                 :         14 :                         }
    1424                 :            :                         // "Else
    1425                 :            :                         else {
    1426                 :            :                             // Deactivate all points"
    1427         [ +  + ]:         54 :                             for (auto child_id : children) {
    1428         [ +  - ]:         48 :                               deactivate_deref_tet_edges(child_id);
    1429                 :            :                             }
    1430         [ +  - ]:          6 :                             tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::skip);
    1431                 :          6 :                             trace_out << "giving up on deref decision. deactivate near 8:4 ntd = 3" << std::endl;
    1432                 :            :                         }
    1433                 :            : 
    1434                 :         20 :                     }
    1435                 :            :                 }
    1436                 :            : 
    1437                 :            :                 // "If nderefine = 4
    1438         [ +  + ]:      13286 :                 else if (num_to_derefine == 4)
    1439                 :            :                 //else if (children.size() == 4)
    1440                 :            :                 {
    1441                 :            :                     // we have a list of (non-parent) nodes that is marked
    1442                 :            :                     // for derefinement. First, determine the nodes that are
    1443                 :            :                     // unmarked for derefinement (or inactive_nodes). Then,
    1444                 :            :                     // determine if these are on a single face.
    1445                 :          2 :                     std::unordered_set<size_t> inactive_node_set;
    1446         [ +  + ]:         14 :                     for (auto npn : non_parent_nodes) {
    1447 [ +  - ][ +  + ]:         12 :                       if (derefine_node_set.count(npn) == 0)
    1448         [ +  - ]:          4 :                         inactive_node_set.insert(npn);
    1449                 :            :                     }
    1450 [ -  + ][ -  - ]:          2 :                     Assert(inactive_node_set.size() == 2, "Incorrectly "
         [ -  - ][ -  - ]
    1451                 :            :                       "sized inactive-node set");
    1452                 :            :                     // Check if the inactive point belong to the same parent
    1453                 :            :                     // face and deactivate the third point on that face
    1454         [ +  - ]:          2 :                     auto same_face = check_same_face(tet_id, inactive_node_set);
    1455                 :            : 
    1456         [ +  - ]:          2 :                     if (same_face.first == true)
    1457                 :            :                     {
    1458                 :            :                         // deactivate the edges associated with same_face.second
    1459         [ +  + ]:         18 :                         for (size_t i = 0; i < children.size(); i++)
    1460                 :            :                         {
    1461         [ +  - ]:         16 :                           edge_list_t edge_list = tet_store.generate_edge_keys(children[i]);
    1462         [ +  + ]:        112 :                           for (size_t k = 0; k < NUM_TET_EDGES; k++)
    1463                 :            :                           {
    1464                 :         96 :                             edge_t edge = edge_list[k];
    1465                 :         96 :                             size_t A = edge.first();
    1466                 :         96 :                             size_t B = edge.second();
    1467 [ +  + ][ +  + ]:         96 :                             if (A == same_face.second || B == same_face.second)
    1468         [ +  - ]:         36 :                               tet_store.edge_store.get(edge).needs_derefining = false;
    1469                 :            :                           }
    1470                 :            :                         }
    1471                 :            : 
    1472                 :            :                         // create a vector of node-array-pairs to mark edges
    1473                 :            :                         // for refinement 1:4
    1474         [ +  - ]:          2 :                         inactive_node_set.insert(same_face.second);
    1475                 :          2 :                         std::vector< std::array< std::size_t, 2 > > ref_edges;
    1476         [ +  + ]:          8 :                         for (auto n:inactive_node_set) {
    1477 [ +  - ][ +  - ]:          6 :                           ref_edges.push_back(node_connectivity.get(n));
    1478                 :            :                         }
    1479                 :            : 
    1480 [ +  - ][ +  - ]:          2 :                         tet_store.edge_store.mark_edges_for_deref_ref(ref_edges);
    1481                 :            : 
    1482                 :            :                         // Accept as 8:4 derefinement
    1483                 :          2 :                         trace_out << "Accept as 8:4" << std::endl;
    1484                 :            :                         //refiner.derefine_eight_to_four(tet_store,  node_connectivity, tet_id);
    1485         [ +  - ]:          2 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::eight_to_four);
    1486                 :          2 :                     }
    1487                 :            :                     // "Else
    1488                 :            :                     else {
    1489                 :            :                         // Deactivate all points"
    1490         [ -  - ]:          0 :                         for (auto child_id : children) {
    1491         [ -  - ]:          0 :                           deactivate_deref_tet_edges(child_id);
    1492                 :            :                         }
    1493         [ -  - ]:          0 :                         tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::skip);
    1494                 :          0 :                         trace_out << "giving up on deref decision. deactivate near 8:4 ntd = 4" << std::endl;
    1495                 :            :                     }
    1496                 :          2 :                 }
    1497                 :            : 
    1498                 :            :                 // "If nderefine = 5
    1499         [ -  + ]:      13284 :                 else if (num_to_derefine == 5)
    1500                 :            :                 {
    1501                 :            :                     // Accept as 8:2 derefine"
    1502                 :          0 :                     trace_out << "Accept as 8:2 " << std::endl;
    1503                 :            :                     //refiner.derefine_eight_to_two(tet_store,  node_connectivity, tet_id);
    1504         [ -  - ]:          0 :                     tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::eight_to_two);
    1505                 :            :                 }
    1506                 :            : 
    1507                 :            :                 // "If nderefine = 6
    1508         [ +  + ]:      13284 :                 else if (num_to_derefine == 6)
    1509                 :            :                 {
    1510                 :            :                     // Accept as 8:1 derefine"
    1511                 :      13166 :                     trace_out << "Accept as 8:1" << std::endl;
    1512                 :            :                     //refiner.derefine_eight_to_one(tet_store,  node_connectivity, tet_id);
    1513         [ +  - ]:      13166 :                     tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::eight_to_one);
    1514                 :            :                 }
    1515                 :            : 
    1516                 :            :                 // "If nderefine = 0
    1517                 :            :                 else {
    1518         [ +  - ]:        118 :                     tet_store.mark_derefinement_decision(tet_id, AMR::Derefinement_Case::skip);
    1519                 :            :                     // Deactivate all points"
    1520         [ +  + ]:       1062 :                     for (auto child_id : children) {
    1521         [ +  - ]:        944 :                       deactivate_deref_tet_edges(child_id);
    1522                 :            :                     }
    1523                 :        118 :                     trace_out << "giving up with no deref decision because nderefine = 0" << std::endl;
    1524                 :            :                 }
    1525         [ +  + ]:      14957 :             }
    1526                 :            : 
    1527                 :            :             // If nothing changed during that round, break
    1528         [ +  + ]:         65 :             if (!tet_store.marked_derefinements.get_state_changed())
    1529                 :            :             {
    1530                 :         40 :                 trace_out << "Terminating loop at iter " << iter << std::endl;
    1531                 :         40 :                 break;
    1532                 :            :             }
    1533                 :         25 :             trace_out << "End iter " << iter << std::endl;
    1534                 :            :             // clear out set of elements considered during this iteration
    1535                 :         25 :             done_deref_marking.clear();
    1536         [ +  + ]:         65 :         }
    1537                 :         40 :         trace_out << "Deref Loop took " << iter << " rounds." << std::endl;
    1538                 :         40 :     }
    1539                 :            : 
    1540                 :            :     // TODO: document
    1541                 :         24 :     void mesh_adapter_t::perform_derefinement()
    1542                 :            :     {
    1543                 :         24 :         trace_out << "Perform deref" << std::endl;
    1544                 :            : 
    1545                 :            :         // Do derefinements
    1546         [ +  + ]:       7638 :         for (const auto& kv : tet_store.tets)
    1547                 :            :         {
    1548                 :       7614 :             size_t tet_id = kv.first;
    1549                 :            :             //size_t parent_id = 0;
    1550                 :            : 
    1551                 :            :             // TODO: Do I really want to loop all tets?
    1552                 :            : 
    1553                 :            :             // TODO: is this doing a double lookup?
    1554 [ +  - ][ +  + ]:       7614 :             if (tet_store.has_derefinement_decision(tet_id))
    1555                 :            :             {
    1556                 :       6562 :                 trace_out << "Do derefine of " << tet_id << std::endl;
    1557                 :            :                 //size_t parent_id = tet_store.get_parent_id(tet_id);
    1558                 :            :                 //trace_out << "Parent = " << parent_id << std::endl;
    1559 [ +  - ][ +  + ]:       6562 :                 switch(tet_store.marked_derefinements.get(tet_id))
         [ -  + ][ -  - ]
                 [ +  - ]
    1560                 :            :                 {
    1561                 :        130 :                     case AMR::Derefinement_Case::two_to_one:
    1562         [ +  - ]:        130 :                         refiner.derefine_two_to_one(tet_store,node_connectivity,tet_id);
    1563                 :        130 :                         trace_out << "Completed derefine 2:1 of " << tet_id << std::endl;
    1564                 :        130 :                         break;
    1565                 :        176 :                     case AMR::Derefinement_Case::four_to_one:
    1566         [ +  - ]:        176 :                         refiner.derefine_four_to_one(tet_store,node_connectivity,tet_id);
    1567                 :        176 :                         trace_out << "Completed derefine 4:1 of " << tet_id << std::endl;
    1568                 :        176 :                         break;
    1569                 :          0 :                     case AMR::Derefinement_Case::four_to_two:
    1570         [ -  - ]:          0 :                         refiner.derefine_four_to_two(tet_store,node_connectivity,tet_id);
    1571                 :          0 :                         trace_out << "Completed derefine 4:2 of " << tet_id << std::endl;
    1572                 :          0 :                         break;
    1573                 :       5974 :                     case AMR::Derefinement_Case::eight_to_one:
    1574         [ +  - ]:       5974 :                         refiner.derefine_eight_to_one(tet_store,node_connectivity,tet_id);
    1575                 :       5974 :                         trace_out << "Completed derefine 8:1 of " << tet_id << std::endl;
    1576                 :       5974 :                         break;
    1577                 :          0 :                     case AMR::Derefinement_Case::eight_to_two:
    1578         [ -  - ]:          0 :                         refiner.derefine_eight_to_two(tet_store,node_connectivity,tet_id);
    1579                 :          0 :                         trace_out << "Completed derefine 8:2 of " << tet_id << std::endl;
    1580                 :          0 :                         break;
    1581                 :          0 :                     case AMR::Derefinement_Case::eight_to_four:
    1582         [ -  - ]:          0 :                         refiner.derefine_eight_to_four(tet_store,node_connectivity,tet_id);
    1583                 :          0 :                         trace_out << "Completed derefine 8:4 of " << tet_id << std::endl;
    1584                 :          0 :                         break;
    1585                 :        282 :                     case AMR::Derefinement_Case::skip:
    1586                 :            :                         // What do we do with skip?
    1587                 :        282 :                         break;
    1588                 :            :                 }
    1589                 :            :                 // Mark tet as not needing derefinement
    1590         [ +  - ]:       6562 :                 tet_store.marked_derefinements.erase(tet_id);
    1591                 :            :             }
    1592                 :            :         }
    1593                 :            : 
    1594                 :         24 :         node_connectivity.print();
    1595                 :         24 :         refiner.delete_intermediates_of_children(tet_store);
    1596                 :         24 :         tet_store.process_delete_list();
    1597                 :         24 :         tet_store.print_node_types();
    1598                 :            : 
    1599                 :         24 :         lock_intermediates();
    1600                 :            : 
    1601         [ +  + ]:     130086 :         for (auto& kv : tet_store.edge_store.edges) {
    1602                 :     130062 :            auto& local = kv.second;
    1603                 :     130062 :            local.needs_derefining = 0;
    1604         [ +  + ]:     130062 :            if (local.needs_refining == 2) local.needs_refining = 0;
    1605                 :            :         }
    1606                 :         24 :     }
    1607                 :            : 
    1608                 :            : }
    1609                 :            : 
    1610                 :            : #if defined(__clang__)
    1611                 :            :   #pragma clang diagnostic pop
    1612                 :            : #endif

Generated by: LCOV version 1.16