Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

cad_geometry_example #29

Closed
wants to merge 25 commits into from
Closed
Show file tree
Hide file tree
Changes from 5 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
# CAD-based Geometry Workflow for Multiphysics Fusion Problems Using OpenMC and MOOSE

This demonstration describes a workflow for modeling fusion problems in OpenMC and MOOSE using a computer aided design (CAD)-based geometry workflow.
It is based on the work published in [!cite](Eltawila2004_PNBC).
meltawila marked this conversation as resolved.
Show resolved Hide resolved

!figures transfers.png
id=transfers
caption=OpenMC and MOOSE Coupling
style=width:60%;margin-left:auto;margin-right:auto

In this example, you'll learn how to:

- Couple OpenMC and MOOSE using Cardinal for fixed source Monte Carlo calculations.
- Use Cardinal to tally values of interest such as tritium production and heating which would be used in MOOSE to solve for the temperature distribution

An extremely simplified tokamak was modeled in CAD and was considered for this example. The meshed geometry was prepared using direct accelerated geometry Monte Carlo (DAGMC) for particle transport, and a volumetric mesh was also prepared to be used in MOOSE’s finite element solver and to tally OpenMC results for heat source distribution and tritium production. Cardinal was used to run OpenMC Monte Carlo particle transport within MOOSE framework. The data transfer system transfered heat source and temperature distribution between OpenMC and MOOSE, with coupling between neutron transport and heat conduction achieved via Picard iteration.

## Generating the meshes

The CAD model was first developed in FUSION360 and was imported into Cubit to assign blocks, materials, and side sets and generate the mesh (tmesh_1.e). A corresponding DAGMC surface mesh (tmesh_1.h5m) was exported directly from the meshed geometry in Cubit (by loading the volumetric meshed geometry in Cubit and exporting a DAGMC surface mesh).

In this example, `tmesh_1.e` is the finite element mesh used in MOOSE on which the heat conduction physics is solved. `tmesh_1.h5m` is the DAGMC surface mesh used for particle transport in OpenMC (which bounds the surfaces between different materials). Cardinal also allows for mesh tallying for tallying OpenMC results directly on the mesh overlayed on the OpenMC geometry which `tmesh_1.e` could be used for as well as an unstructered volume mesh. This could be used by changing the tally type and adding a mesh template (`tally_type = mesh`, `mesh_template = tmesh_1.e`).

!figures mesh_1.png
id=volumetric_mesh
caption=Volumetric mesh [!citep](Eltawila2004_PNBC).
meltawila marked this conversation as resolved.
Show resolved Hide resolved
style=width:60%;margin-left:auto;margin-right:auto

!figures d1.png
id=dagmc
caption=DAGMC surface mesh [!citep](Eltawila2004_PNBC).
meltawila marked this conversation as resolved.
Show resolved Hide resolved
style=width:60%;margin-left:auto;margin-right:auto

## OpenMC

!listing /input_files/model.py language=python

## Cardinal

!listing /input_files/openmc.i

## MOOSE Heat transfer
meltawila marked this conversation as resolved.
Show resolved Hide resolved

!listing /input_files/solid.i

## Results

!figures Temps.png
id=temps
caption=Temperature distribution [!citep](Eltawila2004_PNBC).
meltawila marked this conversation as resolved.
Show resolved Hide resolved
style=width:60%;margin-left:auto;margin-right:auto

!figures tritium_production.png
id=h3production
caption=Tritium production rate density [!citep](Eltawila2004_PNBC).
meltawila marked this conversation as resolved.
Show resolved Hide resolved
style=width:60%;margin-left:auto;margin-right:auto

!table id=results caption=Results summary
| Armor Max. Temp. [K]| 1062.4 |
| First Wall Max. Temp. [K]| 1057.6 |
| Breeder Max. Temp. [K]| 987.4 |
| Heat Source [W] | 2.44 × 10^5 ± 3 × 10^3 |
| Tritium Production [atoms/s] | 4.70 × 10^13 ± 8 × 10^11 |

Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Original file line number Diff line number Diff line change
@@ -0,0 +1,106 @@
import openmc
import math

mat1 = openmc.Material(name="mat1")
mat1.set_density('g/cc', 19.30)
mat1.add_element('W', 1.0)

eurofer = openmc.Material(name="eurofer")
eurofer.add_element('Fe', 0.011 , 'wo')
eurofer.add_element('Al', 0.002 , 'wo')
eurofer.add_element('As', 0.0002 , 'wo')
eurofer.add_element('B', 0.0012 , 'wo')
eurofer.add_element('C', 0.0005 , 'wo')
eurofer.add_element('Co', 0.0005 , 'wo')
eurofer.add_element('Cr', 0.0005 , 'wo')
eurofer.add_element('Cu', 0.00005 , 'wo')
eurofer.add_element('Mn', 0.00005 , 'wo')
eurofer.add_element('Mo', 0.0001 , 'wo')
eurofer.add_element('N', 0.0001 , 'wo')
eurofer.add_element('Nb', 0.00005 , 'wo')
eurofer.add_element('Ni', 0.0003 , 'wo')
eurofer.add_element('O', 0.00005 , 'wo')
eurofer.add_element('P', 0.004 , 'wo')
eurofer.add_element('S', 0.0001 , 'wo')
eurofer.add_element('Sb', 0.09 , 'wo')
eurofer.add_element('Sn', 0.0001 , 'wo')
eurofer.add_element('Si', 0.0011 , 'wo')
eurofer.add_element('Ta', 0.00002 , 'wo')
eurofer.add_element('Ti', 0.0005 , 'wo')
eurofer.add_element('V', 0 , 'wo')
eurofer.add_element('W', 0.0001 , 'wo')
eurofer.add_element('Zr', 0.88698 , 'wo')
eurofer.set_density("g/cm3", 7.798)

Helium = openmc.Material(name="Helium")
Helium.add_element('He', 1.0)
Helium.set_density("kg/m3", 0.166)

mat2 = openmc.Material.mix_materials([eurofer, Helium], [0.65, 0.35], 'ao')

beryllium = openmc.Material(name="beryllium")
beryllium.add_element('Be', 1.0)
beryllium.set_density("g/cm3", 1.85)

Li4SiO4 = openmc.Material(name="Li4SiO4")
Li4SiO4.add_element('Li', 4.0)
Li4SiO4.add_element('Si', 1.0)
Li4SiO4.add_element('O', 4.0)
Li4SiO4.set_density("g/cm3", 2.39)

mat3 = openmc.Material.mix_materials([eurofer, beryllium, Li4SiO4, Helium], [0.1, 0.37, 0.15, 0.38], 'ao')

mats = openmc.Materials([mat1, eurofer, Helium, mat2, beryllium, Li4SiO4, mat3])
mats.export_to_xml()

pz = openmc.Plot()
pz.basis = 'yz'
pz.origin = (0.0, 0.0, 0.0)
pz.width = (200.0, 200.0)
pz.pixels = (500, 500)
pz.color_by = 'material'

px = openmc.Plot()
px.basis = 'xy'
px.origin = (0.0, 0.0, 0.0)
px.width = (200, 200)
px.pixels = (500, 500)
px.color_by = 'material'

plots = openmc.Plots([pz,px])
plots.export_to_xml()



settings = openmc.Settings()
settings.dagmc = True
settings.batches = 100
settings.particles = 10000000
settings.run_mode = "fixed source"

settings.temperature = {'default': 800.0,
'method': 'interpolation',
'range': (294.0, 3000.0),
'tolerance': 1000.0}

source = openmc.Source()

r = openmc.stats.PowerLaw(55, 65, 1.0)
phi = openmc.stats.Uniform(0.0, 2*math.pi)
z = openmc.stats.Discrete([0,], [1.0,])
spatial_dist = openmc.stats.CylindricalIndependent(r, phi, z)

source.angle = openmc.stats.Isotropic()
source.energy = openmc.stats.Discrete([14.08e6], [1.0])
source.space=spatial_dist
settings.source = source
settings.export_to_xml()




dagmc_univ = openmc.DAGMCUniverse(filename='torus7v2t2.h5m')


geometry = openmc.Geometry(root=dagmc_univ)
geometry.export_to_xml()
Original file line number Diff line number Diff line change
@@ -0,0 +1,83 @@
[Mesh]
[file]
type = FileMeshGenerator
file = tmesh_1.e
[]
[]

[AuxVariables]
[cell_temperature]
family = MONOMIAL
order = CONSTANT
[]
[]

[AuxKernels]
[cell_temperature]
type = CellTemperatureAux
variable = cell_temperature
[]
[]

[Problem]
type = OpenMCCellAverageProblem
tally_type = cell
tally_name = 'heat_source H3'
lowest_cell_level = 0
temperature_blocks = '1 2 3'
check_tally_sum = false
source_strength = 1e18
volume_calculation = vol
tally_score = 'heating_local H3_production'
tally_trigger = 'rel_err none'
tally_trigger_threshold = '0.1 0.1'
verbose = true
max_batches = 10
batch_interval = 5
particles = 5000
output = unrelaxed_tally_std_dev
skinner = moab
[]

[UserObjects]
[vol]
type = OpenMCVolumeCalculation
n_samples = 5000
[]
[moab]
type = MoabSkinner
temperature_min = 800
temperature_max = 1100
n_temperature_bins = 100
temperature = temp
build_graveyard = true
[]
[]

[Postprocessors]
[heat_source]
type = ElementIntegralVariablePostprocessor
variable = heat_source
[]
[tritium_production]
type = ElementIntegralVariablePostprocessor
variable = H3
[]
[tritium_RelativeError]
type = TallyRelativeError
tally_score = h3_production
[]
[heat_source_RelativeError]
type = TallyRelativeError
tally_score = heating_local
[]
[]

[Executioner]
type = Transient
[]

[Outputs]
exodus = true
csv = true
[]
Original file line number Diff line number Diff line change
@@ -0,0 +1,112 @@
[Mesh]
[file]
type = FileMeshGenerator
file = tmesh_1.e
[]
[]

[Variables]
[temp]
initial_condition = 800.0
[]
[]

[AuxVariables]
[heat_source]
family = MONOMIAL
order = CONSTANT
[]
[]

[Kernels]
[hc]
type = HeatConduction
variable = temp
[]
[heat]
type = CoupledForce
variable = temp
v = heat_source
[]
[]

[BCs]
[surface]
type = DirichletBC
variable = temp
boundary = 1
value = 800.0
[]
[]

[Materials]
[k_1]
type = GenericConstantMaterial
prop_values = '1.64'
prop_names = 'thermal_conductivity'
block = 'Armour'
[]
[k_2]
type = GenericConstantMaterial
prop_values = '0.45'
prop_names = 'thermal_conductivity'
block = 'FW'
[]
[k_3]
type = GenericConstantMaterial
prop_values = '0.65'
prop_names = 'thermal_conductivity'
block = 'BM'
[]
[]

[Executioner]
type = Transient
nl_abs_tol = 1e-8
num_steps = 2
solve_type = 'NEWTON'
[]

[Outputs]
exodus = true
print_linear_residuals = false
perf_graph = true
[]

[MultiApps]
[openmc]
type = TransientMultiApp
app_type = CardinalApp
input_files = 'openmc.i'
execute_on = timestep_end
[]
[]

[Transfers]
[heat_source_from_openmc]
type = MultiAppGeneralFieldShapeEvaluationTransfer
from_multi_app = openmc
variable = heat_source
source_variable = heat_source
from_postprocessors_to_be_preserved = heat_source
to_postprocessors_to_be_preserved = source_integral
[]
[temp_to_openmc]
type = MultiAppGeneralFieldShapeEvaluationTransfer
to_multi_app = openmc
variable = temp
source_variable = temp
[]
[]

[Postprocessors]
[source_integral]
type = ElementIntegralVariablePostprocessor
variable = heat_source
execute_on = transfer
[]
[max_T]
type = NodalExtremeValue
variable = temp
[]
[]
Binary file not shown.