Commit 4d53528d authored by Conor McCoid's avatar Conor McCoid
Browse files

Tetra: higher dimensions - orphaning algo, practical problems started in notes

parent d05638c7
function [H_new,J_new,T_new,C] = SUB_simplices_Orphan_v1_20211126(H,J,T,A,y,X)
% ORPHAN generates adjacency matrix for intersections
% [H_new,J_new,T_new,C] = Orphan(H,J,T,A,y,X) finds all intersections
% H_new between the n-simplex X and the hyperplanes indexed by T and y.
%
% Inputs:
% H - intersections between X and Y, indexed by J and T;
% J - vertices of X that form the m-face on which the intersections
% of H lie;
% T - hyperplanes of Y where the intersections H lie;
% A - defines which intersections of H are neighbours, ie. the
% adjacency matrix for the convex polyhedron embedded in T;
% y - next hyperplane to intersect the convex hull of H with;
% X - position of the vertices of X.
%
% Outputs:
% H_new - intersections between X and Y indexed by J_new and T_new;
% J_new - vertices of X for the intersections of H_new;
% T_new - T in union with y;
% C - adjacency matrix for H_new.
[n,m]=size(H);
m_max=floor((m+1)/2)*ceil((m+1)/2);
m_new=0; B=false(m,m_max);
T_new=T; T_new(y)=1; J_new=zeros(n,m_max); H_new=zeros(n,m_max);
for i=1:m-1
ind=1:m; ind=ind(A(i,:));
for j=ind
if sign(H(y,i))~=sign(H(y,j))
% calculate intersection
m_new=m_new+1; B([i,j],m_new)=true;
J_new(:,m_new)=J(:,i) | J(:,j);
% H_new(:,m_new)=SUB_simplices_Intersection_v1_?(X,J_new(:,m_new),T_new);
end
end
end
J_new=J_new(:,1:m_new); H_new=H_new(:,1:m_new);
C=false(m_new,m_new);
for i=1:m_new
temp=B(B(:,i),i+1:m_new);
C(i,(i+1):end)=temp(1,:) | temp(2,:);
end
\ No newline at end of file
......@@ -1146,4 +1146,58 @@ Likewise, one can show the calculations in step m, repeated once in this case, a
Since these calculations are the same, the comparison of the signs of the intersections will return the same results.
Thus, both algorithms will agree as to which vertices of $Y$ lie inside $X$.
\subsection{Implementation of the algorithm}
There are a number of practical concerns when it comes to implementing the algorithm.
\subsubsection{Adjacent intersections}
When computing the intersections at a given step $m$ it is important to know which pairs of intersections from the previous step to compare.
Comparing all possible pairs will result in an excessive number of intersections, many of which will need to be discarded immediately as they lie on the inside of $Z$.
Instead, one must determine which intersections at a given step are 'adjacent'.
Two intersections $\vec{h}(J | \Gamma)$ and $\vec{h}(K | \Gamma)$ are said to be adjacent if $J$ and $K$ differ by a single index.
If two intersections are adjacent then both lie on the same $m$--face of $X$ and the edge between them lies on the exterior of $X$.
Knowing only $J$ and $K$ it takes $n+1$ comparisons to determine whether the two intersections are adjacent.
Performing this for all possible pairs could prove prohibitively expensive.
As an alternative, we propose a subroutine to be performed on the adjacency matrix from the previous step.
Take the graph of the set of intersections $\set{\vec{h}(J | \Gamma)}_J$ with edges defined by its adjacency matrix $A$.
Place the nodes of this graph into two groups, one with non-negative values along the $\vec{e}_\gamma$ direction and one with negative values.
Remove the intragroup edges.
The graph is now bipartite.
Place a new node on each edge.
These new nodes are 'children' of the old 'parent' nodes.
Each child has two parents.
Place an edge between all children that share a parent (siblings), then delete the parent nodes.
The adjacency matrix of this graph is the adjacency matrix for the new set of calculated intersections.
To form this adjacency matrix numerically first use the previous adjacency matrix to compute the intersections.
Whenever such an intersection is computed add a column to the adjacency matrix with two edges indicated between this new intersection and the two that were used to calculate it.
This produces the children nodes.
All children must be formed before placing the sibling edges.
Then, the row associated with a given child node is the 'xor' combination of the two rows of its parents.
The final step of orphaning the children is achieved by taking the lower right block of the adjacency matrix.
% add graph of orphaning, tikz layout of adjacency matrix, and pseudocode
\subsubsection{Combinations of hyperplanes}
As we proceed through the algorithm we must continuously add hyperplanes to the intersections.
We must somehow navigate the many combinations of these hyperplanes in an efficient manner.
\subsubsection{Calculation of the intersections}
The actual process of calculating $\vec{h}(J|\Gamma)$ should not be done using the formulas presented here, as these are effectively Cramer's rule and therefore prone to error.
The form presented here is merely to illustrate the connection between the signs of the intersections at various stages.
As this is the most numerically intensive part of the algorithm it is crucial to use a fast and, ideally, accurate method to perform this calculation.
The accuracy of this step will not affect the robustness of the overall algorithm since this relies only on the consistency between signs.
This consistency may (and should) be determined externally to the magnitude of the intersection.
One should consider adapting GMRES to find this solution.
We could also calculate numerators and denominators separately and combine them as needed.
\end{document}
\ No newline at end of file
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment