FissionFragment/src/FFDetectorConstruction.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // -------------------------------------------------------------
0027 //  =============== Begin Documentation Comments ===============
0028 //!
0029 //! \file       FFDetectorConstruction.cc
0030 //! \author     B. Wendt (brycen.linn.wendt@cern.ch)
0031 //! \date       June 06, 2014
0032 //!
0033 //! \brief      Implementation of the FFDetectorConstruction class
0034 //!
0035 //! \details    The model simulated is based off a subcritical assembly design
0036 //!             with 20% enriched meat
0037 //!
0038 //  ================ End Documentation Comments ================
0039 //
0040 //  Modified:
0041 //
0042 //  23-06-14                                              BWendt
0043 //  Fixed issue with the automatic placement of the meat not working. Solution
0044 //      was to use the correct units "inch" in the y-direction as well.
0045 //  Coincidentally eliminated the need to using the 'std::abs()" from the
0046 //      "cmath" library.
0047 //  Implemented method "PlaceFuelPlates"
0048 //
0049 // -------------------------------------------------------------
0050
0051 #include "FFDetectorConstruction.hh"
0052
0053 #include "G4Box.hh"
0054 #include "G4Element.hh"
0055 #include "G4Isotope.hh"
0056 #include "G4LogicalVolume.hh"
0057 #include "G4NistManager.hh"
0058 #include "G4PVPlacement.hh"
0059 #include "G4SystemOfUnits.hh"
0060 #include "G4Tubs.hh"
0061 #include "G4VPhysicalVolume.hh"
0062 #include "globals.hh"
0063
0064 static const G4double inch = 2.54 * cm;
0065
0066 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0067 FFDetectorConstruction::FFDetectorConstruction() : G4VUserDetectorConstruction()
0068 {
0069   DefineMaterials();
0070 }
0071
0072 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0073 G4VPhysicalVolume* FFDetectorConstruction::Construct()
0074 {
0075   G4ThreeVector position;
0076 #ifdef NDEBUG
0077   G4bool const overlapChecking = false;
0078 #else
0079   G4bool const overlapChecking = true;
0080 #endif  // NDEBUG
0081
0082   //
0083   // Create the world
0084   //
0085   const G4double worldSize = 40.0 * inch;
0086   G4Box* const solidWorld = new G4Box("World",  // the name
0087                                       worldSize,  // x size
0088                                       worldSize,  // y size
0089                                       worldSize);  // z size
0090   G4LogicalVolume* const logicalWorld = new G4LogicalVolume(solidWorld,  // the solid volume
0091                                                             fAir,  // the material
0092                                                             solidWorld->GetName());  // the name
0093   // Center at the origin
0094   position.set(0.0, 0.0, 0.0);
0095   G4VPhysicalVolume* const physicalWorld =
0096     new G4PVPlacement(NULL,  // no rotation
0097                       position,  // must be at origin
0098                       logicalWorld,  // the logical volume
0099                       logicalWorld->GetName(),  // the name
0100                       NULL,  // no mother volume
0101                       false,  // no boolean ops
0102                       0,  // copy number
0103                       overlapChecking);  // check for overlaps
0104
0105   //
0106   // Create the graphite pile that the subcritical assembly rests on.
0107   //
0108   const G4double floorH = 30.0 * inch;
0109   const G4ThreeVector floorPosition(0.0, 0.0, 0.0);
0110   G4Box* const solidFloor = new G4Box("Floor",  // the name
0111                                       worldSize,  // x size
0112                                       worldSize,  // y size
0113                                       floorH * 0.5);  // z size
0114   G4LogicalVolume* const logicalFloor = new G4LogicalVolume(solidFloor,  // the solid volume
0115                                                             fGraphite,  // the material
0116                                                             solidFloor->GetName());  // the name
0117   // Shift down so the top is at the origin
0118   position.set(0.0, 0.0, -floorH * 0.5);
0119   new G4PVPlacement(NULL,  // no rotation
0120                     position,  // position
0121                     logicalFloor,  // the logical volume
0122                     logicalFloor->GetName(),  // the name
0123                     logicalWorld,  // the mother volume
0124                     false,  // no boolean ops
0125                     0,  // copy number
0126                     overlapChecking);  // check for overlaps
0127
0128   //
0129   // Create the tank
0130   //
0131   const G4double tankWallThickness = 0.25 * inch;
0132   const G4double tankOR = 18.0 * inch;
0133   const G4double tankH = 39.0 * inch;
0134   G4Tubs* const solidTank = new G4Tubs("Tank_Wall",  // the name
0135                                        0.0,  // inner radius
0136                                        tankOR,  // outer radius
0137                                        tankH * 0.5,  // half height
0138                                        0.0 * deg,  // start angle
0139                                        360.0 * deg);  // end angle
0140   G4LogicalVolume* const logicalTank = new G4LogicalVolume(solidTank,  // the solid volume
0141                                                            fAluminum,  // the material
0142                                                            solidTank->GetName());  // the name
0143   // Shift up so the base is at the origin
0144   position.set(0.0, 0.0, tankH * 0.5);
0145   new G4PVPlacement(NULL,  // no rotation
0146                     position,  // shift up
0147                     logicalTank,  // the logical volume
0148                     logicalTank->GetName(),  // the name
0149                     logicalWorld,  // the mother volume
0150                     false,  // no boolean ops
0151                     0,  // copy number
0152                     overlapChecking);  // check for overlaps
0153   // Top 3 inches are air
0154   const G4double tankAirH = 3.0 * inch;
0155   G4Tubs* const solidTankAir = new G4Tubs("Tank_Air",  // the name
0156                                           0.0,  // inner radius
0157                                           tankOR - tankWallThickness,  // outer radius
0158                                           tankAirH * 0.5,  // half height
0159                                           0.0 * deg,  // start angle
0160                                           360.0 * deg);  // end angle
0161   G4LogicalVolume* const logicalTankAir = new G4LogicalVolume(solidTankAir,  // the solid volume
0162                                                               fAir,  // the material
0163                                                               solidTankAir->GetName());  // the name
0164   // Shift up so that the top of the air is the same as the top of the tank
0165   position.set(0.0, 0.0, (tankH - tankAirH) * 0.5);
0166   new G4PVPlacement(NULL,  // no rotation
0167                     position,  // shift ip
0168                     logicalTankAir,  // the logical volume
0169                     logicalTankAir->GetName(),  // the name
0170                     logicalTank,  // the mother volume
0171                     false,  // no boolean ops
0172                     0,  // copy number
0173                     overlapChecking);  // check for overlaps
0174   // Fill remaining area with water
0175   const G4double tankH2OH = (tankH - (tankAirH + tankWallThickness));
0176   G4Tubs* const solidTankH2O = new G4Tubs("Tank_H2O",  // the name
0177                                           0.0,  // inner radius
0178                                           tankOR - tankWallThickness,  // outer radius
0179                                           tankH2OH * 0.5,  // half height
0180                                           0.0 * deg,  // start angle
0181                                           360.0 * deg);  // end angle
0182   G4LogicalVolume* const logicalTankH2O = new G4LogicalVolume(solidTankH2O,  // the solid volume
0183                                                               fAluminum,  // the material
0184                                                               solidTankH2O->GetName());  // the name
0185   // Shift up so that the top of the water is at the bottom of the air
0186   const G4double centerOfH2O = (tankH - tankH2OH) * 0.5 - tankAirH;
0187   position.set(0.0, 0.0, centerOfH2O);
0188   new G4PVPlacement(NULL,  // no rotation
0189                     position,  // shift to origin
0190                     logicalTankH2O,  // the logical volume
0191                     logicalTankH2O->GetName(),  // the name
0192                     logicalTank,  // the mother volume
0193                     false,  // no boolean ops
0194                     0,  // copy number
0195                     overlapChecking);  // check for overlaps
0196
0197   //
0198   // Fuel plates
0199   //
0200   const G4double plateX = 3.0 * inch;
0201   const G4double plateY = 0.08 * inch;
0202   const G4double plateZ = 26.0 * inch;
0203   const G4double meatX = 2.75 * inch;
0204   const G4double meatY = 0.04 * inch;
0205   const G4double meatZ = 24.0 * inch;
0206   const G4double xSpacing = 5.0 * inch;
0207   const G4double ySpacing = 0.3 * inch;
0208   const G4double plateRadius = 12.0 * inch;
0209   // Define the aluminim claddiing
0210   G4Box* const solidPlate = new G4Box("Plate_Cladding",  // the name
0211                                       plateX * 0.5,  // x size
0212                                       plateY * 0.5,  // y size
0213                                       plateZ * 0.5);  // z size
0214   G4LogicalVolume* const logicalPlate = new G4LogicalVolume(solidPlate,  // the solid volume
0215                                                             fAluminum,  // the material
0216                                                             solidPlate->GetName());  // the name
0217   // Place the meat inside the cladding
0218   G4Box* const solidMeat = new G4Box("Plate_Meat",  // the name
0219                                      meatX * 0.5,  // x size
0220                                      meatY * 0.5,  // y size
0221                                      meatZ * 0.5);  // z size
0222   G4LogicalVolume* const logicalMeat = new G4LogicalVolume(solidMeat,  // the solid volume
0223                                                            fUO2_20E,  // the material
0224                                                            solidMeat->GetName());  // the name
0225   // The meat goes into the exact center of the plate
0226   position.set(0.0, 0.0, 0.0);
0227   new G4PVPlacement(NULL,  // no rotation
0228                     position,  // position
0229                     logicalMeat,  // the logical volume
0230                     logicalMeat->GetName(),  // the name
0231                     logicalPlate,  // the mother volume
0232                     false,  // no boolean ops
0233                     0,  // copy number
0234                     overlapChecking);  // check for overlaps
0235   // The plate will be centered in the z-direction within the water
0236   // Simulate a subcritical assembly loading within a radius of 12 inches
0237   bool placeMe;
0238
0239   position.setZ(0.0);
0240   fCopyNumber = 0;
0241   for (double x = 0.0; x <= plateRadius; x += xSpacing) {
0242     // 5 rows of plates
0243     for (double y = 0.0; y <= plateRadius; y += ySpacing) {
0244       placeMe = false;
0245
0246       // Fuel plate must be completely within the radius to be placed
0247       if (std::sqrt(x * x + y * y) < plateRadius) {
0248         // Leave a 1 inch radius opening in the middle for the neutron
0249         // source
0250         if (std::sqrt(x * x + y * y) > 1.0 * inch) {
0251           placeMe = true;
0252         }
0253       }
0254
0255       if (placeMe) {
0256         PlaceFuelPlate(x, y, logicalPlate, logicalTankH2O);
0257         PlaceFuelPlate(x, -y, logicalPlate, logicalTankH2O);
0258         if (x > 0.0) {
0259           PlaceFuelPlate(-x, y, logicalPlate, logicalTankH2O);
0260           PlaceFuelPlate(-x, -y, logicalPlate, logicalTankH2O);
0261         }
0262       }
0263     }
0264   }
0265   G4cout << fCopyNumber << " plates were added to the subcritical assembly" << G4endl;
0266
0267   //
0268   // Neutron Source
0269   //
0270   // TODO create the AmBe material in DefineMaterials() and use it here
0271   //      For now steel is used, but the logical volume is used in the
0272   //      PrimaryGeneratorAction to know where to emit the neutrons from
0273   const G4double sourceH = 2 * inch;
0274   const G4double sourceR = 0.2 * inch;
0275   G4Tubs* const solidSource = new G4Tubs("NeutronSource",  // the name
0276                                          0.0,  // inner radius
0277                                          sourceR,  // outer radius
0278                                          sourceH * 0.5,  // half height
0279                                          0.0 * deg,  // start angle
0280                                          360.0 * deg);  // end angle
0281   G4LogicalVolume* const logicalSource = new G4LogicalVolume(solidSource,  // the solid volume
0282                                                              fStainlessSteel,  // the material
0283                                                              solidSource->GetName());  // the name
0284   // Place in the exact center of the water tank
0285   position.set(0.0, 0.0, 0.0);
0286   new G4PVPlacement(NULL,  // no rotation
0287                     position,  // shift to origin
0288                     logicalSource,  // the logical volume
0289                     logicalSource->GetName(),  // the name
0290                     logicalTankH2O,  // the mother volume
0291                     false,  // no boolean ops
0292                     0,  // copy number
0293                     overlapChecking);  // check for overlaps
0294
0295   //
0296   // Detector Tower
0297   //
0298   const G4double polyS = 3.0 * inch;
0299   const G4double polyH = 18.0 * inch;
0300   G4Box* const solidPoly = new G4Box("Poly",  // the name
0301                                      polyS,  // x size
0302                                      polyS,  // y size
0303                                      polyH);  // z size
0304   G4LogicalVolume* const logicalPoly = new G4LogicalVolume(solidPoly,  // the solid volume
0305                                                            fPolyethylene,  // the material
0306                                                            solidPoly->GetName());  // the name
0307   // The polyethylene detector tower goes just outside the tank at 45 deg
0308   G4double radiusToPolyCenter = (tankOR / std::sqrt(2.0)) + std::sqrt(2.0) * polyS;
0309   position.set(-radiusToPolyCenter, radiusToPolyCenter, polyH);
0310   new G4PVPlacement(NULL,  // no rotation
0311                     position,  // position
0312                     logicalPoly,  // the logical volume
0313                     logicalPoly->GetName(),  // the name
0314                     logicalWorld,  // the mother volume
0315                     false,  // no boolean ops
0316                     0,  // copy number
0317                     overlapChecking);  // check for overlaps
0318   // Create the detector shell
0319   G4double shellR = 0.3 * inch;
0320   G4double shellH = 6.5 * inch;
0321   G4Tubs* const solidShell = new G4Tubs("Detector_Shell",  // the name
0322                                         0.0,  // inner radius
0323                                         shellR,  // outer radius
0324                                         shellH * 0.5,  // half height
0325                                         0.0 * deg,  // start angle
0326                                         360.0 * deg);  // end angle
0327   G4LogicalVolume* const logicalShell = new G4LogicalVolume(solidShell,  // the solid volume
0328                                                             fStainlessSteel,  // the material
0329                                                             solidShell->GetName());  // the name
0330   // Place in the exact center of the polyethylene tower
0331   position.set(0.0, 0.0, 0.0);
0332   new G4PVPlacement(NULL,  // no rotation
0333                     position,  // shift to origin
0334                     logicalShell,  // the logical volume
0335                     logicalShell->GetName(),  // the name
0336                     logicalPoly,  // the mother volume
0337                     false,  // no boolean ops
0338                     0,  // copy number
0339                     overlapChecking);  // check for overlaps
0340   // Create the BF3 detector
0341   G4double BF3R = 0.2 * inch;
0342   G4double BF3H = 6.0 * inch;
0343   G4Tubs* const solidBF3 = new G4Tubs("Detector_BF3_Core",  // the name
0344                                       0.0,  // inner radius
0345                                       BF3R,  // outer radius
0346                                       BF3H * 0.5,  // half height
0347                                       0.0 * deg,  // start angle
0348                                       360.0 * deg);  // end angle
0349   G4LogicalVolume* const logicalBF3 = new G4LogicalVolume(solidBF3,  // the solid volume
0350                                                           fBF3_96E,  // the material
0351                                                           solidBF3->GetName());  // the name
0352   // Place in the exact center of the shell
0353   position.set(0.0, 0.0, 0.0);
0354   new G4PVPlacement(NULL,  // no rotation
0355                     position,  // shift to origin
0356                     logicalBF3,  // the logical volume
0357                     logicalBF3->GetName(),  // the name
0358                     logicalShell,  // the mother volume
0359                     false,  // no boolean ops
0360                     0,  // copy number
0361                     overlapChecking);  // check for overlaps
0362
0363   return physicalWorld;
0364 }
0365
0366 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0367 void FFDetectorConstruction::DefineMaterials(void)
0368 {
0369   static G4NistManager* const nist = G4NistManager::Instance();
0370
0371   fAir = nist->FindOrBuildMaterial("G4_AIR");
0372   fAluminum = nist->FindOrBuildMaterial("G4_Al");
0373   fGraphite = nist->FindOrBuildMaterial("G4_GRAPHITE");
0374   fPolyethylene = nist->FindOrBuildMaterial("G4_POLYETHYLENE");
0375   fStainlessSteel = nist->FindOrBuildMaterial("G4_STAINLESS-STEEL");
0376   fWater = nist->FindOrBuildMaterial("G4_WATER");
0377
0378   /*// List available materials
0379   std::vector< G4String > materials = nist->GetNistMaterialNames();
0380   for(unsigned int i = 0;
0381       i < materials.size();
0382       ++i)
0383   {
0384       G4cout << materials[i] << G4endl;
0385   }*/
0386
0387   //
0388   // Define the 20% enriched UO2
0389   //
0390   // First we need to start by creating the isotopes
0391   G4double const U235Enrichment = 0.2;
0392   G4double const U238Enrichment = 0.8;
0393   G4Isotope* const iU235 = new G4Isotope("iU235",  // name
0394                                          92,  // ZZZ
0395                                          235,  // AAA
0396                                          235.053930 * (g / mole));  // molecular weight
0397   G4Isotope* const iU238 = new G4Isotope("iU238",  // name
0398                                          92,  // ZZZ
0399                                          238,  // AAA
0400                                          238.050788 * (g / mole));  // molecular weight
0401   // Now create the elements and add the isotopes
0402   G4Element* const U235 = new G4Element("U235",  // name
0403                                         "U235",  // symbol
0404                                         1);  // number of isotopes
0405   U235->AddIsotope(iU235,  // isotope
0406                    1.0);  // abundance
0407   G4Element* const U238 = new G4Element("U238",  // name
0408                                         "U238",  // symbol
0409                                         1);  // number of isotopes
0410   U238->AddIsotope(iU238,  // isotope
0411                    1.0);  // abundance
0412   G4Element* const oxygen = nist->FindOrBuildElement("O");
0413   // Calculate the mass fractions
0414   const G4double UO2MolecularWeight =
0415     U235->GetA() * U235Enrichment + U238->GetA() * U238Enrichment + oxygen->GetA() * 2;
0416   const G4double U235MassFraction = (U235->GetA() * U235Enrichment) / UO2MolecularWeight;
0417   const G4double U238MassFraction = (U238->GetA() * U238Enrichment) / UO2MolecularWeight;
0418   const G4double oxygenMassFraction = (oxygen->GetA() * 2) / UO2MolecularWeight;
0419   // create the material and add the elements
0420   fUO2_20E = new G4Material("UO2_20E",  // name
0421                             10.97 * (g / cm3),  // density
0422                             3);  // number of components
0423   fUO2_20E->AddElement(U235,  // element
0424                        U235MassFraction);  // mass fraction
0425   fUO2_20E->AddElement(U238,  // element
0426                        U238MassFraction);  // mass fraction
0427   fUO2_20E->AddElement(oxygen,  // element
0428                        oxygenMassFraction);  // mass fraction
0429
0430   //
0431   // Define the BF3
0432   //
0433   // The BF3 is B-10 enriched
0434   // http://www.orau.org/ptp/collection/proportional%20counters/bf3info.htm
0435   G4double const B10Enrichment = 0.96;
0436   G4double const B11Enrichment = 0.04;
0437   G4Isotope* const iB10 = new G4Isotope("iB10",  // name
0438                                         5,  // ZZZ
0439                                         10,  // AAA
0440                                         10.0129370 * (g / mole));  // molecular weight
0441   G4Isotope* const iB11 = new G4Isotope("iB11",  // name
0442                                         5,  // ZZZ
0443                                         11,  // AAA
0444                                         11.0093054 * (g / mole));  // molecular weight
0445   // Now create the elements and add the isotopes
0446   G4Element* const B10 = new G4Element("B10",  // name
0447                                        "B10",  // symbol
0448                                        1);  // number of isotopes
0449   B10->AddIsotope(iB10,  // isotope
0450                   1.0);  // abundance
0451   G4Element* const B11 = new G4Element("B11",  // name
0452                                        "B11",  // symbol
0453                                        1);  // number of isotopes
0454   B11->AddIsotope(iB11,  // isotope
0455                   1.0);  // abundance
0456   G4Element* const flouride = nist->FindOrBuildElement("F");
0457   // Calculate the mass fractions
0458   const G4double BF3MolecularWeight =
0459     B10->GetA() * B10Enrichment + B11->GetA() * B11Enrichment + flouride->GetA() * 3;
0460   const G4double B10MassFraction = (B10->GetA() * B10Enrichment) / BF3MolecularWeight;
0461   const G4double B11MassFraction = (B11->GetA() * B11Enrichment) / BF3MolecularWeight;
0462   const G4double flourideMassFraction = (flouride->GetA() * 3) / BF3MolecularWeight;
0463   // create the material and add the elements
0464   fBF3_96E = new G4Material("BF3_96E",  // name
0465                             2.5 * (kg / m3),  // density
0466                             3);  // number of components
0467   fBF3_96E->AddElement(B10,  // element
0468                        B10MassFraction);  // mass fraction
0469   fBF3_96E->AddElement(B11,  // element
0470                        B11MassFraction);  // mass fraction
0471   fBF3_96E->AddElement(flouride,  // element
0472                        flourideMassFraction);  // mass fraction
0473 }
0474
0475 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0476 void FFDetectorConstruction::PlaceFuelPlate(double x, double y,
0477                                             G4LogicalVolume* const myLogicalVolume,
0478                                             G4LogicalVolume* const parentLogicalVolume)
0479 {
0480   G4ThreeVector position(x, y);
0481   std::ostringstream copyName;
0482   copyName << "Plate@Location    X:" << std::setprecision(2) << x / inch << "    Y:" << y / inch;
0483
0484   new G4PVPlacement(NULL,  // no rotation
0485                     position,  // position
0486                     myLogicalVolume,  // the logical volume
0487                     copyName.str(),  // the name
0488                     parentLogicalVolume,  // the mother volume
0489                     false,  // no boolean ops
0490                     fCopyNumber++,  // copy number
0491                     true);  // check for overlaps
0492 }
0493
0494 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0495 FFDetectorConstruction::~FFDetectorConstruction()
0496 {
0497   // Nothing here
0498 }