# DAGMC/Cubit Workflow Issues for Cylinders in OpenMC Plot

**URL:** <https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485>\
**Category:** User Support\
**Created:** [June 30, 2024, 5:36pm UTC](https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485 "2024-06-30T17:36:01Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![asunkari](https://avatars.discourse-cdn.com/v4/letter/a/838e76/32.png) [@asunkari](https://openmc.discourse.group/u/asunkari)\
**Post date:** [June 30, 2024, 5:36pm UTC](https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485/1 "2024-06-30T17:36:02Z")

</div>

_Let me first start off by saying I’m new to DAGMC. I’ve used GEANT4 for most of my work and I’ve imported CAD files using GDML. OpenMC’s DAGMC/h5m format is new to me._

_I have a couple of questions about merging and just general workflow guidance when it comes to importing STEP files and meshing them into a DAGMC format. I’ve gone through the webinars and tutorials for the DAGMC workflow, but I’m still running into some problems. I’ve outlined a thorough walkthrough of my workflow, so if anyone could please point out where the issue lies, it would be greatly appreciated!_

I am currently attempting to generate a rough “strawman” geometry of a molten salt reactor core. It’s very simplistic, at this point it is 3 cylinders. The three cylinders are coaxial, fitting precisely within one another. The innermost cylinder is the fuel, the middle cylinder is the moderator, and the outer cylinder is the container. As mentioned, they are coaxial and fit precisely within each other, so the inner surface(s) of each geometry is shared with the outer surface of the cylinder that fits within it. The geometry is shown in **Figure A.**

I first imprinted and merged all the relevant geometries/surfaces. This is where I run into problems. I first run the standard

```auto
imprint body all
merge body all

```

commands to automatically merge any shared surfaces. Here is the output:

```auto
Preserving undo information...done
Group imprint finished.
Updated volume(s): 2, 3
IMPRINT completed.
Finished Command: imprint body all

```

```auto
Consolidated 0 pair of surfaces
Curve 1 and 10 consolidated
Curve 2 and 9 consolidated
Curve 11 and 14 consolidated
Curve 12 and 13 consolidated
Curve 6 and 8 consolidated
Consolidated 5 curves
Consolidated 0 pairs of vertices
Finished Command: merge body all

```

Notice how none of the surfaces are actually merged. This is confirmed when I try to mesh the geometry, as I receive a warning that there are overlapping surfaces. The command

`find overlap all`

confirms this notion:

```auto
Finding surface overlap...
Found 6 overlapping surface pairs (added to group 'surf_overlap')
Surface 1 and 15 overlap
Surface 2 and 14 overlap
Surface 3 and 13 overlap
Surface 7 and 11 overlap
Surface 8 and 12 overlap
Surface 9 and 10 overlap
Finished Command: find overlap all

```

Inspecting these surfaces visually confirms that these are indeed surfaces that _should’ve been_ merged with the imprint and merge commands. For example, surfaces 3 and 13 (shared between the innermost and middle cylinders) do indeed represent the same surface when checked visually (as indicated by the orange circle in **Figure B** ). One can even see the artefacts from the outermost and middle container in **Figure B.**

In an attempt to solve this problem, I force merged all those overlapping surfaces.

```auto
merge surface 1 15 force
merge surface 2 14 force
…
merge surface 9 10 force

```

This indeed allowed me to mesh the geometry. However, I ran into more problems when trying to export the DAGMC file:

```auto
Found 7 entities of dimension 0
Found 7 entities of dimension 1
Found 9 entities of dimension 2
Found 3 entities of dimension 3
Surface 9 has forward sense with multiple volume 2 and volume 3
Error storing surface senses: 16
ERROR: Command Failed.

```

For multiple surfaces, I got this error regarding the forward sense of certain surfaces. After trying to do some research about the error online, I assumed that I must reverse it:

```auto
reverse surface 9
…
reverse surface 2

```

After reversing these surfaces, I attempt to mesh again, but get the overlapping surfaces error. Lo and behold, it’s the latest surface that I reversed:

```auto
Finding surface overlap...
Found 1 overlapping surface pairs
Surface 19 and 2 overlap
Finished Command: find overlap all

```

At this point I decide the force merge surfaces 2 and 19, which allows me to mesh the object. At this point I can also export the object as a DAGMC file. I didn’t mention this earlier, but the rest of the workflow (creating materials, assigning materials to blocks, creating reflecting boundary surfaces) have been completed.

I don’t know if my workflow is accurate, because **my OpenMC slice plot of the geometry looks completely wrong. I get some strange slices with the material colors all over the place, with the streak lines of the geometry meshing showing as well.** In addition, the slice coordinates seem to only work from between around 0.5 to -0.1, otherwise a blank image is outputted. This is strange, because the geometry is around 60 cm tall. It extends 10 cm in the +y direction from the origin, and extends 50 cm in the -y direction from the origin. **These can be seen in Figure C.**

**However, if I use the mbconvert tool from MOAB and convert the .h5m file into an STL file, the geometry looks exactly how I want it to look like.**

I have uploaded my OpenMC code, STEP files, and Cubit File here: [Filebin | yifmddkdvqm0to0q](https://filebin.net/yifmddkdvqm0to0q)  
_Note that the numbers of the surfaces may have changed, but the workflow is still the same when looking at the file._

 ![image](https://global.discourse-cdn.com/free1/uploads/openmc/original/2X/0/03b81fedb1e50dbe0e7d21193cceba5d09aa9c95.jpeg)

**Edit: added files**

---

<div class="post-metadata">

**Author:** ![asunkari](https://avatars.discourse-cdn.com/v4/letter/a/838e76/32.png) [@asunkari](https://openmc.discourse.group/u/asunkari)\
**Post date:** [June 30, 2024, 8:23pm UTC](https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485/2 "2024-06-30T20:23:00Z")

</div>

I have read this post by Dr. Shriwise regarding artefacts when plotting the slices in OpenMC:

> <https://github.com/openmc-dev/openmc/issues/2660>
>
> \*\*TL;DR:\*\* Artifacts sometimes appear in the plotter in regions with unmerged su…rfaces in DAGMC models. The order in which we check volumes/cells during the \`find\_cell\` operation matters in these cases and we should defer checking DAGMC's implicit complement until the end.
> 
> \---
> \## The problem
> 
> "Streaks" sometimes appear in plotting images of DAGMC models. As in our CSG geometry, this can occur if the slicing plane is on the boundary of a region, but these artifacts also appear throughout volumes in DAGMC universes. This is an explanation as to why this occurs and a proposal for a fix.
> 
> !\[image\](https://github.com/openmc-dev/openmc/assets/4563941/2cdc817b-e611-46e8-bb1e-2165ecc67533)
> 
> Image from the Copenhagen Atomics open \[MSRE model\](https://github.com/openmsr/msre)
> 
> 
> I can confirm that there are unmerged, coincident surfaces in this model in this region that are causing these problems. It's possible that other models have different causes for these artifacts, but here's my explanation of why it occurs for unmerged, coincident surfaces.
> 
> \## DAGMC model description
> 
> Let's say a model is created with the following two volumes, 1 and 2.
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/bed4820e-d7d3-48d4-b07a-7412e083400e"\>
> 
> 
> Ideally, there would be a single surface, surface 1, between the two volumes. It would have a "forward" sense with respect to volume 1 and a "reverse" sense with respect to volume 2. This sense is used to determine whether a particle crossing is entering or exiting one of the volumes.
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/0e540d7b-6562-4627-80c3-a5368d82ff23"\>
> 
> 
> This shared surface is generated as part of the "imprint and merge" process in the CAD engine. Part of this process is to combine, or merge, coincident surfaces into one -- resulting in surface 1 in this case. In the case that this process fails, the model has two coincident surfaces, surfaces 1 and 2, between volumes 1 and 2. In this case, we get the following:
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/c4bcd327-0e5e-407f-af21-3948b57b8134"\>
> 
> 
> The surfaces are shown with some separation to ease visualization, but they are numerically coincident with one another. Surface 1 now has a "forward" sense with respect to volume 1 and nothing on the other side. Surface 2 has a "forward" sense with respect to volume 2 as well and nothing on the other side. When building the \[implicit complement\](https://svalinn.github.io/DAGMC/contribute/devtheory.html#implicit-complement) in DAGMC, the programatically constructed space claims the "nothing" side of these surfaces. A key thing to remember is that \_both\_ of these surfaces are now part of the IPC volume generated by 
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/1c07b283-13e5-4bb5-b332-1a72f90704f1"\>
> 
> 
> At first glance, it might appear that the implicit complement occupies the space between the two unmerged surfaces, but this isn't quite how it works in practice. 
> 
> \### The implicit complement from the perspective of volume 1:
> 
> In this case the implicit complement is bounded by surface 1.
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/7ca5f588-9862-4e1d-bb12-a4909629fdf9"\>
> 
> \### Implicit complement from the perspective of volume 2:
> 
> In this case the implicit complement is bounded by surface 2.
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/54b30567-87be-4a70-9545-a96a4edb5089"\>
> 
> \## How this affects the plotting process:
> 
> In OpenMC we check each volume in a DAGMC universe for point containment to establish the position of a particle. To perform this check, we fire a ray in any given direction from the particle's location in an attempt to find an \_exiting\_ intersection with the triangles of the DAGMC volume.
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/0f2d9054-c7f1-49cd-ba3c-9fe1eb0075c2"\>
> 
> 
> As we get into the subtleties of the problem, we'll recall that these surfaces are coincident and represent them as such in the image:
> 
> \<img width="677" alt="image" src="https://github.com/openmc-dev/openmc/assets/4563941/4fb6c845-d9e5-47fe-96f4-44dd4adcba57"\>
> 
> While looping over volumes searching for a successful point containment check, it's possible that the implicit complement may be checked before volume 1 and 2. This is where the problem occurs. Because both surfaces 1 and 2 are part of the implicit complement and are coincident, a nearest intersection with either surface is valid and results in two outcomes:
> 
> 1. The particle ray intersects with surface 1: Based on the image of the IPC from the perspective of volume 1, this is considered to be an \_entering\_ intersection as if the particle were moving from volume 1 to the IPC volume and therefore the point is considered to be outside the IPC. The search loop continues and eventually identifies that the particle is in volume 1 (only the triangles of surface 2 are part of this search, so there's no ambiguity in that case)
> 2. The particle ray intersects with surface 2: Based on the image of the IPC from the perspective of volume 2, this is considered to be an exiting intersection as if the particle were moving from the IPC to volume 2. The particle is considered to be in the IPC and the search loop stops \*\*resulting in the artifacts we see in the plotter\*\*.
> 
> \### A proposed solution
> 
> Instead of including the IPC volume in the initial search for a containing cell, we can first check all of the intentionally created volumes in the CAD model and fall back to a check on the IPC at the end. In this way we still check all cells in the DAGMC universe, but we do not inadvertently exit the loop too soon and miss a particle placement. An implementation of this fix can be found \[here\](https://github.com/pshriwise/openmc/tree/dagmc-find-cell-updates) and does seem to remove these artifacts. Note there is some special behavior depending on whether or not the \`DAGUniverse\` is the root universe of the geometry or not.
> 
> !\[image\](https://github.com/openmc-dev/openmc/assets/4563941/7871beaa-b12a-4c07-8b2d-780f7b3fa55d)
> 
> 
> These changes requires the \`DAGUniverse::find\_cell\` method to take ownership of the search loop (before we could rely on the \`Universe::find\_cell\` routine). That was nice, but I don't see a problem with having a universe containing a fundamentally different geometry type have its own search loop.
> 
> \---
> 
> I thought I'd tackle a couple of questions that could come to mind reading the post above.
> \### Why does this occur in some locations of a volume but not others? Why is it direction-dependent?
> 
> The intersection with coincident surfaces of triangles is going to come down to traversal order of the bounding volume hierarchy and floating point arithmetic during the triangle intersection calculation -- slight changes in the angle and position will result in selecting triangles from one surface or the other of the coincident surfaces. From the perspective of the ray tracer, a nearest intersection with either surface is valid.
> 
> \### How might this adversely affect transport (after a particle is sourced)?
> 
> As long as the particle starts in the correct volume, the unmerged surfaces shouldn't have a negative effect during transport due to how intersections are filtered. During point containment checks, DAGMC searches for the nearest \_exiting or entering\_ intersection. During transport, DAGMC performs a query to find the nearest \_exiting\_ intersection. So in the example above as a particle moves out of volume 1 it will:
> 
> 1. Intersect with surface 1
> 2. Cross into the IPC volume
> 3. Intersect with surface 2 at a zero distance (intersections with surface 1 are filtered b/c they would be entering intersections)
> 4. Cross into volume 2
> 
> So while there are additional steps involved, the unmerged surfaces shouldn't cause a problem with transport.
> 
> \### Why didn't we see this in testing?
> 
> The test models were created as analogs of eigenvalue problems. As a result the geometries are simple and don't contain unmerged surfaces. Going forward, it's worth considering adding some models with unmerged surfaces to ensure we can handle these cases as expected -- though this is arguably something that should be part of DAGMC testing as opposed to a downstream application.

I think this may be related to the issue, as I had issues with the forward and reverse sense surfaces in Cubit. My OpenMC version is 0.15.0, so I should have the latest changes that reflect the improved plotting methodology. I can run the teapot example just fine and get a good plot:

![image](https://global.discourse-cdn.com/free1/uploads/openmc/original/2X/7/72d793ff9c1f2378513b8b44eb0e1ebd583872b5.png)

So I’m thinking it has to do with those meshing issues I had while working with Cubit and preparing the model for DAGMC.

---

<div class="post-metadata">

**Author:** ![asunkari](https://avatars.discourse-cdn.com/v4/letter/a/838e76/32.png) [@asunkari](https://openmc.discourse.group/u/asunkari)\
**Post date:** [July 1, 2024, 12:44am UTC](https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485/3 "2024-07-01T00:44:34Z")

</div>

Some more work on this topic:

I’ve made a simpler version of the geometry where the ends aren’t “capped” and the cylinders dont sit inside each other like Russian dolls. Instead, the geometry is “open faced” on both sides, like the pincell example.

This geometry seems to merge and mesh with no issues. Indeed, it was always the ends of each cylinder which were causing issues previously. The plot definitely looks more correct now. Additionally, I also realized that the issue with the “slice” layer value was because it seems that my geometry must’ve had a mismatch when exporting in terms of centimeters.

I still don’t know what the issue with the more complex geometry is. Furthermore, even after applying a reflecting boundary to all exterior surfaces of the model, I still see this “void” region (white) when plotting my geometry.

 ![image](https://global.discourse-cdn.com/free1/uploads/openmc/original/2X/6/6469cf317da74e9df5b4677e8f0d86c096c7120a.png)

---

<div class="post-metadata">

**Author:** ![nlinden](https://yyz2.discourse-cdn.com/free1/user_avatar/openmc.discourse.group/nlinden/32/2897_2.png) [@nlinden](https://openmc.discourse.group/u/nlinden)\
**Post date:** [July 1, 2024, 11:47am UTC](https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485/4 "2024-07-01T11:47:28Z")

</div>

Regarding your original problem, there seem to be a problem with the geometry contained in your step, as it seems Cubit fails because of surface senses.

I recreated your case in “pure Cubit” and it works fine:

```auto
reset
create Cylinder height 0.6 radius 0.15 
create Cylinder height 0.5 radius 0.1 
create Cylinder height 0.4 radius 0.05

subtract vol 2 from vol 1 keep_tool
subtract vol 3 from vol 2 keep_tool

create material "Inner" property_group "CUBIT-ABAQUS" 
create material "Outer" property_group "CUBIT-ABAQUS" 
create material "Inter" property_group "CUBIT-ABAQUS" 

set duplicate block elements off
undo group begin
block 1 add volume 3
block 1 name "Inner"
undo group end

set duplicate block elements off
undo group begin
block 2 add volume 4
block 2 name "Outer"
undo group end

set duplicate block elements off
undo group begin
block 3 add volume 5
block 3 name "Inter"
undo group end

block 1 material 'Inner'
block 2 material 'Outer'
block 3 material 'Inter'

merge tolerance 0.0005
imprint volume all
merge volume all

set trimesher coarse on ratio 100 angle 5
surface all scheme trimesh
mesh surface all

export cf_dagmc "C:/Users/NicolasLINDEN/Documents/cylinder.h5m" overwrite

```

How was your step file created? The problem probably stems from that.

From experience, imprinting and merging is the most error inducing part of the procedure for creating DAGMC geometries. I usually try to leverage the [implicit complement](https://svalinn.github.io/DAGMC/usersguide/codes/openmc.html#implicit-complement-materials) feature as much as I can, and I usally don’t have to do any I&M at all.

In your case, you could simply remove the intermediate cylinder and create an implicit material to replace it, this would eliminate the need for merging surfaces.

Another thing that I usually do and that may be useful to you, is that I always insert DAGMC geometries in OpenMC CSG cells, rather than directly creating a `openmc.Geometry` object from the DAGMCUniverse.

---

<div class="post-metadata">

**Author:** ![asunkari](https://avatars.discourse-cdn.com/v4/letter/a/838e76/32.png) [@asunkari](https://openmc.discourse.group/u/asunkari)\
**Post date:** [July 1, 2024, 8:51pm UTC](https://openmc.discourse.group/t/dagmc-cubit-workflow-issues-for-cylinders-in-openmc-plot/4485/5 "2024-07-01T20:51:09Z")

</div>

Thanks for the detailed help, @nlinden !

My STEP file was exported from my workspace in OnShape. I’ve tried exporting individual parts, part studios, and assemblies, but it still shows up the same way. I haven’t tried “pure cubit” as you suggested, but that’s probably a cleaner approach since the geometry isn’t that complicated.

I do wonder, however, if the extrusion direction in Onshape makes a difference? The only reason why I could think that the surface faces the wrong way is because I extruded it in the “opposite” direction when creating the geometry in the CAD program.

For now I think I’ll use the “Pure Cubit” geometry, or even just the “built-in” OpenMC geometry for simplicity while I try to figure out the CAD model. I’ll be sure to try the implicit complement feature for the intermediate cylinder, too, and update the thread on how that goes.

Thanks again for the help!

**EDIT:** Here is the link to the OnShape document, if anyone is interested in exporting the STEP file to see if they have better luck. The file must be exported in centimeters. [Onshape](https://cad.onshape.com/documents/7f449c1388a7ab813aa3e27a/w/8ac99a4f1b1c7d966d060ae8/e/0dc1ca07d7b17630f3b78089?renderMode=0&uiState=66831709cffd54186b48e843)
