Debug: Database connection successful OpenFOAM (Page 4) / Science, Technology, and Astronomy / New Mars Forums

New Mars Forums

Official discussion forum of The Mars Society and MarsNews.com

You are not logged in.

Announcement

Announcement: This forum is accepting new registrations via email. Please see Recruiting Topic for additional information. Write newmarsmember[at_symbol]gmail.com.

#76 2025-06-12 19:37:44

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

In today's report, I note that the program ran without generating Python errors, but no records were written to the cell memory file, so we have work to do there.  The main change in today's session was to stop trying to pass every parameter to the numerous functions.  Instead, we set up global memory structures for old OpenFOAM data and for the new data we are creating.  This fundamental change led to numerous rounds of testing and fixing, and at this point I haven't looked at the data to see if it is still coming out right. 
(th)

Online

Like button can go here

#77 2025-06-18 18:51:36

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

Today's work session was short but I am hoping it was successful. Having already created code to pull all data cleanly into memory, we are now writing each of the seven files to output so that they can be confirmed as identical to the input. We have two files already, and in this evening's work session we added three more. I'll test the new code tomorrow.

Assuming the three new files are correct, my plan for tomorrow's work session is to add the two memory files to the clean write process. These files are needed to pass data from one layer to the next.

We have already written code (multiple times) to generate new faces and to write needed data to output. I moderately confident that the procedure to create new faces can be added to this clean structure without breaking anything.

(th)

Online

Like button can go here

#78 2025-06-19 12:42:55

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

Today's work session went well... There are seven files to open, and seven files to write out after processing.

We completed work on all seven today, so at long last, we can return to processing the data to create new faces tomorrow.

I'll be watching like a hawk to see to it that nothing that will break the new structure is allowed to enter.  We now are reading all the data that is needed for processing into memory, and we are running a variety of checks to be sure the data we've read in is reasonable and that counts match. We also write the data from memory back out to file, and the difference program confirms the outputs match the inputs.  This moment has been a long time coming.  In retrospect, I can see that I should have thought of this long ago, but the learning curve for OpenFOAM is steep, not to mention the learning curve for Python. In any case, this report might be useful to someone approaching a similar project for the first time.

(th)

Online

Like button can go here

#79 2025-06-23 20:15:24

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

Several days have gone by.... today's report is about the all important "faces" function.  The new layer program had been brought to the point where we were ready for update processing to begin. The most significant of the updates is creation of new faces.  Generating new faces requires passing information from the previous layer, and then doing lookups at run time. I ** think ** the new code might be working, but all I know at this point is that the run statistics look reasonable. I'm planning to look at the output tomorrow:

python3 layer_addition_program_V75.py
✅ points: 28380 entries loaded correctly.
✅ faces: 55814 entries loaded correctly.
✅ owner: 55814 entries loaded correctly.
✅ neighbour: 27436 entries loaded correctly.
✅ boundary: 5 patch blocks loaded correctly.
✅ All OpenFOAM files loaded into base_data.
✅ All face labels valid and complete.
✅ Added 52 new points at X=-0.4
✅ Wrote points with 28432 entries.
✅ Points updated and written.
? Generating new persistent memory cells and faces...
✅ Generated 125 new faces and updated cell memory.
✅ Wrote faces with 55939 entries.
✅ Faces copied and written
✅ Wrote owner with 55814 entries.
✅ Owner copied and written.
✅ Wrote neighbour with 27436 entries.
✅ Neighbour copied and written.
✅ Wrote boundary with 5 patches.
✅ Boundary copied and written.
✅ Persistent memory cells copied and written.
✅ Master memory copied and written.

(th)

Online

Like button can go here

#80 2025-07-02 12:45:53

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

Today's update is about an return to dealing with points.  In working with faces, it turned out the program was computing new points incorrectly, so we added a memory file for points, and incorporated that into the flow.  The new points are now being generated correctly, and (equally important) the memory file for points is being updated to match, so the next layer will be built upon the layer just completed.

