54 unsigned targetNodeIndex,
60 assert(DIM == 2 || DIM == 3);
65 bool current_node_contained = !containing_elements.empty();
66 bool target_node_contained = !new_location_containing_elements.empty();
69 assert(new_location_containing_elements.size() < 2);
71 if (!current_node_contained && !target_node_contained)
73 EXCEPTION(
"At least one of the current node or target node must be in an element.");
76 if (current_node_contained && target_node_contained)
78 if (*(new_location_containing_elements.begin()) == *(containing_elements.begin()))
80 EXCEPTION(
"The current node and target node must not be in the same element.");
85 unsigned neighbours_in_same_element_as_current_node = 0;
86 unsigned neighbours_in_same_element_as_target_node = 0;
87 std::set<unsigned> target_neighbouring_node_indices = rCellPopulation.
rGetMesh().GetVonNeumannNeighbouringNodeIndices(targetNodeIndex);
88 for (std::set<unsigned>::iterator iter = target_neighbouring_node_indices.begin();
89 iter != target_neighbouring_node_indices.end();
92 std::set<unsigned> neighbouring_node_containing_elements = rCellPopulation.
rGetMesh().GetNode(*iter)->rGetContainingElementIndices();
95 assert(neighbouring_node_containing_elements.size() < 2);
97 bool neighbouring_node_contained = !neighbouring_node_containing_elements.empty();
99 if (neighbouring_node_contained && target_node_contained)
101 unsigned neighbour_element = (*neighbouring_node_containing_elements.begin());
102 unsigned target_element = (*new_location_containing_elements.begin());
103 if (target_element == neighbour_element)
105 neighbours_in_same_element_as_target_node++;
108 if (neighbouring_node_contained && current_node_contained)
110 unsigned neighbour_element = (*neighbouring_node_containing_elements.begin());
111 unsigned current_element = (*containing_elements.begin());
112 if (current_element == neighbour_element)
114 neighbours_in_same_element_as_current_node++;
121 assert(neighbours_in_same_element_as_current_node<5);
122 assert(neighbours_in_same_element_as_target_node<5);
124 double change_in_surface_area[5] = {4,2,0,-2,-4};
126 if (current_node_contained)
128 unsigned current_element = (*containing_elements.begin());
129 double current_surface_area = rCellPopulation.
rGetMesh().GetSurfaceAreaOfElement(current_element);
130 double current_surface_area_difference = current_surface_area - mMatureCellTargetSurfaceArea;
131 double current_surface_area_difference_after_switch = current_surface_area_difference + change_in_surface_area[neighbours_in_same_element_as_current_node];
133 delta_H += mDeformationEnergyParameter*(current_surface_area_difference_after_switch*current_surface_area_difference_after_switch - current_surface_area_difference*current_surface_area_difference);
135 if (target_node_contained)
137 unsigned target_element = (*new_location_containing_elements.begin());
138 double target_surface_area = rCellPopulation.
rGetMesh().GetSurfaceAreaOfElement(target_element);
139 double target_surface_area_difference = target_surface_area - mMatureCellTargetSurfaceArea;
140 double target_surface_area_difference_after_switch = target_surface_area_difference - change_in_surface_area[neighbours_in_same_element_as_target_node];
142 delta_H += mDeformationEnergyParameter*(target_surface_area_difference_after_switch*target_surface_area_difference_after_switch - target_surface_area_difference*target_surface_area_difference);
147 assert(neighbours_in_same_element_as_current_node < 7);
148 assert(neighbours_in_same_element_as_target_node < 7);
150 double change_in_surface_area[7] = {6,4,2,0,-2,-4,-6};
152 if (current_node_contained)
154 unsigned current_element = (*containing_elements.begin());
155 double current_surface_area = rCellPopulation.
rGetMesh().GetSurfaceAreaOfElement(current_element);
156 double current_surface_area_difference = current_surface_area - mMatureCellTargetSurfaceArea;
157 double current_surface_area_difference_after_switch = current_surface_area_difference + change_in_surface_area[neighbours_in_same_element_as_current_node];
159 delta_H += mDeformationEnergyParameter*(current_surface_area_difference_after_switch*current_surface_area_difference_after_switch - current_surface_area_difference*current_surface_area_difference);
161 if (target_node_contained)
163 unsigned target_element = (*new_location_containing_elements.begin());
164 double target_surface_area = rCellPopulation.
rGetMesh().GetSurfaceAreaOfElement(target_element);
165 double target_surface_area_difference = target_surface_area - mMatureCellTargetSurfaceArea;
166 double target_surface_area_difference_after_switch = target_surface_area_difference - change_in_surface_area[neighbours_in_same_element_as_target_node];
168 delta_H += mDeformationEnergyParameter*(target_surface_area_difference_after_switch*target_surface_area_difference_after_switch - target_surface_area_difference*target_surface_area_difference);