crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-26 21:25:36 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Allow materials to be created/assigned from C API bindings. Fix a few errors.

Browse source
This commit is contained in:
Paul Romano 2017-09-27 15:44:40 -05:00
parent a1eb39d8fc
commit bc84687378
8 changed files with 219 additions and 15 deletions
Download patch file Download diff file Expand all files Collapse all files

38
openmc/capi/cell.py
View file

@ -1,12 +1,14 @@
from collections import Mapping
from collections import Mapping, Iterable
from ctypes import c_int, c_int32, c_double, c_char_p, POINTER
from weakref import WeakValueDictionary
import numpy as np
from numpy.ctypeslib import as_array
from . import _dll
from .core import _View
from .error import _error_handler
from .material import MaterialView
__all__ = ['CellView', 'cells']
@ -14,6 +16,12 @@ __all__ = ['CellView', 'cells']
_dll.openmc_cell_get_id.argtypes = [c_int32, POINTER(c_int32)]
_dll.openmc_cell_get_id.restype = c_int
_dll.openmc_cell_get_id.errcheck = _error_handler
_dll.openmc_cell_get_fill.argtypes = [
c_int32, POINTER(c_int), POINTER(POINTER(c_int32)), POINTER(c_int32)]
_dll.openmc_cell_set_fill.argtypes = [
c_int32, c_int, c_int32, POINTER(c_int32)]
_dll.openmc_cell_set_fill.restype = c_int
_dll.openmc_cell_set_fill.errcheck = _error_handler
_dll.openmc_cell_set_temperature.argtypes = [
c_int32, c_double, POINTER(c_int32)]
_dll.openmc_cell_set_temperature.restype = c_int
@ -55,6 +63,34 @@ class CellView(_View):
_dll.openmc_cell_get_id(self._index, cell_id)
return cell_id.value
@property
def fill(self):
fill_type = c_int()
indices = POINTER(c_int32)()
n = c_int32()
_dll.openmc_cell_get_fill(self._index, fill_type, indices, n)
if fill_type.value == 1:
if n.value > 1:
return [MaterialView(i) for i in indices[:n.value]]
else:
return MaterialView(indices[0])
else:
raise NotImplementedError
@fill.setter
def fill(self, fill):
if isinstance(fill, Iterable):
n = len(fill)
indices = (c_int*n)(*(m._index for m in fill))
_dll.openmc_cell_set_fill(self._index, 1, 1, indices)
elif isinstance(fill, MaterialView):
materials = [fill]
indices = (c_int*1)(fill._index)
_dll.openmc_cell_set_fill(self._index, 1, 1, indices)
else:
raise NotImplementedError
def set_temperature(self, T, instance=None):
"""Set the temperature of a cell

2
openmc/capi/error.py
View file

@ -62,7 +62,7 @@ def _error_handler(err, func, args):
elif err == errcode('e_invalid_type'):
raise InvalidTypeError(msg)
if err == errcode('e_invalid_id'):
raise IDError(msg)
raise InvalidIDError(msg)
elif err == errcode('e_geometry'):
raise GeometryError(msg)
elif err == errcode('e_data'):

7
openmc/capi/filter.py
View file

@ -144,11 +144,16 @@ class MaterialFilterView(FilterView):
_dll.openmc_material_filter_set_bins(self._index, n, bins)
@classmethod
def new(cls, bins=None):
def new(cls, bins=None, filter_id=None):
# Determine filter ID to assign
if filter_id is None:
filter_id = max(filters) + 1
index = c_int32()
_dll.openmc_extend_filters(1, index, None)
_dll.openmc_filter_set_type(index, b'material')
f = cls(index.value)
f.id = filter_id
if bins is not None:
f.bins = bins
return f

24
openmc/capi/material.py
View file

@ -13,6 +13,9 @@ from .error import _error_handler
__all__ = ['MaterialView', 'materials']
# Material functions
_dll.openmc_extend_materials.argtypes = [c_int32, POINTER(c_int32), POINTER(c_int32)]
_dll.openmc_extend_materials.restype = c_int
_dll.openmc_extend_materials.errcheck = _error_handler
_dll.openmc_get_material_index.argtypes = [c_int32, POINTER(c_int32)]
_dll.openmc_get_material_index.restype = c_int
_dll.openmc_get_material_index.errcheck = _error_handler
@ -35,6 +38,9 @@ _dll.openmc_material_set_densities.argtypes = [
c_int32, c_int, POINTER(c_char_p), POINTER(c_double)]
_dll.openmc_material_set_densities.restype = c_int
_dll.openmc_material_set_densities.errcheck = _error_handler
_dll.openmc_material_set_id.argtypes = [c_int32, c_int32]
_dll.openmc_material_set_id.restype = c_int
_dll.openmc_material_set_id.errcheck = _error_handler
class MaterialView(_View):
@ -67,6 +73,10 @@ class MaterialView(_View):
_dll.openmc_material_get_id(self._index, mat_id)
return mat_id.value
@id.setter
def id(self, mat_id):
_dll.openmc_material_set_id(self._index, mat_id)
@property
def nuclides(self):
return self._get_densities()[0]
@ -100,7 +110,7 @@ class MaterialView(_View):
density_array = as_array(densities, (n.value,))
return nuclide_list, density_array
def add_nuclide(name, density):
def add_nuclide(self, name, density):
"""Add a nuclide to a material.
Parameters
@ -113,6 +123,18 @@ class MaterialView(_View):
"""
_dll.openmc_material_add_nuclide(self._index, name.encode(), density)
@classmethod
def new(cls, material_id=None):
# Determine ID to assign
if material_id is None:
material_id = max(materials) + 1
index = c_int32()
_dll.openmc_extend_materials(1, index, None)
mat = cls(index.value)
mat.id = material_id
return mat
def set_density(self, density):
"""Set density of a material.