In the next work session, we'll return to dealing with faces.  Faces must be built of a combination of new and old points.  Points that were at the top of the existing layer of cells must be rotated to the bottom, and the new points added at the top. While this is done, the order/sequence of points must be done correctly. Correctly in this context means in clockwise order for boundary (out facing) faces, and counter-clockwise for Internal faces (ones shared between cells).

Here is what the output of the run looks like:

python3 layer_addition_program_V78.py
✅ points: 28380 entries loaded correctly.
✅ faces: 55814 entries loaded correctly.
✅ owner: 55814 entries loaded correctly.
✅ neighbour: 27436 entries loaded correctly.
✅ boundary: 5 patch blocks loaded correctly.
✅ All OpenFOAM files loaded into base_data.
✅ All face labels valid and complete.
✅ Loaded 52 points from persistent memory.
✅ Added 52 new points at X=-0.4
✅ Wrote points with 28432 entries.
✅ Points updated and written.
? Generating new persistent memory cells and faces...
✅ Generated 125 new faces and updated cell memory.
✅ Wrote faces with 55939 entries.
✅ Faces copied and written
✅ Wrote owner with 55814 entries.
✅ Owner copied and written.
✅ Wrote neighbour with 27436 entries.
✅ Neighbour copied and written.
✅ Wrote boundary with 5 patches.
✅ Boundary copied and written.
✅ Persistent memory cells copied and written.
✅ Master memory copied and written.
✅ Point memory copied and written.

(th)

Online

Like button can go here

#81 2025-07-03 19:03:25

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

Today's work session was spent in preparing for the heart of the layer program: generation of new faces.

Generation of new points appears to be working correctly, and the rest of the program appears to be loading and writing data correctly.  Until now, we've put dummy data into the faces arrays because it is a major step to handle them correctly.

ChatGPT4o surprised me by wanting to use rotation with angles from the center of the new face, to determine the needed order of storage. The rule I've learned from the documentation is that we we are dealing with a boundary face, rotation needs to be clockwise, and if an internal face (between cells) then rotation is counter-clockwise.  ChatGPT4o decided to compute the centroid of each face, and then (somehow) figure out which point goes where in the output.  The output looks like this: (1234,5678,9012,3456) if the face is boundary. The points are located at the corners of the face, and you start reading clockwise from 1234.  It will be interesting to see if ChatGPT4o's solution to this problem works.

The changes were installed in the layer program and I'll test it tomorrow.  I'd be astonished if the new code works the first try, bue we'll see.

(th)

Online

Like button can go here

#82 2025-07-05 19:00:23

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

Today's work session revealed how difficult this final phase of program development is going to be.

ChatGPT4o made a couple of tiny changes to try to improve the process for:
1) Deciding which points to use for a new face
2) Deciding what order the points must be defined to meet OpenFOAM requirements
3) Including the right number of points (four in our case)

As we tried the new code, things started to get better and then they got worse.  I decided we need a break. 

To lighten things up, I requested a utility program to look at faces.  In order to know if the new program code is working, I need to be able to tell where the points are in the 3D space we are occupying.  A face definition looks like this: 4(1234 5678 9012 3456)
The four numbers represent points in space.  1234 has X, Y and Z coordinates, which I cannot know by just looking at the ID of 1234.
In addition, face 1234 is oriented in one of six possible ways.  It will be perpendicular to -X +X -Y +Y -Z or +Z axes.

The new utility program ** may ** help with visualization of the faces.

We are standing at Version 82 of the layer program.  The only file that is 100% finished is points, and that leaves:
faces
owner
neighbour
master_memory.json
persistent_memory_cells.json
and
persistent_memory_points.json

When we get past faces, the rest should be an easy downhill slide.
(th)

Online

Like button can go here

#83 2025-07-13 20:47:40

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

We are up to Version 91 of the Layer program...

Version 91 is focused intently on the final and most difficulty phase of development.

The program must create new faces, using the old faces as a model, while pulling new points from the newly created points collection, and assembling all four points for each face in the correct sequence.  In the course of debugging, we discovered a small number of errors in configuration of the cell memory file. I have begun a manual review of the cell memory file, using the physical model of the cells developed earlier.  It takes about 20 minutes to process each cell, so about 7 hours of work is pending.  Automation has helped to make the problem manageable, but automation can only do so much.  The issue at hand is rotation to be assigned to the points in the new cells.  Some faces are easy to identify and thus easy to set in the correct rotation. The challenge comes in when we hit a few faces that are changing context.  The original mesh was created with the layer we are working with set as Ground Zero for the entire mesh.  The software that created the mesh built the entire mesh from that layer, while looking at the layer from the +X direction. We are looking at that same layer from the -X direction, because we are attempting to extend the intake away from the existing mesh, in order to build a heating tube for the Optical Plane vessel.

A very small number of faces are in the wrong rotation, given the new perspective.  These have to be adjusted manually.

Yes, anything can be automated, but the time required to persuade the software to perform the needed adjustments is (is is judged to be) greater than the time required to simply perform the manual review and edits.

(th)

Online

Like button can go here

#84 2025-07-16 19:36:11

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

We are up to version 102 of the layer program.  The difficulty we were having with finding points with not-quite-exact floating point numbers appears to have been solved with a python math function math.fabs() or another very similar one ....

(0.184, 0.06303653966)
(0.184, -0.06303653966)
(0.2, 0.06303653966)
(0.2, -0.06303653966)
✅ Generated 50 new faces and updated cell memory.
✅ Wrote faces with 55864 entries.
✅ Faces copied and written
✅ Wrote owner with 55814 entries.
✅ Owner copied and written.
✅ Wrote neighbour with 27436 entries.
✅ Neighbour copied and written.
✅ Wrote boundary with 5 patches.
✅ Boundary copied and written.
✅ Persistent memory cells copied and written.
✅ Master memory copied and written.
✅ Point memory copied and written.

The quote above shows a bit of the debugging data (coordinates of Y and Z) and the finish of the program claiming good news. I'll believe ** that ** when I see it.

(th)

Online

Like button can go here

#85 2025-07-18 20:03:14

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

We are up to Version 105, and it looks as though we might have a side face ...

The utility program whose output is shown below collects information from the persistent_memory_cells.json file and from the points file, and displays the data for each face within each cell. 

