Skip to content

Commit d88d662

Browse files
authored
Merge pull request #56 from nasa/lb-and-ub
Lb and ub
2 parents b9b5e4b + 305312f commit d88d662

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

59 files changed

+163988
-1452
lines changed

.gitignore

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,8 @@
66
# vs code folders
77
.vscode/
88
*.p3d
9+
# Allow test data p3d files
10+
!python/tests/data/cross_plane_pair.p3d
911
*.x
1012
*.lock
1113
python/mesh.pickle

README.md

Lines changed: 28 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -67,6 +67,34 @@ Never use jupyter notebook unless you are done with your .py file and want to de
6767
| Christopher Keokot | Intern | Fall 2021 | Plot3D ReactJS GUI | |
6868
| Tim Beach | Advisor | Fall 2021 | Block splitting help | |
6969

70+
# Connectivity & Orientation
71+
72+
`connectivity_fast()` finds matching faces between blocks and returns orientation data for each match. When two faces lie on different constant planes (e.g., a K-face connects to a J-face), the varying axes differ and an axis swap may be needed. This is encoded as a `permutation_index` (0–7) using 3 bits:
73+
74+
| Bit | Flag | Meaning |
75+
|-----|------|---------|
76+
| 0 | `u_reversed` | Flip first varying axis |
77+
| 1 | `v_reversed` | Flip second varying axis |
78+
| 2 | `swapped` | Transpose u and v axes |
79+
80+
Each face match includes an `orientation` dict:
81+
82+
```python
83+
face_matches, outer_faces = connectivity_fast(blocks)
84+
# Each match in face_matches has:
85+
# match['orientation'] = {
86+
# 'permutation_index': int, # 0-7
87+
# 'plane': str, # 'in-plane' or 'cross-plane'
88+
# 'permutation_matrix': list # 2x2 signed permutation matrix
89+
# }
90+
```
91+
92+
The 8 permutation matrices are also available as `plot3d.PERMUTATION_MATRICES`.
93+
94+
# Documentation
95+
- [Face Orientation: Cross-Plane Connectivity](docs/notes/unverified_connectivity_findings.md) — why cross-plane face connections need orientation flags beyond lb/ub, and how the 8-permutation system works
96+
- [Presentation (PowerPoint)](docs/notes/unverified_connectivity_findings.pptx) — visual walkthrough of 2D→3D diagonal combinatorics
97+
7098
# Rust Version
7199
The rust version is available here. This version of the code is useful for creating CFD solvers in rust [Plot3D-RS](https://github.com/pjuangph/plot3d-rs) It has most of the functionality of the python version except for some of the GlennHT pre-processing code.
72100

docs/_static/block-split_paraview_plot.py

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -155,30 +155,30 @@ def CreateSubset(block_source,voi:List[int],name:str,opacity:float=1):
155155
# Add Plots for Matched Faces
156156
if f['block1']['block_index'] == b or f['block2']['block_index'] == b :
157157
if f['block1']['block_index'] == b:
158-
voi = [f['block1']['IMIN'], f['block1']['IMAX'], f['block1']['JMIN'], f['block1']['JMAX'],f['block1']['KMIN'], f['block1']['KMAX']]
158+
voi = [f['block1']['lb'][0], f['block1']['ub'][0], f['block1']['lb'][1], f['block1']['ub'][1],f['block1']['lb'][2], f['block1']['ub'][2]]
159159
else:
160-
voi = [f['block2']['IMIN'], f['block2']['IMAX'], f['block2']['JMIN'], f['block2']['JMAX'],f['block2']['KMIN'], f['block2']['KMAX']]
160+
voi = [f['block2']['lb'][0], f['block2']['ub'][0], f['block2']['lb'][1], f['block2']['ub'][1],f['block2']['lb'][2], f['block2']['ub'][2]]
161161
CreateSubset(block_source, voi, name='match '+str(match_indx))
162162

163163
# Plot the outer faces
164164
for surface_indx, o in enumerate(outer_faces):
165165
# Add Plots for Outer Faces
166166
if o['block_index'] == b:
167-
voi = [o['IMIN'], o['IMAX'], o['JMIN'], o['JMAX'],o['KMIN'], o['KMAX']]
167+
voi = [o['lb'][0], o['ub'][0], o['lb'][1], o['ub'][1],o['lb'][2], o['ub'][2]]
168168
CreateSubset(block_source, voi, name='outer_face '+str(surface_indx+1),opacity=0.2)
169169

170170
# Plot the outer faces
171171
for periodic_indx, p in enumerate(periodic_faces):
172172
# Add Plots for Outer Faces
173173
if p['block1']['block_index'] == b and p['block2']['block_index'] == b: # Periodicity within the block
174-
voi = [p['block1']['IMIN'], p['block1']['IMAX'], p['block1']['JMIN'], p['block1']['JMAX'],p['block1']['KMIN'], p['block1']['KMAX']]
174+
voi = [p['block1']['lb'][0], p['block1']['ub'][0], p['block1']['lb'][1], p['block1']['ub'][1],p['block1']['lb'][2], p['block1']['ub'][2]]
175175
CreateSubset(block_source, voi, name='periodic '+str(periodic_indx))
176-
voi = [p['block2']['IMIN'], p['block2']['IMAX'], p['block2']['JMIN'], p['block2']['JMAX'],p['block2']['KMIN'], p['block2']['KMAX']]
176+
voi = [p['block2']['lb'][0], p['block2']['ub'][0], p['block2']['lb'][1], p['block2']['ub'][1],p['block2']['lb'][2], p['block2']['ub'][2]]
177177
CreateSubset(block_source, voi, name='periodic '+str(periodic_indx))
178178

179179
elif p['block1']['block_index'] == b or p['block2']['block_index'] == b: # Periodicity from block to block
180180
if p['block1']['block_index'] == b:
181-
voi = [p['block1']['IMIN'], p['block1']['IMAX'], p['block1']['JMIN'], p['block1']['JMAX'],p['block1']['KMIN'], p['block1']['KMAX']]
181+
voi = [p['block1']['lb'][0], p['block1']['ub'][0], p['block1']['lb'][1], p['block1']['ub'][1],p['block1']['lb'][2], p['block1']['ub'][2]]
182182
else:
183-
voi = [p['block2']['IMIN'], p['block2']['IMAX'], p['block2']['JMIN'], p['block2']['JMAX'],p['block2']['KMIN'], p['block2']['KMAX']]
183+
voi = [p['block2']['lb'][0], p['block2']['ub'][0], p['block2']['lb'][1], p['block2']['ub'][1],p['block2']['lb'][2], p['block2']['ub'][2]]
184184
CreateSubset(block_source, voi, name='periodic '+str(periodic_indx))

docs/_static/paraview_code_part2.py

Lines changed: 8 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -113,30 +113,30 @@ def CreateSubset(block_source,voi:List[int],name:str,opacity:float=1):
113113
# Add Plots for Matched Faces
114114
if f['block1']['block_index'] == b or f['block2']['block_index'] == b :
115115
if f['block1']['block_index'] == b:
116-
voi = [f['block1']['IMIN'], f['block1']['IMAX'], f['block1']['JMIN'], f['block1']['JMAX'],f['block1']['KMIN'], f['block1']['KMAX']]
116+
voi = [f['block1']['lb'][0], f['block1']['ub'][0], f['block1']['lb'][1], f['block1']['ub'][1],f['block1']['lb'][2], f['block1']['ub'][2]]
117117
else:
118-
voi = [f['block2']['IMIN'], f['block2']['IMAX'], f['block2']['JMIN'], f['block2']['JMAX'],f['block2']['KMIN'], f['block2']['KMAX']]
118+
voi = [f['block2']['lb'][0], f['block2']['ub'][0], f['block2']['lb'][1], f['block2']['ub'][1],f['block2']['lb'][2], f['block2']['ub'][2]]
119119
CreateSubset(block_source, voi, name='match '+str(match_indx))
120120

121121
# Plot the outer faces
122122
for o in outer_faces:
123123
# Add Plots for Outer Faces
124124
if o['block_index'] == b:
125-
voi = [o['IMIN'], o['IMAX'], o['JMIN'], o['JMAX'],o['KMIN'], o['KMAX']]
125+
voi = [o['lb'][0], o['ub'][0], o['lb'][1], o['ub'][1],o['lb'][2], o['ub'][2]]
126126
CreateSubset(block_source, voi, name='outer_face '+str(surface_indx),opacity=0.2)
127127
surface_indx +=1
128128

129129
for periodic_indx, p in enumerate(periodic_faces):
130130
# Add Plots for Outer Faces
131131
if p['block1']['block_index'] == b and p['block2']['block_index'] == b: # Periodicity within the block
132-
voi = [p['block1']['IMIN'], p['block1']['IMAX'], p['block1']['JMIN'], p['block1']['JMAX'],p['block1']['KMIN'], p['block1']['KMAX']]
132+
voi = [p['block1']['lb'][0], p['block1']['ub'][0], p['block1']['lb'][1], p['block1']['ub'][1],p['block1']['lb'][2], p['block1']['ub'][2]]
133133
CreateSubset(block_source, voi, name='periodic '+str(periodic_indx))
134-
voi = [p['block2']['IMIN'], p['block2']['IMAX'], p['block2']['JMIN'], p['block2']['JMAX'],p['block2']['KMIN'], p['block2']['KMAX']]
134+
voi = [p['block2']['lb'][0], p['block2']['ub'][0], p['block2']['lb'][1], p['block2']['ub'][1],p['block2']['lb'][2], p['block2']['ub'][2]]
135135
CreateSubset(block_source, voi, name='periodic '+str(periodic_indx))
136136

137-
elif p['block1']['block_index'] == b or p['block2']['block_index'] == b: # Periodicity from block to block
137+
elif p['block1']['block_index'] == b or p['block2']['block_index'] == b: # Periodicity from block to block
138138
if p['block1']['block_index'] == b:
139-
voi = [p['block1']['IMIN'], p['block1']['IMAX'], p['block1']['JMIN'], p['block1']['JMAX'],p['block1']['KMIN'], p['block1']['KMAX']]
139+
voi = [p['block1']['lb'][0], p['block1']['ub'][0], p['block1']['lb'][1], p['block1']['ub'][1],p['block1']['lb'][2], p['block1']['ub'][2]]
140140
else:
141-
voi = [p['block2']['IMIN'], p['block2']['IMAX'], p['block2']['JMIN'], p['block2']['JMAX'],p['block2']['KMIN'], p['block2']['KMAX']]
141+
voi = [p['block2']['lb'][0], p['block2']['ub'][0], p['block2']['lb'][1], p['block2']['ub'][1],p['block2']['lb'][2], p['block2']['ub'][2]]
142142
CreateSubset(block_source, voi, name='periodic '+str(periodic_indx))

docs/notes/connectivity.rst

Lines changed: 161 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -40,6 +40,167 @@ This is the matching face within the omesh
4040
:figclass: align-center
4141

4242

43+
Face Orientation: Cross-Plane Connectivity
44+
*********************************************
45+
46+
When two connected faces share the **same constant axis** (e.g., both K-constant), the two varying axes are identical and only 4 direction combinations exist. These are fully encodable by the ``lb/ub`` diagonal convention.
47+
48+
When faces lie on **different constant planes** (e.g., a K-face connects to a J-face), the varying axes differ. Face A varies in (i, j) while Face B varies in (i, k). This introduces an additional degree of freedom — which axis maps to which:
49+
50+
- **No swap**: A(i,j) maps to B(i,k) — 4 direction combos (encodable by lb/ub)
51+
- **Swapped**: A(i,j) maps to B(k,i) — 4 direction combos (need orientation)
52+
53+
**4 × 2 = 8 total permutations. lb/ub only encodes 4.**
54+
55+
Permutation Matrix System
56+
~~~~~~~~~~~~~~~~~~~~~~~~~~
57+
58+
Both the Python and Rust implementations encode all 8 orientations as a 3-bit ``permutation_index`` into an array of 2×2 signed matrices (``PERMUTATION_MATRICES``):
59+
60+
.. code-block:: text
61+
62+
permutation_index = u_reversed | (v_reversed << 1) | (swapped << 2)
63+
64+
.. list-table:: The 8 Permutation Matrices
65+
:header-rows: 1
66+
:widths: 8 8 8 8 8 20 20
67+
68+
* - Index
69+
- Binary
70+
- u_rev
71+
- v_rev
72+
- swap
73+
- Matrix
74+
- Effect
75+
* - 0
76+
- ``000``
77+
- no
78+
- no
79+
- no
80+
- ``[[ 1, 0],[ 0, 1]]``
81+
- identity
82+
* - 1
83+
- ``001``
84+
- yes
85+
- no
86+
- no
87+
- ``[[-1, 0],[ 0, 1]]``
88+
- flip u
89+
* - 2
90+
- ``010``
91+
- no
92+
- yes
93+
- no
94+
- ``[[ 1, 0],[ 0,-1]]``
95+
- flip v
96+
* - 3
97+
- ``011``
98+
- yes
99+
- yes
100+
- no
101+
- ``[[-1, 0],[ 0,-1]]``
102+
- flip both
103+
* - 4
104+
- ``100``
105+
- no
106+
- no
107+
- yes
108+
- ``[[ 0, 1],[ 1, 0]]``
109+
- transpose
110+
* - 5
111+
- ``101``
112+
- yes
113+
- no
114+
- yes
115+
- ``[[ 0,-1],[ 1, 0]]``
116+
- transpose + flip u
117+
* - 6
118+
- ``110``
119+
- no
120+
- yes
121+
- yes
122+
- ``[[ 0, 1],[-1, 0]]``
123+
- transpose + flip v
124+
* - 7
125+
- ``111``
126+
- yes
127+
- yes
128+
- yes
129+
- ``[[ 0,-1],[-1, 0]]``
130+
- transpose + both
131+
132+
In **Python**, the ``PERMUTATION_MATRICES`` constant is defined in ``connectivity.py`` and ``_orient_vec_to_permutation()`` converts the legacy orientation vector to a ``permutation_index``.
133+
134+
In **Rust**, the ``PERMUTATION_MATRICES`` constant is in ``face_record.rs`` and each ``FaceMatch`` carries an ``Orientation { permutation_index, plane }`` set by ``verify_connectivity``.
135+
136+
The ``u`` and ``v`` names are abstract — they map to concrete i/j/k axes depending on which axis is constant:
137+
138+
.. list-table:: u/v to axis mapping
139+
:header-rows: 1
140+
:widths: 30 35 35
141+
142+
* - Constant axis
143+
- u (outer loop)
144+
- v (inner loop)
145+
* - I-constant
146+
- j
147+
- k
148+
* - J-constant
149+
- i
150+
- k
151+
* - K-constant
152+
- i
153+
- j
154+
155+
So ``u_reversed: true`` means the outer-loop axis runs opposite direction on block2 vs block1, and ``v_reversed: true`` means the inner-loop axis runs opposite.
156+
157+
Legacy Orientation Vector
158+
~~~~~~~~~~~~~~~~~~~~~~~~~~
159+
160+
The older Python ``_compute_orientation()`` in ``connectivity.py`` produces an orientation vector that maps each face1 axis to a face2 axis:
161+
162+
.. code-block:: python
163+
164+
# orientation = [d1_maps_to, d2_maps_to, d3_maps_to]
165+
# 1-indexed: 1=I, 2=J, 3=K
166+
#
167+
# Example: orientation = [2, 1, 3]
168+
# face1 I-axis -> face2 J-axis
169+
# face1 J-axis -> face2 I-axis
170+
# face1 K-axis -> face2 K-axis
171+
172+
Direction (forward/reverse) is encoded in the ``lb/ub`` values, not in the orientation vector. The ``_orient_vec_to_permutation()`` function converts this vector to a ``permutation_index`` for the unified system. The verification modules (``verify_connectivity`` / ``verify_periodicity``) test all 8 permutations to confirm the match.
173+
174+
For a detailed analysis with diagrams, see the `Root Cause Analysis <https://github.com/nasa/Plot3D_utilities/blob/main/docs/notes/unverified_connectivity_findings.md>`_ document.
175+
176+
177+
Directed Diagonal for GHT_CONN Export
178+
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
179+
180+
The GlennHT connectivity file (``.ght_conn``) uses a **directed diagonal** convention where each face is specified by two corners ``(IMIN, JMIN, KMIN)`` and ``(IMAX, JMAX, KMAX)``. For cross-plane matches, the traversal direction is encoded by allowing the "max" corner to be less than the "min" corner on reversed axes. For example:
181+
182+
.. code-block:: text
183+
184+
1 1 1 1 25 409 1
185+
2 409 1 1 1 1 25
186+
187+
Here block 2's i-axis runs 409 → 1 (reversed), encoding the cross-plane orientation without a separate permutation matrix.
188+
189+
When connectivity is computed in **Python** (via ``connectivity_fast``), the ``lb``/``ub`` values already encode the directed diagonal from the point-match traversal order.
190+
191+
When connectivity is computed in **Rust** (via the ``connectivity-finder`` binary in grid-packed), the JSON output uses ascending ``lo``/``hi`` bounds with a ``permutation_index`` (0-7). Use ``reconstruct_directed_diagonal()`` to convert these to directed ``lb``/``ub`` before exporting:
192+
193+
.. code-block:: python
194+
195+
from plot3d.connectivity import reconstruct_directed_diagonal
196+
197+
# face_match has ascending lb/ub + permutation_index from Rust JSON
198+
directed_match = reconstruct_directed_diagonal(face_match)
199+
# directed_match now has directed lb/ub and permutation_index = -1
200+
201+
The ``export_to_glennht_conn()`` function applies this reconstruction automatically, so callers do not need to call it explicitly.
202+
203+
43204
Plotting Connectivity using Paraview
44205
****************************************
45206
This example shows how you can take a mesh created Numeca Autogrid and plot the connecitivity using paraview. With this information, anyone should be able to comprehend the format and export it to whatever solver.

0 commit comments

Comments
 (0)