fusion-energy / paramak

Create parametric 3D fusion reactor CAD and neutronics models
https://paramak.readthedocs.io/en/main/
MIT License
58 stars 17 forks source link

Dagmc export of shapes producted from the EU DEMO example #306

Open Adambar24 opened 1 year ago

Adambar24 commented 1 year ago

Hi, I tired to export_dagmc_h5m the Reactor example Eudemofrom2105paper. The export could not complete due to a 'no matching part error' in the brep_part_finder function. Example code below:

import paramak

my_reactor = paramak.EuDemoFrom2015PaperDiagram( rotation_angle=180.0, number_of_tf_coils=16, )

my_reactor.export_dagmc_h5m(filename='demotest.h5m')

I tired copying the shapes of the EU demo reactor class and combining them into another reactor class separately i,e, i copied the coodinates points of the blanket, divertor, vac vessel, vac vessel interior and tf coil casting. When I instance a reactor class(without the vac vessel interior, since this is only used to create the vac vessel shape) with these shapes and try to export_dagmc_h5m. I get the same error. Though if i instance the reactor class without the vac vessel it will create a h5m file and also if i instance the reactor class with only the vac vessel shape it will create a h5m. The problem only arises with all shapes combined. Example these two individual reactors would create two .h5m files.

my_reactor = paramak.Reactor([blanket, divertor,tf_coil_casing]) vac_vessel = paramak.Reactor([vac_vessel]) my_reactor.export_dagmc_h5m(filename='test.h5m', min_mesh_size = 5, max_mesh_size = 20) vac_vessel.export_dagmc_h5m(filename='test2.h5m', min_mesh_size = 5, max_mesh_size = 20)

Whereas:

my_reactor = paramak.Reactor([blanket, divertor,tf_coil_casing,vac_vessel]) my_reactor.export_dagmc_h5m(filename='test.h5m', min_mesh_size = 5, max_mesh_size = 20)' Would not create one h5m file

shimwell commented 1 year ago

Sorry about this, tracking the volumes through the workflow is tricky stuff. I have a longer term plan to attach meta data to the brep file so that we can track volumes and won't have to search for them using matching parts. This is a while away from being complete but it my current focus to really fixing this problem instead of just patching it up. I can take a look at preproducing this error and see if there are any easy fixes, thanks for reporting and making it easy to reproduce

Adambar24 commented 1 year ago

Not a problem, I am just mainly trying to gauge which method is best to use in the future. Since I want to start building some more complex reactors that are likely to reuse some of the parts from the EuDemo reactor. I will try with the cad_to_h5m method to see if I can get that working. If not I will use the stl_to_h5m. The problem with the later is that I dont yet understand how to insure my geometry is not overlapping since I havent tried this method yet. Also it seems the most manual approach so I will try this last.

Though in the future it would be better using dagmc_to_h5m since its open source rather than relying on cubit

shimwell commented 1 year ago

I think I've found the bug. Not fixed but I know what is going wrong.

For some reason the Brep export in the cad-to-dagmc package is not working for this geometry. Oddly enough the paramak.export_brep() method does work. They are quite similar so I should be able to track it down

So perhaps in the meantime you could export the Brep using the paramak then use Brep-to-h5m to export it as a h5m.

Material assignment will be a bit tricky till I can get the cad-to-dagmc brep export function working for this geometry.

Adambar24 commented 1 year ago

Hi, thanks I tried the work around and the paramak.export_brep() method. I then used this file in the Brep-to-h5m to create the h5m file as suggested. Here its easier to add tags to the materials (see below) so its easier to check them later.

brep_to_h5m(
    brep_filename='test.brep',
    material_tags=[
        'Blanket',
        'divertor',
        'tf coil casing',
        'vacuum vessel',
    ],
    h5m_filename='dagmc.h5m',
    min_mesh_size= 5,
    max_mesh_size = 20,
    mesh_algorithm = 1,
)

After which you can use the property material_names found in the openmc.DAGMCUniverse class to print the material names in the set in the .h5m file.

dag_univ = openmc.DAGMCUniverse(filename='dagmc.h5m')
vac_surf = openmc.Sphere(r=10000, surface_id=9999, boundary_type="vacuum")

print(dag_univ.material_names)

This allows you to assign the material names easily for the openmc simulation. I thought I would be explicit if somebody was having similar problems until the dagmc to h5m is fixed. Thanks for your help

shimwell commented 1 year ago

Thanks @Adambar24 and congrats for completing the work around. Posting the example is much appreciated, I'm sure others will find it useful.

by the way I added an automatic bounded_universe method to openmc a while back that you can use to save one line if you are interested.

~vac_surf = openmc.Sphere(r=10000, surface_id=9999, boundary_type="vacuum")~

dag_univ = openmc.DAGMCUniverse(filename='dagmc.h5m').bounded_universe()

I shall get back to you on this thread when I make progress on fixing the workflow