-
Notifications
You must be signed in to change notification settings - Fork 7
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
Adding Kobayashi Dogleg model #25
base: master
Are you sure you want to change the base?
Conversation
Removing large files Fixing Problem_II New results with modified Problem II
Hi Ethan, we spoke a few months ago about this project coming along this summer. Any comments or suggested changes you might have will be handled by one or both of us. I wanted Nina to issue the pull request, though, since she did most of the hard work, and I just formatted it to conform with your other models. |
def _parse_args(): | ||
"""Parse and return commandline arguments""" | ||
parser = argparse.ArgumentParser() | ||
parser.add_argument("-b", "--batches", type=int, default=20) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
parser.add_argument("-b", "--batches", type=int, default=20) | |
parser.add_argument("-b", "--batches", type=int, default=100) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
We usually leave 100 batches as default just for convention.
parser = argparse.ArgumentParser() | ||
parser.add_argument("-b", "--batches", type=int, default=20) | ||
parser.add_argument("-p", "--particles", type=int, default=int(1e7)) | ||
parser.add_argument("-s", "--threads", type=int,default=5) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
parser.add_argument("-s", "--threads", type=int,default=5) | |
parser.add_argument("-s", "--threads", type=int) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Just by convention, we prefer not to specify any default number of threads in the other benchmarks.
Source = openmc.Material(name='source_mat') | ||
Source.set_density('macro', 1.0) | ||
Source.add_macroscopic('source') | ||
|
||
Absorber = openmc.Material(name='absorber_mat') | ||
Absorber.set_density('macro', 1.0) | ||
Absorber.add_macroscopic('absorber') | ||
|
||
Void = openmc.Material(name='void_mat') | ||
Void.set_density('macro', 1.0) | ||
Void.add_macroscopic('void') | ||
|
||
# Instantiate a Materials collection and export to XML | ||
materials = openmc.Materials([Source,Absorber,Void]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Source = openmc.Material(name='source_mat') | |
Source.set_density('macro', 1.0) | |
Source.add_macroscopic('source') | |
Absorber = openmc.Material(name='absorber_mat') | |
Absorber.set_density('macro', 1.0) | |
Absorber.add_macroscopic('absorber') | |
Void = openmc.Material(name='void_mat') | |
Void.set_density('macro', 1.0) | |
Void.add_macroscopic('void') | |
# Instantiate a Materials collection and export to XML | |
materials = openmc.Materials([Source,Absorber,Void]) | |
source = openmc.Material(name='source_mat') | |
source.set_density('macro', 1.0) | |
source.add_macroscopic('source') | |
absorber = openmc.Material(name='absorber_mat') | |
absorber.set_density('macro', 1.0) | |
absorber.add_macroscopic('absorber') | |
void = openmc.Material(name='void_mat') | |
void.set_density('macro', 1.0) | |
void.add_macroscopic('void') | |
# Instantiate a Materials collection and export to XML | |
materials = openmc.Materials([source, absorber, void]) |
source_cell = openmc.Cell(fill = Source, | ||
region = box2, | ||
name = 'source region') | ||
|
||
void_cell = openmc.Cell(fill = Void, | ||
region = (segment1 | segment2 | -segment3 | -segment4), | ||
name = 'void region') | ||
|
||
absorber_cell = openmc.Cell(fill = Absorber, | ||
region = box0 & ~void_cell.region & ~source_cell.region, | ||
name = 'absorber region') | ||
|
||
|
||
# Fill the Cell with the Material | ||
Ucell.fill = Void |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
source_cell = openmc.Cell(fill = Source, | |
region = box2, | |
name = 'source region') | |
void_cell = openmc.Cell(fill = Void, | |
region = (segment1 | segment2 | -segment3 | -segment4), | |
name = 'void region') | |
absorber_cell = openmc.Cell(fill = Absorber, | |
region = box0 & ~void_cell.region & ~source_cell.region, | |
name = 'absorber region') | |
# Fill the Cell with the Material | |
Ucell.fill = Void | |
source_cell = openmc.Cell(fill = source, | |
region = box2, | |
name = 'source_cell') | |
void_cell = openmc.Cell(fill = void, | |
region = (segment1 | segment2 | -segment3 | -segment4), | |
name = 'void_cell') | |
absorber_cell = openmc.Cell(fill = absorber, | |
region = box0 & ~void_cell.region & ~source_cell.region, | |
name = 'absorber_cell') | |
# Fill the Cell with the Material | |
ucell.fill = void |
Ucell = openmc.Cell(name='base universe', region = box) | ||
|
||
# domain: the domain for spatial homogenization aka | ||
# averaging the cross sections over volume and flux spectrum | ||
|
||
total = mgxs.TotalXS(domain=Ucell, energy_groups=groups) | ||
absorption = mgxs.AbsorptionXS(domain=Ucell, energy_groups=groups) | ||
scattering = mgxs.ScatterXS(domain=Ucell, energy_groups=groups) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm not convinced by starting with a capital letter for naming the objects.
Ucell = openmc.Cell(name='base universe', region = box) | |
# domain: the domain for spatial homogenization aka | |
# averaging the cross sections over volume and flux spectrum | |
total = mgxs.TotalXS(domain=Ucell, energy_groups=groups) | |
absorption = mgxs.AbsorptionXS(domain=Ucell, energy_groups=groups) | |
scattering = mgxs.ScatterXS(domain=Ucell, energy_groups=groups) | |
ucell = openmc.Cell(name='base_universe', region = box) | |
# domain: the domain for spatial homogenization aka | |
# averaging the cross sections over volume and flux spectrum | |
total = mgxs.TotalXS(domain=ucell, energy_groups=groups) | |
absorption = mgxs.AbsorptionXS(domain=ucell, energy_groups=groups) | |
scattering = mgxs.ScatterXS(domain=ucell, energy_groups=groups) |
# generate weight windows tallies | ||
try: | ||
p = subprocess.Popen(["python3", "openmc_model.py", "--Problem_I","--weight_windows"]) | ||
except: | ||
p = subprocess.Popen(["python", "openmc_model.py", "--Problem_I","--weight_windows"]) | ||
|
||
# wait for the simulation to finish | ||
p.wait() | ||
|
||
# run simulation without scattering | ||
try: | ||
p = subprocess.Popen(["python3", "openmc_model.py", "--Problem_I"]) | ||
except: | ||
p = subprocess.Popen(["python", "openmc_model.py", "--Problem_I"]) | ||
|
||
# wait for the simulation to finish | ||
p.wait() | ||
|
||
try: | ||
p = subprocess.Popen(["python3", "openmc_model.py", "--weight_windows"]) | ||
except: | ||
p = subprocess.Popen(["python", "openmc_model.py", "--weight_windows"]) | ||
|
||
# wait for the simulation to finish | ||
p.wait() | ||
|
||
# run simulation with scattering | ||
try: | ||
p = subprocess.Popen(["python3", "openmc_model.py"]) | ||
except: | ||
p = subprocess.Popen(["python", "openmc_model.py"]) | ||
|
||
# wait for the simulation to finish | ||
p.wait() |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I would not default to a ww generation run. Especially for the run_and_store
file. I would rather provide already some ww files and leave the ww generation as optional. Or does this benchmark specifically focus on ww generation and performances?
# generate weight windows tallies | |
try: | |
p = subprocess.Popen(["python3", "openmc_model.py", "--Problem_I","--weight_windows"]) | |
except: | |
p = subprocess.Popen(["python", "openmc_model.py", "--Problem_I","--weight_windows"]) | |
# wait for the simulation to finish | |
p.wait() | |
# run simulation without scattering | |
try: | |
p = subprocess.Popen(["python3", "openmc_model.py", "--Problem_I"]) | |
except: | |
p = subprocess.Popen(["python", "openmc_model.py", "--Problem_I"]) | |
# wait for the simulation to finish | |
p.wait() | |
try: | |
p = subprocess.Popen(["python3", "openmc_model.py", "--weight_windows"]) | |
except: | |
p = subprocess.Popen(["python", "openmc_model.py", "--weight_windows"]) | |
# wait for the simulation to finish | |
p.wait() | |
# run simulation with scattering | |
try: | |
p = subprocess.Popen(["python3", "openmc_model.py"]) | |
except: | |
p = subprocess.Popen(["python", "openmc_model.py"]) | |
# wait for the simulation to finish | |
p.wait() | |
# run simulation without scattering | |
try: | |
p = subprocess.Popen(["python3", "openmc_model.py", "--Problem_I"]) | |
except: | |
p = subprocess.Popen(["python", "openmc_model.py", "--Problem_I"]) | |
# wait for the simulation to finish | |
p.wait() | |
# run simulation with scattering | |
try: | |
p = subprocess.Popen(["python3", "openmc_model.py"]) | |
except: | |
p = subprocess.Popen(["python", "openmc_model.py"]) | |
# wait for the simulation to finish | |
p.wait() |
print("Reading Statepoint file") | ||
# read statepoint file | ||
ProblemI_file = ofb.ResultsFromOpenmc('statepoint.20.h5', 'Problem_I') | ||
# read statepoint file | ||
ProblemII_file = ofb.ResultsFromOpenmc('statepoint.20.h5', 'Problem_II') |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
print("Reading Statepoint file") | |
# read statepoint file | |
ProblemI_file = ofb.ResultsFromOpenmc('statepoint.20.h5', 'Problem_I') | |
# read statepoint file | |
ProblemII_file = ofb.ResultsFromOpenmc('statepoint.20.h5', 'Problem_II') | |
print("Reading Statepoint file") | |
# read statepoint file | |
problemI_file = ofb.ResultsFromOpenmc('statepoint.100.h5', 'Problem_I') | |
# read statepoint file | |
problemII_file = ofb.ResultsFromOpenmc('statepoint.100.h5', 'Problem_II') |
Hi @ninadol and @sjcrs, thank you so much for contributing to this! The Kobayashi benchmark is a great addition to the repository! I'll leave a couple of comments here that we think should help better homogenize this benchmark with the other ones, if you want to address them: Style:
I would like to open a discussion about using the magic method in the default run. I believe we can either provide the ww file already and foget about it or maybe having the ww generation run as optional with the Lastly, and most important, did you happen to have the time to explore why openmc is giving results that are so different from mcnp? Thank you again for the great work! |
Hi @SteSeg , thank you for your review and for your suggestions. As you said, these are mostly for style conformity and therefore we have no objection to making these changes. The goal of @ninadol 's work was to benchmark the weight windows generation using the MAGIC method, and the goal with this problem was to demonstrate its utility. I had thought about structuring this so that it first produced an analog run with large error bars before moving on to the MAGIC method weight windows generation, but was concerned about producing too many comparisons. I think it's reasonable to provide a ww file and optinally generate a new one via Finally, I have spoken to @jtramm about the differences between the results here, because i know he has been looking at the Kobayashi dogleg for a different variance reduction feature. Movement on this PR will encourage us to meet and puzzle this one out. |
Actually, I'd like to correct myself, |
Thank you for the clarification. I think the best way to go is to provide an original ww file, then add the ww generation as optional and have just two notebooks. One for the actual results analysis ( Anyways, I think we should dig more in the results differences, usually openmc and mcnp results are not so far apart. I'll run the benchmark soon just to add another couple of eyes to search for the problem. |
I worked with Dr. Stephen Coleman at Radiasoft LLC to construct the Kobayashi Dogleg Benchmark for problem 3i and 3ii. Stephen formatted my initial model for this package. The relevant plots are included in the post-processing notebook.
This is novel as a benchmark because it utilizes the new MAGIC Method implementation in OpenMC. This is also novel because it is outside the Sinbad database. Stephen will connect below, and to help with any feedback you might have for us.