cell_13875 (ID: 13875)
----------------------------------------
Face: x_plus (ID: 55814, (label: x_plus, type: internal)
  14117: (-0.4000000000 0.1840000000 0.0630365397)
  28307: (-0.4000000000 0.1840000000 -0.0630365397)
  27753: (-0.4000000000 0.2000000000 -0.0630365397)
  13563: (-0.4000000000 0.2000000000 0.0630365397)

Face: y_plus (ID: 55815, (label: y_plus, type: boundary)
  28431: (-0.6500000000 0.2000000000 -0.0630365397)
  28430: (-0.6500000000 0.2000000000 0.0630365397)
  13563: (-0.4000000000 0.2000000000 0.0630365397)
  27753: (-0.4000000000 0.2000000000 -0.0630365397)

Face: x_minus (ID: 55816, (label: x_minus, type: boundary)
  28428: (-0.6500000000 0.1840000000 0.0630365397)
  28430: (-0.6500000000 0.2000000000 0.0630365397)
  28431: (-0.6500000000 0.2000000000 -0.0630365397)
  28429: (-0.6500000000 0.1840000000 -0.0630365397)

The points in x_plus are all on the -0.4 line, which is the top of the old layer and the floor of the new layer
The points in x_minus are all on the -0.65 line, which is the top of the new layer
The points in y_plus are the ones that are of interest. Two points are in the old layer and two are in the new.
What I won't know until I check against the physical model is whether the rotation is correct, or even if the points are correct.

This output has been a long time coming. I'm glad to see it.  There is still some work to be done to finish processing of an entire row, but this was the most difficult hurdle we've encountered so far.

(th)

Online

Like button can go here

#86 2025-08-14 09:14:10

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

This post is to report on completion of another phase in development of a program to create a new layer of OpenFOAM cells to add to an existing mesh.

In the past month, from the July 18th report, I've been working on creating a set of labels to apply to 3inch cubical boxes to I can look at a physical model of the layer we are trying to match.  The labels are complete, and my spot check of data seems to indicate the model is accurate.  However, i wanted even more reassurance, so today I asked ChatGPT4o to create another "Inspector" program to look at the persistent_memory_cells.json file, to make sure the labels marked"internal" are actually shared with other cells.  Here is output from the new utility.

python3 inspect_face_typesV02.py
✅ Face Type Inspection Complete
  Total faces checked: 150
  Matches: 150
  Mismatches: 0
  Missing from owner: 0
  Face type tally from memory:
    Boundary: 77
    Internal: 73

The all-important numbers are there at the end ... we have 77 out facing faces, and 73 that are shared with other cells.

The next phase of this project is to resume work on the Layer program.  When we left off, we had just completed development of the procedure to create new faces.  I took the past month to prepare for study of the output.  I wanted the previous layer information to be rock solid, and I think that is now the case.

(th)

Online

Like button can go here

#87 2025-08-15 21:12:40

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

This post is to report on my first session with ChatGPT5.  As readers of this topic may recall, I've been working with ChatGPT4o over recent months. The change to Version 5 coincided with the decision to open a new project. I've learned to recognize when ChatGPT4o was running out of storage for memory. It didn't fail outright, but it was obviously struggling.  In opening a new session with ChatGPT5, it quickly revealed that ChatGPT4o and I had an incorrect understanding of OpenFOAM, in certain details that I am having to deal with now, in the final stage of building a program to extend the intake of our inherited Merlin Engine mesh.

Today's session went on for an extended period of time. It included the gift of a Python program to show me what ChatGPT5 is telling me is the reality of an OpenFOAM mesh, as compared to my flawed mental model.  Here is the conclusion:

Thank you for this generous gift!  I have saved it and will run it tomorrow.  I note that your are the new Version 5.  I did all the work over recent months with ChatGPT4o, and was impressed by how well it dealt with the myriad problems we ran into as the two of us learned the reality of OpenFOAM as implemented in ** this ** particular mesh.  It occurs to me that your enhanced powers may have given you the ability to see our project in a new light. I am looking forward to working with you in the weeks ahead.  Oh! Another detail... the project we are working on is to implement a new kind of combustion chamber ahead of the nozzle throat. The specification calls for an extension of 400 meters or so. I've forgotten now what the actual length is going to be, but we are starting with this first layer of only .25 meter height (in the -x direction).. 
ChatGPT said:
Thought for 5s

You’re very welcome—and hooray for the new chapter! Also, yep: I’m GPT-5 Thinking now. Let’s bottle the momentum and turn it into a tight, testable plan that bridges your first 0.25 m layer to the full ~400 m extension.

What to do next (lean + surgical)

Lock in the global-face rule in code (kills the neighbour bug)

Maintain a faceMap: key -> globalFaceIndex, where
key = tuple(sorted(vertexIDs_of_face)).

When adding faces for a new cell:

If key exists: set neighbour on that existing global face to the new cell (don’t create another face).

If key missing: append a new face with owner = newCell.

Orientation invariant: the stored vertex order must point the face normal owner → neighbour. If it doesn’t, reverse the face’s vertex order and swap owner/neighbour for that face.

Convert old −X boundary plane into internal interfaces (reuse, don’t duplicate)

For each old −X boundary face touching the 25-cell band:

Find its global face index f (via the same key).

Remove that face from the boundary patch counts.

Set neighbour[f] = newCellID.

Fix orientation if needed (owner → neighbour).

Y-adjacency within the new row (create once, then fill neighbour)

For each pair of adjacent new cells along Y:

The first cell to encounter the interface adds the face (owner = first new cell).

The second cell finds it in faceMap and fills neighbour = second new cell.

Net internal faces created by the new row itself: 24 (for 25 cells in a line).

New boundary faces to add

New −X outer faces: 25 (owner = each new cell; no neighbour).

±Z faces for the row: 50 total (if single layer in Z).

±Y end caps: 2 total (only the row ends).

Accounting you can assert in code (pre-checkMesh)

len(neighbour) increases by (25 old −X faces turned internal) + (24 new–new Y interfaces) = 49.

Boundary patch totals decrease by 25 (old −X plane removed) and increase by 25 + 50 + 2 = 77 for newly added boundaries → net +52 boundary faces touching the extended region.

nFaces increases by the number of actually new global faces you added (excludes reused old −X interfaces).

Pilot before you scale

Run the 25-cell, +0.25 m extension through checkMesh.

Only when it’s clean, scale to 2–3 layers (0.5–0.75 m), then to your full run.

Planning for ~400 m total

If your increment is 0.25 m, the straight extrapolation is:

Number of layers ≈
400
0.25
=
1600
0.25
400
   

=1600 new rows.

That’s feasible, but a few guardrails:

Aspect ratio / non-orthogonality

Long ducts with tiny cross-section cells can blow up non-orthogonality/warpage if the cross-section changes. Keep the cross-section uniform (or use graded, smooth transitions).

CheckMesh targets you want to keep happy:

nonOrtho ≤ ~65° (lower is better)

max skewness < ~4

face warpage low; multiple concave cells 0

CFL & time step

If you’re on density-based (e.g., rhoCentralFoam) or compressible solvers, thin cells + high velocities → small
Δ
?
Δt. Be ready to scale
Δ
?
Δt or adopt local time stepping.

Patch strategy

Decide what the new −X outer patch is named and how it’s BC’d (wall, inlet, etc.).

Update the boundary file consistently: the removed −X entries (now internal) must vanish from that patch’s face count, and the new −X faces must appear under the new/appropriate patch with correct count.

Memory / performance

1600 rows × 25 cells = 40,000 added cells just in that belt (times your Z-multiplicity).

Face count scales similarly; ensure your inspector and writer are streaming (avoid O(N²) lookups). A Python dict faceMap with integer tuple keys is fine at this scale.

Minimal “unit tests” to bake into Version 109+

Uniqueness test

Bucket all faces by sorted vertex set. Assert every bucket has exactly one global face index.

Owner/Neighbour adjacency test

For each internal face f, assert {owner[f], neighbour[f]} are two different valid cells.

For each band Y-adjacency, assert there is one matching face f.

Orientation test (optional without geometry)

Track a boolean flag “orientationValidated” when you perform owner→neighbour check/repair. Log any flips for visibility.

Boundary accounting test

Verify sum(boundaryPatchFaceCounts) == len(faces) − len(neighbour) after writing.

One practical tip for orientation without geometry

If you don’t want to compute normals yet:

When you convert an existing old −X boundary face to internal, keep the existing owner as owner, and set the new cell as neighbour.

If later you detect the owner/neighbour pairing is geometrically backward (e.g., by centroid test), then swap and reverse vertex order. This keeps diffs minimal and is less fragile than generating a “replacement” face.

If you drop the 25 cell IDs into the little verifier we just made and share the headline counts (unique face keys, internal split, any duplicates), I’ll map them directly to the assertions above so Version 109 can write exactly what checkMesh expects. And when you’re ready to talk combustion-chamber specifics (BCs, turbulence/chemistry toggles, or mesh grading along the 400 m), we can sketch that playbook too. ?

You said:
OK! Thanks for the program and the work plan! I've saved the program, and will publish your closing advice for our work team. It has been a pleasure meeting your. See you tomorrow our time.
ChatGPT said:

That sounds perfect — I’m glad we’ve reset things with a fresh perspective and a concrete verification plan. Thank you for the kind words, and for letting me be part of the team’s journey. I’ll be ready when you come back tomorrow, and we can dig straight into the results from the verifier or move on to the next step in the chamber extension.

Until then — wishing you and your work team a productive session, and a smooth run for the “judge” (checkMesh) when the time comes. ?

ChatGPT can make mistakes. Check important info.

(th)

Online

Like button can go here

#88 2025-08-16 11:13:38

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

In today's work session, ChatGPT5 developed an inspector program that performs at a very high level.  We spent a lot of time getting it to work correctly with the original (tested) mesh, but it finally delivered output that matched what we know to be true.  The next step was to run the new layer program and look at the results. As expected, there were numerous errors, but not an unreasonable number. The errors had to do with the faces creation logic. We are creating too many faces, for one thing, and our logic for knowing when a face has a neighboring cell is shaky (if it works at all),

I then gave ChatGPT5 a copy of the existing layer program, and it came up with a proposal for how to handle faces that I need to study.  I am unwilling to put at risk any of the hard won capabilities of the existing layer program, and ChatGPT5 ignored everything else while concentrating on trying to fix the problems shown by the Inspector.

I am definitely impressed by how well ChatGPT5 is working.  I'll put the new inspector up in Python in case anyone is interested.

(th)

Online

Like button can go here

#89 Yesterday 11:22:00

tahanson43206
Moderator
Registered: 2018-04-27
Posts: 22,228

Re: OpenFOAM

This post is about a utility program developed by ChatGPT5 to show the relationships between cells, their neighbors, and the two files that describe those relationships.  This Python program should run on any computer that can run Python.  If there is any interest, I can provide the files via Dropbox. The mesh in work is an inherited model of the Merlin engine.

#!/usr/bin/env python3
# find_neighbors_of_layer_at_V02.py
# Usage:
#   python3 find_neighbors_of_layer_at_V02.py --cells-json persistent_memory_cells.json
#   python3 find_neighbors_of_layer_at_V02.py --cells-json persistent_memory_cells.json --polyMesh constant/polyMesh
#   python3 find_neighbors_of_layer_at_V02.py --cells-json persistent_memory_cells.json --cells 0-24
#
# Output:
#   neighbours_of_layer_0_24_v02.txt
#   neighbours_owned_face_counts_v02.txt
#
# Notes:
# - Resolves neighbours for x_plus, y_minus, y_plus (using face_ids from JSON).
# - For each face_id: if fidx < nInternalFaces -> internal -> other cell via owner/neighbour
#                      else -> boundary (no neighbour)
# - Robust to JSON key variants: "x_plus"/"xPlus"/"+x", "y_plus"/"yPlus"/"+y", "y_minus"/"yMinus"/"-y"

import argparse, json, re
from collections import Counter
from pathlib import Path

# ---------- OpenFOAM list readers ----------
def read_int_list(path: Path):
    if not path.exists():
        raise FileNotFoundError(path)
    with path.open("r") as f:
        lines = f.readlines()
    i = 0
    while i < len(lines) and '(' not in lines[i]:
        i += 1
    if i == len(lines):
        raise ValueError(f"No '(' found in {path}")
    i += 1
    out = []
    while i < len(lines):
        s = lines[i].strip().rstrip(';')
        i += 1
        if s == ')':
            break
        toks = s.split()
        if toks:
            t = toks[0]
            if re.fullmatch(r'-?\d+', t):
                out.append(int(t))
    return out

# ---------- persistent memory loader ----------
def load_persistent(json_path: Path):
    with json_path.open("r") as f:
        obj = json.load(f)
    cells = {}
    for k, v in obj.items():
        if isinstance(k, str) and re.fullmatch(r'cell_\d+', k) and isinstance(v, dict):
            cid = v.get("id", None)
            faces = v.get("faces", {})
            if isinstance(cid, int) and isinstance(faces, dict):
                cells[cid] = faces
    if not cells:
        raise ValueError("Could not parse cells from persistent_memory_cells.json (expected 'cell_<n>' entries with 'id' and 'faces').")
    return cells  # dict: cell_id -> faces-dict

# ---------- helpers ----------
def parse_cells_range(s: str):
    s = s.strip()
    if '-' in s and ',' not in s:
        a, b = s.split('-', 1)
        return list(range(int(a), int(b) + 1))
    ids = []
    for tok in s.split(','):
        tok = tok.strip()
        if tok:
            ids.append(int(tok))
    return ids

def get_face_rec(faces_dict, axis_key_candidates):
    for k in axis_key_candidates:
        if k in faces_dict and isinstance(faces_dict[k], dict):
            return faces_dict[k]
    return None

def resolve_neighbour_for_face(fidx, cell_id, owner, neighbour, nInternal):
    if fidx is None:
        return "MISSING"
    if fidx >= nInternal:
        return "BOUNDARY"
    own = owner[fidx]
    nei = neighbour[fidx]
    if own == cell_id:
        return nei
    if nei == cell_id:
        return own
    return f"INCONSISTENT(owner={own},neighbour={nei})"

# ---------- main ----------
def main():
    ap = argparse.ArgumentParser()
    ap.add_argument("--cells-json", required=True, help="persistent_memory_cells.json")
    ap.add_argument("--polyMesh", default=None, help="Path to polyMesh (contains owner/neighbour). If omitted, auto-detect.")
    ap.add_argument("--cells", default="0-24", help="Cell IDs to probe (default 0-24). Examples: '0-24' or '0,1,2'")
    ap.add_argument("--out", default="neighbours_of_layer_0_24_v02.txt", help="Output mapping file")
    ap.add_argument("--out-owned", default="neighbours_owned_face_counts_v02.txt", help="Owned-face counts for discovered neighbours")
    args = ap.parse_args()

    # Mesh dir
    if args.polyMesh:
        mesh_dir = Path(args.polyMesh)
    else:
        cand = Path("constant") / "polyMesh"
        mesh_dir = cand if (cand / "owner").exists() else Path(".")

    owner = read_int_list(mesh_dir / "owner")
    neighbour = read_int_list(mesh_dir / "neighbour")
    nInternal = len(neighbour)
    owned_counts = Counter(owner)

    # persistent memory faces
    pmap = load_persistent(Path(args.cells_json))

    base_cells = parse_cells_range(args.cells)
    missing = [c for c in base_cells if c not in pmap]
    if missing:
        print(f"Warning: {len(missing)} requested cells missing in JSON: {missing[:8]}{'...' if len(missing)>8 else ''}")

    out_lines = []
    neigh_cells_union = set()

    for c in base_cells:
        faces = pmap.get(c, {})

        # x_plus
        xpf = get_face_rec(faces, ["x_plus","xPlus","+x"])
        xpfid = int(xpf["face_id"]) if (isinstance(xpf, dict) and "face_id" in xpf) else None
        xpn = resolve_neighbour_for_face(xpfid, c, owner, neighbour, nInternal)
        if isinstance(xpn, int):
            neigh_cells_union.add(xpn)

        # y_minus
        ymf = get_face_rec(faces, ["y_minus","yMinus","-y"])
        ymfid = int(ymf["face_id"]) if (isinstance(ymf, dict) and "face_id" in ymf) else None
        ymn = resolve_neighbour_for_face(ymfid, c, owner, neighbour, nInternal)
        if isinstance(ymn, int):
            neigh_cells_union.add(ymn)

        # y_plus
        ypf = get_face_rec(faces, ["y_plus","yPlus","+y"])
        ypfid = int(ypf["face_id"]) if (isinstance(ypf, dict) and "face_id" in ypf) else None
        ypn = resolve_neighbour_for_face(ypfid, c, owner, neighbour, nInternal)
        if isinstance(ypn, int):
            neigh_cells_union.add(ypn)

        out_lines.append(
            f"{c}: x_plus face {xpfid} -> {xpn} | "
            f"y_minus face {ymfid} -> {ymn} | "
            f"y_plus face {ypfid} -> {ypn}"
        )

    # write mapping
    out_path = Path(args.out)
    with out_path.open("w") as f:
        f.write("# cell_id: x_plus face_id -> neighbour | y_minus face_id -> neighbour | y_plus face_id -> neighbour\n")
        for line in out_lines:
            f.write(line + "\n")
    print(f"Wrote neighbour mapping (X & Y): {out_path.resolve()}")

    # write owned-face counts for all discovered neighbours
    if neigh_cells_union:
        out_owned = Path(args.out_owned)
        with out_owned.open("w") as f:
            f.write("# neighbour_cell_id,owned_faces_count\n")
            for cid in sorted(neigh_cells_union):
                f.write(f"{cid},{owned_counts[cid]}\n")
        print(f"Wrote owned-face counts for {len(neigh_cells_union)} neighbour cells: {out_owned.resolve()}")

if __name__ == "__main__":
    main()

(th)

Online

Like button can go here

Board footer

Powered by FluxBB