Message: Re: Recording (n, gamma) reactions, and misc. newbie questions  Not Logged In (login) 
That went a lot easier than I expected. I actually have an
implementation already thanks to your help! But there may be issues.
Just to make sure I'm not confusing concepts, though and thinking I'm measuring things that I'm not, could you let me know if this looks right?
void SteppingAction::UserSteppingAction(const G4Step* step) { if (step>GetTrack()>GetGlobalTime() / microsecond > 1.0) return;
const auto volume = step>GetPreStepPoint()>GetTouchableHandle()>GetVolume()>GetLogicalVolume()>GetName(); ; fEventAction>data[volume].energy += step>GetTotalEnergyDeposit();
if (step>GetTrack()>GetDefinition()>GetParticleName() == "neutron") { ++fEventAction>data[volume].interactions; for (auto iter : *step>GetSecondary()) { const auto def = iter>GetDefinition(); const auto name = def>GetParticleName(); if (name == "neutron") ++fEventAction>data[volume].multiplications; else { if (def>GetAtomicNumber() >= 1) ++fEventAction>data[volume].captures[name]; } } } } So, that is, only look at events that occur in the first microsecond; register the amount of energy deposited in each object in the scene; and for each neutron, 1) increment the number of neutron interactions; 2) accumulate data on each new secondary particle created in this step. For each secondary, if it's a neutron, we increment the number of neutron multiplication events; and for each (n, gamma) reaction, we increment the number of the new resultant element. The scene set up in the detector is a world box containing a beryllium cylinder:
auto world_shape = new G4Box("_world", 5.1*m, 5.1*m, 5.1*m); auto world = new G4LogicalVolume(world_shape, air, "world", 0, 0, 0); auto world_p = new G4PVPlacement( 0, G4ThreeVector(), world, "world_p", 0, false, 0);
auto multiplier_shape = new G4Tubs("_multiplier", 0*cm, 200*cm, 2*m, 0*deg, 360*deg); auto multiplier = new G4LogicalVolume(multiplier_shape, beryllium, "multiplier", 0,0,0); new G4PVPlacement( new G4RotationMatrix(90*deg, 0*deg, 0*deg), G4ThreeVector(4*m,0*m,0*m), multiplier, "multiplier_p", world, false, 0, true); And the cylinder is bombarded by neutrons with a mean energy of 35 MeV: /gps/particle neutron /gps/time 1 ns /gps/pos/type Plane /gps/pos/shape Circle /gps/pos/centre 5 0 0 m /gps/pos/halfx 5 cm /gps/pos/halfy 5 cm /gps/ang/type beam1d /gps/ang/rot1 0 0 1 /gps/ang/rot2 0 1 0 /gps/ang/sigma_r 0. deg /gps/ene/type Gauss /gps/ene/ezero 35 MeV /gps/ene/sigma 8 MeV /gps/ene/min 1 eV /gps/ene/max 2000 MeV /run/beamOn 10000 I get the following results:
 multiplier Energy: 49700.1 MeV Interactions: 363865 Multiplications: 10861 Captures: Be9: 1818971 He6: 1087 Li6: 5 Li7: 153 Li8: 28 N14: 303 O16: 87 alpha: 9034 deuteron: 32 triton: 160
world Energy: 616.553 MeV Interactions: 13573 Multiplications: 43 Captures: Ar40: 4 B10: 1 B11: 38 Be9: 6693 C12: 6 C13: 10 C14: 3 Li7: 2 N14: 774 N15: 5 O16: 220 O17: 1 alpha: 49 deuteron: 6 proton: 6 triton: 1  This looked mostly plausible in general. I assumed that the Be9 secondaries created in the beryllium and to a lesser extend the world were kicked off of the cylinder by the impact energy. I was surprised to see no Be10, but thought maybe the cross section is just too small versus other possible reactions. There's some random captures in the air... the neutrons and alphas are similar in number, which would make sense if both of the 2n are treated as new particles. So in general, it looked mostly plausible. But what got me was that the yield looked rather low, I'd expect more (n, 2n) reactions. LIFE for example has a useful multiplication factor of 1.8, with only 17MeV neutrons. So I simplified it: /gps/ene/type Gauss /gps/ene/ezero 17 MeV /gps/ene/sigma 0 MeV /run/beamOn 10000 Which gave me:
 multiplier Energy: 8446.55 MeV Interactions: 280122 Multiplications: 2 Captures: Be10: 2 Be9: 1056507 N14: 160 O16: 179
world Energy: 15.4674 MeV Interactions: 13767 Multiplications: 1 Captures: Be9: 9166 C14: 1 N14: 102 N15: 2 O16: 68 proton: 1  Only two multiplications? It's like the (n,2n) reaction doesn't work with monoenergetic 17MeV neutrons, despite the fact that it has a rather large cross section ... And even when I put the GPS in the middle of a 25 meter block of beryllium in a 50meter world, and nothing affects "world" anymore, the energy deposited in the beryllium never exceeds 10033.1 MeV:
multiplier Energy: 10033.1 MeV Interactions: 431651 Multiplications: 1 Transmutation: Be10: 5 Be9: 1726106 ... despite the fact that the GPS should be shooting out 10000 17MeV neutrons, aka 170000 MeV total, not counting nuclear reactions. So there's something I'm not getting here... One thing I forgot to ask last time: I figured out that I can set a solid's temperature by initializing it like one initializes a gas, but with kStateSolid: auto beryllium = new G4Material("Beryllium", 4., 9.012182*g/mole, 1.8480*g/cm3, kStateSolid, 1300*kelvin, 1*atmosphere);
But is it possible to do anything more complex than that  for example, temperature gradients? And does Geant4 take into account heat flow at all? There's not much in the examples I see dealing with heat, which surprises me. ? I hope I haven't been bothering you with these newbie questions; I'm just trying to get the hang of things here.
 kv, Karen

Inline Depth:  Outline Depth:  Add message: 
to: 