10
openmc/capi/tally.py
View file

@ -140,11 +140,17 @@ class TallyView(_View):
_dll.openmc_tally_set_scores(self._index, len(scores), scores_)
@classmethod
def new(cls):
def new(cls, tally_id=None):
# Determine ID to assign
if tally_id is None:
tally_id = max(tallies) + 1
index = c_int32()
_dll.openmc_extend_tallies(1, index, None)
_dll.openmc_tally_set_type(index, b'generic')
return cls(index.value)
tally = cls(index.value)
tally.id = tally_id
return tally
class _TallyMapping(Mapping):

4
src/api.F90
View file

@ -34,10 +34,13 @@ module openmc_api
private
public :: openmc_calculate_volumes
public :: openmc_cell_get_id
public :: openmc_cell_get_fill
public :: openmc_cell_set_fill
public :: openmc_cell_set_temperature
public :: openmc_energy_filter_get_bins
public :: openmc_energy_filter_set_bins
public :: openmc_extend_filters
public :: openmc_extend_materials
public :: openmc_extend_tallies
public :: openmc_filter_get_id
public :: openmc_filter_get_type
@ -59,6 +62,7 @@ module openmc_api
public :: openmc_material_get_densities
public :: openmc_material_set_density
public :: openmc_material_set_densities
public :: openmc_material_set_id
public :: openmc_material_filter_get_bins
public :: openmc_material_filter_set_bins
public :: openmc_mesh_filter_set_mesh

78
src/geometry_header.F90
View file

@ -6,7 +6,7 @@ module geometry_header
use constants, only: HALF, TWO, THREE, INFINITY, K_BOLTZMANN, &
MATERIAL_VOID, NONE
use dict_header, only: DictCharInt, DictIntInt
use material_header, only: Material, materials, material_dict
use material_header, only: Material, materials, material_dict, n_materials
use nuclide_header
use sab_header
use stl_vector, only: VectorReal
@ -447,6 +447,33 @@ contains
end function openmc_get_cell_index
function openmc_cell_get_fill(index, type, indices, n) result(err) bind(C)
integer(C_INT32_T), value, intent(in) :: index
integer(C_INT), intent(out) :: type
integer(C_INT32_T), intent(out) :: n
type(C_PTR), intent(out) :: indices
integer(C_INT) :: err
err = 0
if (index >= 1 .and. index <= size(cells)) then
associate (c => cells(index))
type = c % type
select case (type)
case (FILL_MATERIAL)
n = size(c % material)
indices = C_LOC(c % material(1))
case (FILL_UNIVERSE, FILL_LATTICE)
n = 1
indices = C_LOC(c % fill)
end select
end associate
else
err = E_OUT_OF_BOUNDS
call set_errmsg("Index in cells array is out of bounds.")
end if
end function openmc_cell_get_fill
function openmc_cell_get_id(index, id) result(err) bind(C)
! Return the ID of a cell
integer(C_INT32_T), value :: index
@ -463,9 +490,56 @@ contains
end function openmc_cell_get_id
function openmc_cell_set_fill(index, type, n, indices) result(err) bind(C)
! Set the fill for a fill
integer(C_INT32_T), value, intent(in) :: index ! index in cells
integer(C_INT), value, intent(in) :: type
integer(c_INT32_T), value, intent(in) :: n
integer(C_INT32_T), intent(in) :: indices(n)
integer(C_INT) :: err
integer :: i, j
err = 0
if (index >= 1 .and. index <= size(cells)) then
associate (c => cells(index))
select case (type)
case (FILL_MATERIAL)
if (allocated(c % material)) deallocate(c % material)
allocate(c % material(n))
c % type = FILL_MATERIAL
do i = 1, n
j = indices(i)
if (j == 0) then
c % material(i) = MATERIAL_VOID
else
if (j >= 1 .and. j <= n_materials) then
c % material(i) = j
else
err = E_OUT_OF_BOUNDS
call set_errmsg("Index " // trim(to_str(j)) // " in the &
&materials array is out of bounds.")
end if
end if
end do
case (FILL_UNIVERSE)
c % type = FILL_UNIVERSE
case (FILL_LATTICE)
c % type = FILL_LATTICE
end select
end associate
else
err = E_OUT_OF_BOUNDS
call set_errmsg("Index in cells array is out of bounds.")
end if
end function openmc_cell_set_fill
function openmc_cell_set_temperature(index, T, instance) result(err) bind(C)
! Set the temperature of a cell
integer(C_INT32_T), value, intent(in) :: index ! cell index in cells
integer(C_INT32_T), value, intent(in) :: index ! index in cells
real(C_DOUBLE), value, intent(in) :: T ! temperature
integer(C_INT32_T), optional, intent(in) :: instance ! cell instance

71
src/material_header.F90
View file

@ -14,12 +14,14 @@ module material_header
private
public :: free_memory_material
public :: openmc_extend_materials
public :: openmc_get_material_index
public :: openmc_material_add_nuclide
public :: openmc_material_get_id
public :: openmc_material_get_densities
public :: openmc_material_set_density
public :: openmc_material_set_densities
public :: openmc_material_set_id
!===============================================================================
! MATERIAL describes a material by its constituent nuclides
@ -28,7 +30,7 @@ module material_header
type, public :: Material
integer :: id ! unique identifier
character(len=104) :: name = "" ! User-defined name
integer :: n_nuclides ! number of nuclides
integer :: n_nuclides = 0 ! number of nuclides
integer, allocatable :: nuclide(:) ! index in nuclides array
real(8) :: density ! total atom density in atom/b-cm
real(8), allocatable :: atom_density(:) ! nuclide atom density in atom/b-cm
@ -237,6 +239,37 @@ contains
! C API FUNCTIONS
!===============================================================================
function openmc_extend_materials(n, index_start, index_end) result(err) bind(C)
! Extend the materials array by n elements
integer(C_INT32_T), value, intent(in) :: n
integer(C_INT32_T), optional, intent(out) :: index_start
integer(C_INT32_T), optional, intent(out) :: index_end
integer(C_INT) :: err
type(Material), allocatable :: temp(:) ! temporary materials array
if (n_materials == 0) then
! Allocate materials array
allocate(materials(n))
else
! Allocate materials array with increased size
allocate(temp(n_materials + n))
! Move original materials to temporary array
temp(1:n_materials) = materials(:)
! Move allocation from temporary array
call move_alloc(FROM=temp, TO=materials)
end if
! Return indices in materials array
if (present(index_start)) index_start = n_materials + 1
if (present(index_end)) index_end = n_materials + n
n_materials = n_materials + n
err = 0
end function openmc_extend_materials
function openmc_get_material_index(id, index) result(err) bind(C)
! Returns the index in the materials array of a material with a given ID
integer(C_INT32_T), value :: id
@ -269,6 +302,7 @@ contains
real(8) :: awr
integer, allocatable :: new_nuclide(:)
real(8), allocatable :: new_density(:)
logical, allocatable :: new_p0(:)
character(:), allocatable :: name_
name_ = to_f_string(name)
@ -277,7 +311,7 @@ contains
if (index >= 1 .and. index <= size(materials)) then
associate (m => materials(index))
! Check if nuclide is already in material
do j = 1, size(m % nuclide)
do j = 1, m % n_nuclides
k = m % nuclide(j)
if (nuclides(k) % name == name_) then
awr = nuclides(k) % awr
@ -296,15 +330,20 @@ contains
if (err == 0) then
! Extend arrays
n = size(m % nuclide)
n = m % n_nuclides
allocate(new_nuclide(n + 1))
new_nuclide(1:n) = m % nuclide
if (n > 0) new_nuclide(1:n) = m % nuclide
call move_alloc(FROM=new_nuclide, TO=m % nuclide)
allocate(new_density(n + 1))
new_density(1:n) = m % atom_density
if (n > 0) new_density(1:n) = m % atom_density
call move_alloc(FROM=new_density, TO=m % atom_density)
allocate(new_p0(n + 1))
if (n > 0) new_p0(1:n) = m % p0
new_p0(n + 1) = .false.
call move_alloc(FROM=new_p0, TO=m % p0)
! Append new nuclide/density
k = nuclide_dict % get_key(to_lower(name_))
m % nuclide(n + 1) = k
@ -368,6 +407,23 @@ contains
end function openmc_material_get_id
function openmc_material_set_id(index, id) result(err) bind(C)
! Set the ID of a material
integer(C_INT32_T), value, intent(in) :: index
integer(C_INT32_T), value, intent(in) :: id
integer(C_INT) :: err
if (index >= 1 .and. index <= n_materials) then
materials(index) % id = id
call material_dict % add_key(id, index)
err = 0
else
err = E_OUT_OF_BOUNDS
call set_errmsg("Index in materials array is out of bounds.")
end if
end function openmc_material_set_id
function openmc_material_set_density(index, density) result(err) bind(C)
! Set the total density of a material in atom/b-cm
integer(C_INT32_T), value, intent(in) :: index
@ -396,14 +452,15 @@ contains
integer(C_INT) :: err
integer :: i
integer :: stat
character(C_CHAR), pointer :: string(:)
character(len=:, kind=C_CHAR), allocatable :: name_
if (index >= 1 .and. index <= size(materials)) then
associate (m => materials(index))
! If nuclide/density arrays are not correct size, reallocate
if (n /= size(m % nuclide)) then
deallocate(m % nuclide, m % atom_density, m % p0)
if (n /= m % n_nuclides) then
deallocate(m % nuclide, m % atom_density, m % p0, STAT=stat)
allocate(m % nuclide(n), m % atom_density(n), m % p0(n))
end if