EGSnrc
EGSnrc copied to clipboard
nrc-cnrc
Geometry error in CD geometry with conestack
It appears that one of the changes we've made (it was probably me) to fix other geometry errors has introduced a new one for CD geometries, that arises in this very simple input file. It only seems to occur immediately with a phase-space source, for some reason. I haven't narrowed down the exact commit where this occurs yet, but it does disappear in v2017 and appear in v2018.
Edit: I have narrowed it down, see comment below.
Here is the original reddit thread where this was reported: https://www.reddit.com/r/EGSnrc/comments/kv1b0o/geometry_errors_in_egs_chamber_using_iaea_phase/
##################################################################################
##################################################################################
#
# Chamber in air at 0,0,100
# Varian 6X Truebeam (v2) phase Space source is at 26.7cm
# Also uses EGS collimated source for comparison
#
##################################################################################
##################################################################################
#################################################################
# MC transport parameters
#################################################################
:start MC transport parameter:
Global Ecut = 0.521
Global Pcut = 0.01
:stop MC transport parameter:
#################################################################
# Random number generator seeds
#################################################################
:start rng definition:
initial seeds = 91 2556
:stop rng definition:
###############################################################################
# Media definitions
###############################################################################
:start media definition:
ae = 0.521
ap = 0.01
ue = 50.511
up = 50
:start dryair:
density correction file = air_dry_nearsealevel
:stop dryair:
:start polystyrene:
density correction file = polystyrene
:stop polystyrene:
:start c552:
density correction file = c-552air-equivalentplastic
:stop c552:
:stop media definition:
############################################################################
# Source Definition
############################################################################
:start source definition:
:start source:
name = phsps
library = iaea_phsp_source
iaea phase space file = /home/rwtownson/Downloads/Varian_TrueBeam6MV_02
particle type = all
cutout = -1.33 1.33 -1.33 1.33 #gives 10x10 @100cm
recycle photons = 0
recycle electrons = 0
:stop source:
:start source:
name = photons
library = egs_collimated_source
charge = 0
:start source shape:
type = point
position = 0 0 200
:stop source shape:
:start target shape:
library = egs_rectangle
rectangle = -5 -5 5 5
:stop target shape:
:start spectrum:
type = monoenergetic
energy = 6
:stop spectrum:
:stop source:
:start source: # translate/rotate collimated source to mimic phasespace source position and beam direction
name = trans_photons
library = egs_transformed_source
source name = photons
:start transformation:
translation = 0 0 200
rotation vector = 0 0 -1
:stop transformation:
:stop source:
#simulation source = trans_photons # choose one or the other
simulation source = phsps
:stop source definition:
############################################################################
# Build the geometry
############################################################################
:start geometry definition:
# create chamber tip as sphere #
:start geometry:
name = chamber_tip
library = egs_spheres
midpoint = 0.0 0.0 100.0
radii = 0.05 0.3 0.34 # units cm
:start media input:
media = c552 dryair
set medium = 0 0
set medium = 1 1
set medium = 2 0
:stop media input:
:stop geometry:
# create chamber body as conestack
:start geometry:
name = chamber_body
library = egs_cones
type = EGS_ConeStack
axis = 0.5 0.0 100.0 -1 0 0 # ie x-axis, first module starts at 0.5, note extra 5mm to provide overlap
# unsheathed electrode
:start layer:
thickness = 1.1
top radii = 0.05 0.3 0.34
bottom radii = 0.05 0.3 0.34
media = c552 dryair c552
:stop layer:
# sheathed electrode
:start layer:
thickness = 0.23
top radii = 0.06 0.3 0.34
bottom radii = 0.06 0.3 0.34
media = c552 dryair c552
:stop layer:
# chamber body start
:start layer:
thickness = 0.5
top radii = 0.34
bottom radii = 0.34
media = polystyrene
:stop layer:
# next chamber body piece
:start layer:
thickness = 1.3
top radii = 0.34
bottom radii = 0.34
media = polystyrene
:stop layer:
:stop geometry:
# define chopping planes for final chamber composite geometry
:start geometry:
name = cd_planes
library = egs_planes
type = EGS_Xplanes
positions = -4.0 0.0 4.0
:stop geometry:
# create the composite geometry chamber from tip and body
:start geometry:
name = chamber
library = egs_cdgeometry
base geometry = cd_planes
set geometry = 0 chamber_body
set geometry = 1 chamber_tip
:stop geometry:
# create air column
:start geometry:
name = air_column
library = egs_box
box size = 50 50 400
:start media input:
media = dryair
:stop media input:
:stop geometry:
# inscribe chamber into air column
:start geometry:
name = chamber_in_air
library = egs_genvelope
base geometry = air_column
inscribed geometries = chamber
:stop geometry:
simulation geometry = chamber_in_air
:stop geometry definition:
################################################################################
# Run Control
################################################################################
:start run control:
ncase = 1e6
nbatch = 10
geometry error limit = 20
:stop run control:
#################################################################
# Variance reduction - NB these are egs_chamber specific
#################################################################
:start variance reduction:
# photon splitting = 100 # commented out, not used
cs enhancement = 0 # turned off
# :start range rejection: # commented out, not used
rejection = 256
Esave = 0.521
cavity geometry = air_column
rejection range medium = dryair
:stop range rejection:
:stop variance reduction:
#################################################################
# Scoring options - NB these are egs_chamber specific
#################################################################
:start scoring options:
:start calculation geometry:
geometry name = chamber_in_air
cavity regions = 2 5 14 # chamber cavity regions
cavity mass = 0.1
cavity geometry = chamber # for RR, not used
enhance regions = 2 5 # for xcse, not used
enhancement = -128 # use -ve to set all regions, not used
:stop calculation geometry:
:stop scoring options:
############################################################################
### AUSGAB OBJECTS
############################################################################
:start ausgab object definition:
:start ausgab object:
name = tracks
library = egs_track_scoring
score photons = yes
score electrons = yes
score positrons = yes
start scoring = 0
stop scoring = 5000
buffer size = 20000
:stop ausgab object:
:stop ausgab object definition:
Apparently adding edep_local to struct EGS_Epcont in egs_interface2.h was the cause of this. Commenting out one or more of the EGS_Epcont parameters in egs_interface2.h makes the geometry errors go away! @mainegra @ftessier any idea what's going on here?? It seems that the number of parameters in EGS_Epcont plays some role, but the parameters from the interface match what you find in egsnrc.macros. Is there some step in the interface process that has been missed?
The troublesome commit: https://github.com/nrc-cnrc/EGSnrc/commit/85d4b7d6a6344400a4cb5f6e74df255a9ebda2d6
Does that change not require a modification on the mortran back end as well, so the structs remain aligned? Oh wait, it is there in egsnrc.macros.
Woah, there is an unrelated bug with the order of the variables Epcont->vstep and Epcont->tvstep, compared to the Mortran side! Fortunately, these variables are never used on the egs++, but I will fix that in passing.
I can't figure it out. The structs seems correct (modulo the inversion of vstep and tvstep). So my only possible conclusion at this point in time is that the code is correct, and by making it correct it is revealing an existing bug elsewhere? That's a stretch I know. Other things to consider:
- are you testing from a fresh clone?
- is it necessary to reconfigure the code to get the interface linked on the Mortran side?
- are there compiler issues or flags that cause a misalignment of the structs between g++ and fortran? unlikely.
- are both fortran and egs++ compiled in double presision?
It was not a fresh clone in my case, but it was for the users who found the bug. I don't believe it's necessary to reconfigure, all you need to do is recompile the application for the interface to be updated. I'm using the default compiler flags except with C++11, but the other users were using the defaults only.
I think you could be right that this part of the code is fine and just highlighting a different bug.. I'll have to try going back to yet older versions or (sob) actually stepping through the geometry error to see what is happening.