62 double heterotypic_boundary_length = 0.0;
63 double total_shared_edges_length = 0.0;
64 double num_heterotypic_pairs = 0.0;
65 double total_num_pairs = 0.0;
73 cell_iter != pCellPopulation->
End();
80 bool cell_is_labelled = cell_iter->template HasCellProperty<CellLabel>();
86 for (std::set<unsigned>::iterator neighbour_iter = neighbour_indices.begin();
87 neighbour_iter != neighbour_indices.end();
91 unsigned neighbour_index = *neighbour_iter;
98 total_shared_edges_length += edge_length;
99 total_num_pairs += 1.0;
103 bool neighbour_is_labelled = p_neighbour_cell->template HasCellProperty<CellLabel>();
106 if (cell_is_labelled != neighbour_is_labelled)
109 heterotypic_boundary_length += edge_length;
110 num_heterotypic_pairs += 1.0;
117 heterotypic_boundary_length *= 0.5;
118 total_shared_edges_length *= 0.5;
121 num_heterotypic_pairs *= 0.5;
122 total_num_pairs *= 0.5;
124 *this->mpOutStream << heterotypic_boundary_length <<
"\t" << total_shared_edges_length <<
"\t" << num_heterotypic_pairs <<
"\t" << total_num_pairs;
131 double heterotypic_boundary_length = 0.0;
132 double total_shared_edges_length = 0.0;
133 double num_heterotypic_pairs = 0.0;
134 double total_num_pairs = 0.0;
138 cell_iter != pCellPopulation->
End();
145 bool cell_is_labelled = cell_iter->template HasCellProperty<CellLabel>();
148 std::set<unsigned> neighbour_node_indices = pCellPopulation->
rGetMesh().GetVonNeumannNeighbouringNodeIndices(index);
151 for (std::set<unsigned>::iterator neighbour_iter = neighbour_node_indices.begin();
152 neighbour_iter != neighbour_node_indices.end();
156 double edge_length = 1.0;
158 unsigned neighbour_index = *neighbour_iter;
164 total_shared_edges_length += edge_length;
165 total_num_pairs += 1.0;
168 bool neighbour_is_labelled = p_neighbour_cell->template HasCellProperty<CellLabel>();
171 if (cell_is_labelled != neighbour_is_labelled)
174 heterotypic_boundary_length += edge_length;
175 num_heterotypic_pairs += 1.0;
182 heterotypic_boundary_length *= 0.5;
183 total_shared_edges_length *= 0.5;
186 num_heterotypic_pairs *= 0.5;
187 total_num_pairs *= 0.5;
189 *this->mpOutStream << heterotypic_boundary_length <<
"\t" << total_shared_edges_length <<
"\t" << num_heterotypic_pairs <<
"\t" << total_num_pairs;
197 pCellPopulation->
Update();
200 double heterotypic_boundary_length = 0.0;
201 double total_shared_edges_length = 0.0;
202 double num_heterotypic_pairs = 0.0;
203 double total_num_pairs = 0.0;
207 cell_iter != pCellPopulation->
End();
211 bool cell_is_labelled = cell_iter->template HasCellProperty<CellLabel>();
220 if (!neighbour_indices.empty())
223 for (std::set<unsigned>::iterator neighbour_iter = neighbour_indices.begin();
224 neighbour_iter != neighbour_indices.end();
228 double neighbour_radius = pCellPopulation->
GetNode(*neighbour_iter)->
GetRadius();
232 double sum_of_radii = node_radius + neighbour_radius;
235 if (separation < sum_of_radii)
238 double a = node_radius;
239 double b = neighbour_radius;
240 double c = separation;
241 double s = 0.5*(a + b + c);
242 double A = sqrt(s*(s-a)*(s-b)*(s-c));
243 double edge_length = 4.0*A/c;
245 total_shared_edges_length += edge_length;
246 total_num_pairs += 1.0;
250 bool neighbour_is_labelled = p_neighbour_cell->template HasCellProperty<CellLabel>();
253 if (cell_is_labelled != neighbour_is_labelled)
256 heterotypic_boundary_length += edge_length;
257 num_heterotypic_pairs += 1.0;
265 heterotypic_boundary_length *= 0.5;
266 total_shared_edges_length *= 0.5;
269 num_heterotypic_pairs *= 0.5;
270 total_num_pairs *= 0.5;
272 *this->mpOutStream << heterotypic_boundary_length <<
"\t" << total_shared_edges_length <<
"\t" << num_heterotypic_pairs <<
"\t" << total_num_pairs;
281 double heterotypic_boundary_length = 0.0;
282 double total_shared_edges_length = 0.0;
283 double num_heterotypic_pairs = 0.0;
284 double total_num_pairs = 0.0;
288 cell_iter != pCellPopulation->
End();
292 bool cell_is_labelled = cell_iter->template HasCellProperty<CellLabel>();
295 unsigned num_nodes_in_elem = pCellPopulation->
rGetMesh().GetElement(elem_index)->GetNumNodes();
298 for (
unsigned local_index=0; local_index<num_nodes_in_elem; local_index++)
301 unsigned global_index = pCellPopulation->
rGetMesh().GetElement(elem_index)->GetNodeGlobalIndex(local_index);
302 std::set<unsigned> neighbour_node_indices = pCellPopulation->
rGetMesh().GetVonNeumannNeighbouringNodeIndices(global_index);
305 for (std::set<unsigned>::iterator neighbour_iter = neighbour_node_indices.begin();
306 neighbour_iter != neighbour_node_indices.end();
312 if (neighbour_elem_indices.size() == 1)
314 unsigned neigbour_elem_index = *(neighbour_elem_indices.begin());
316 if (neigbour_elem_index != elem_index)
319 total_shared_edges_length += 1.0;
323 bool neighbour_is_labelled = p_neighbour_cell->template HasCellProperty<CellLabel>();
326 if (cell_is_labelled != neighbour_is_labelled)
329 heterotypic_boundary_length += 1.0;
343 for (std::set<unsigned>::iterator neighbour_iter = neighbour_node_indices.begin();
344 neighbour_iter != neighbour_node_indices.end();
347 total_num_pairs += 1.0;
350 bool neighbour_is_labelled = p_neighbour_cell->template HasCellProperty<CellLabel>();
353 if (cell_is_labelled != neighbour_is_labelled)
356 num_heterotypic_pairs += 1.0;
362 heterotypic_boundary_length *= 0.5;
363 total_shared_edges_length *= 0.5;
366 num_heterotypic_pairs *= 0.5;
367 total_num_pairs *= 0.5;
369 *this->mpOutStream << heterotypic_boundary_length <<
"\t" << total_shared_edges_length <<
"\t" << num_heterotypic_pairs <<
"\t" << total_num_pairs;
377 pCellPopulation->
Update();
380 double heterotypic_boundary_length = 0.0;
381 double total_shared_edges_length = 0.0;
382 double num_heterotypic_pairs = 0.0;
383 double total_num_pairs = 0.0;
387 cell_iter != pCellPopulation->
End();
391 bool cell_is_labelled = cell_iter->template HasCellProperty<CellLabel>();
398 for (std::set<unsigned>::iterator neighbour_iter = neighbour_elem_indices.begin();
399 neighbour_iter != neighbour_elem_indices.end();
403 unsigned neighbour_index = *neighbour_iter;
406 total_shared_edges_length += edge_length;
407 total_num_pairs += 1.0;
411 bool neighbour_is_labelled = p_neighbour_cell->template HasCellProperty<CellLabel>();
414 if (cell_is_labelled != neighbour_is_labelled)
417 heterotypic_boundary_length += edge_length;
418 num_heterotypic_pairs += 1.0;
424 heterotypic_boundary_length *= 0.5;
425 total_shared_edges_length *= 0.5;
428 num_heterotypic_pairs *= 0.5;
429 total_num_pairs *= 0.5;
431 *this->mpOutStream << heterotypic_boundary_length <<
"\t" << total_shared_edges_length <<
"\t" << num_heterotypic_pairs <<
"\t" << total_num_pairs;
437 auto& r_mesh = pCellPopulation->
rGetMesh();
440 double heterotypic_boundary_length = 0.0;
441 double total_shared_edges_length = 0.0;
444 std::set<std::pair<unsigned, unsigned>> total_elem_pairs;
445 std::set<std::pair<unsigned, unsigned>> heterotypic_elem_pairs;
448 const bool update_diagram =
true;
452 const std::vector<unsigned>& node_idx_to_voronoi_cell_id_map = r_mesh.GetVoronoiCellIdsIndexedByNodeIndex();
454 for (
const auto& p_node : r_mesh.rGetNodes())
456 const unsigned this_node_idx = p_node->GetIndex();
457 const unsigned this_elem_idx = *p_node->ContainingElementsBegin();
460 const unsigned voronoi_cell_id = node_idx_to_voronoi_cell_id_map[this_node_idx];
461 const auto& voronoi_cell = vd.cells()[voronoi_cell_id];
464 auto p_edge = voronoi_cell.incident_edge();
467 if (p_edge->is_finite())
470 const unsigned twin_node_idx = p_edge->twin()->cell()->color();
471 const unsigned twin_elem_idx = *r_mesh.GetNode(twin_node_idx)->ContainingElementsBegin();
474 if (this_elem_idx != twin_elem_idx)
477 const unsigned lo_idx = std::min(this_elem_idx, twin_elem_idx);
478 const unsigned hi_idx = std::max(this_elem_idx, twin_elem_idx);
479 total_elem_pairs.insert(std::make_pair(lo_idx, hi_idx));
481 const bool this_is_labelled = pCellPopulation->
GetCellUsingLocationIndex(this_elem_idx)->template HasCellProperty<CellLabel>();
482 const bool twin_is_labelled = pCellPopulation->
GetCellUsingLocationIndex(twin_elem_idx)->template HasCellProperty<CellLabel>();
484 const double edge_length = r_mesh.CalculateLengthOfVoronoiEdge(*p_edge);
486 total_shared_edges_length += edge_length;
488 if (this_is_labelled != twin_is_labelled)
490 heterotypic_boundary_length += edge_length;
491 heterotypic_elem_pairs.insert(std::make_pair(lo_idx, hi_idx));
496 p_edge = p_edge->next();
497 }
while (p_edge != voronoi_cell.incident_edge());
501 heterotypic_boundary_length *= 0.5;
502 total_shared_edges_length *= 0.5;
504 *this->mpOutStream << heterotypic_boundary_length <<
'\t'
505 << total_shared_edges_length <<
'\t'
506 << heterotypic_elem_pairs.size() <<
'\t'
507 << total_elem_pairs.size();