diff --git a/README.md b/README.md index ed60d72..5bd3802 100644 --- a/README.md +++ b/README.md @@ -137,7 +137,7 @@ $$

2. 质子 - 与电子保持一致，额外考虑幂律平移模型$j=j_0(E-\beta)^{-\alpha}$，同电子进行`BIC test`，最终选用幂律平移模型$j_0=559.76377,\ \alpha=2.40873,\ \beta=-1.53424$。 + 与电子保持一致，幂律模型修改为$j=j_0(E-\beta)^{-\alpha}$，进行`BIC test`，最终选用幂律模型$j_0=559.76377,\ \alpha=2.40873,\ \beta=-1.53424$。
@@ -151,19 +151,19 @@ $$

## 辐射源设置 -使用`Geant4`的`General Particle Source`（`GPS`）生成源，考虑在`World`球面上均匀随机发射粒子（各向同性），保证入射方向为球内，则球面上的总通量为： +使用`Geant4`的`General Particle Source`（`GPS`）生成源，考虑在`World`球面上均匀随机发射粒子（各向同性），保证入射方向为球内，则球面上的计数率为： $$ F = A\int_{E_a}^{E_b} \phi(E)\mathrm{d}E,\ A=4\pi R^2 $$ -`GCR`的单位为$\mathrm{m^{-2}sr^{-1}s^{-1}}$，乘以$4\pi A$即可得到通量；俘获辐射的单位为$\mathrm{cm^{-2}s^{-1}}$，乘以$A$即可得到通量，最终结果如下表所述。 +`GCR`的单位为$\mathrm{m^{-2}sr^{-1}s^{-1}}$，乘以$4\pi A$即可得到计数率；俘获辐射的单位为$\mathrm{cm^{-2}s^{-1}}$，乘以$A$即可得到计数率。经过检验，模型中使用的外壳（3层，2mm 铝 + 10mm 芳纶 + 5mm 铝）能够抵挡能量在$50\ \mathrm{MeV}$以下的质子和能量在$100\ \mathrm{keV}$以下的电子，最终结果如下表所述。
| Partical | Model | $E_{min}(\mathrm{MeV})$ | $E_{max}(\mathrm{MeV})$ | Flux($\mathrm{s^{-1}}$) | |:--------:|:--------------------------------------------------------------------:|:-----------------------:|:-----------------------:|:-----------------------:| -| $e^{-}$ | $j=j_0e^{-\frac{E}{E_0}}$ | 0.04 | 10 | 1.91227E+12 | -| `p` | $j=j_0(E-\beta)^{-\alpha}$ | 0.1 | 1000 | 5623595363 | +| $e^{-}$ | $j=j_0e^{-\frac{E}{E_0}}$ | 0.1 | 10 | 9.23242E+11 | +| `p` | $j=j_0(E-\beta)^{-\alpha}$ | 50 | 1000 | 42853994.15 | | `H` | $\frac{C_1\beta^{\alpha-1}}{pc}\left(\frac{pc}{pc+C_2}\right)^{C_3}$ | 200 | 100000 | 12854310.92 | | `He` | $C_1e^{-C_2E}(1-e^{-C_3E+C_4})$ | 200 | 100000 | 1718828.12 | | `Li` | $C_1e^{-C_2E}(1-e^{-C_3E+C_5})$ | 200 | 100000 | 12982.36284 | @@ -195,6 +195,41 @@ $$
+考虑到`GCR`和俘获辐射计数率相差太大，分为两部分分别进行模拟，同时按比例进行放缩。 + +### GCR +电子计数率与质子计数率之比约为`21544`，设置质子的源强为`1`，模拟`923285`个入射事件（数值上为$F\times10^{-6}\ \mathrm{s}$），因此最终的剂量需要乘以$10^6\ \mathrm{s^{-1}}\times30\ \mathrm{d}\times86400\ \mathrm{s/d}$倍。 +```macro +/gps/source/add 21544 +/gps/particle e- +/gps/pos/type Surface +/gps/pos/shape Sphere +/gps/pos/radius 15 m +/gps/ene/type Exp +/gps/ene/min 100 keV +/gps/ene/max 10 MeV +/gps/ene/ezero 0.0824 +/gps/ang/type iso +/gps/ang/maxtheta 90 deg + +/gps/source/add 1 +/gps/particle proton +/gps/pos/type Surface +/gps/pos/shape Sphere +/gps/pos/radius 15 m +/gps/ene/type Pow +/gps/ene/min 200 MeV +/gps/ene/max 100 GeV +/gps/ene/alpha 2.40873 +/gps/ang/type iso +/gps/ang/maxtheta 90 deg + +/run/beamOn 923285 +``` + +### 俘获辐射 +假定`Co`的计数率为1，设定其他粒子的计数率，模拟`707013`个入射事件，因此最终的剂量需要乘以$20.9\ \mathrm{s^{-1}}\times30\ \mathrm{d}\times86400\ \mathrm{s/d}$倍。 + ## 空间站建模 1. 尺寸与分区[^6]
diff --git a/debug.mac b/debug.mac new file mode 100644 index 0000000..6ade0ff --- /dev/null +++ b/debug.mac @@ -0,0 +1,28 @@ +# 多线程设置 +/run/numberOfThreads 1 + +# verbose +/control/verbose 1 +/run/verbose 0 +/event/verbose 0 +/tracking/verbose 0 +/gps/verbose 0 + +# 初始化 +/run/initialize + +# 电子射程 +/gps/particle e- +/gps/energy 100 keV +/gps/position 0 0 5 m +/gps/ang/type iso +/gps/ang/maxphi 60 deg +/gps/ang/maxtheta 20 deg + +# 质子射程 +/gps/particle proton +/gps/energy 50 MeV +/gps/position 0 0 5 m +/gps/ang/type iso +/gps/ang/maxphi 60 deg +/gps/ang/maxtheta 20 deg diff --git a/electron.mac b/electron.mac index 78ba6b9..6a4376b 100644 --- a/electron.mac +++ b/electron.mac @@ -13,7 +13,7 @@ /gps/source/clear # Trapped electron -/gps/source/add 1 +/gps/source/add 21544 /gps/particle e- /gps/pos/type Surface /gps/pos/shape Sphere @@ -25,5 +25,19 @@ /gps/ang/type iso /gps/ang/maxtheta 90 deg +/gps/source/add 1 +/gps/particle proton +/gps/pos/type Surface +/gps/pos/shape Sphere +/gps/pos/radius 15 m +/gps/ene/type Pow +/gps/ene/min 200 MeV +/gps/ene/max 100 GeV +/gps/ene/alpha 2.40873 +/gps/ang/type iso +/gps/ang/maxtheta 90 deg + +/run/beamOn 923285 + # list -/gps/source/list \ No newline at end of file +/gps/source/list diff --git a/src/DetectorConstruction.cpp b/src/DetectorConstruction.cpp index f662299..deedaa8 100644 --- a/src/DetectorConstruction.cpp +++ b/src/DetectorConstruction.cpp @@ -229,6 +229,19 @@ G4VPhysicalVolume* DetectorConstruction::Construct() { G4LogicalVolume* logicWorld = new G4LogicalVolume(solidWorld, G4Material::GetMaterial("Vacuum"), "World"); G4VPhysicalVolume* physicsWorld = new G4PVPlacement(0, G4ThreeVector(), logicWorld, "World", 0, false, 0, true); + // Debug for particle range + // G4Box* solidPlane1 = new G4Box("Plane1", 10 * m, 10 * m, 2 * mm); + // G4LogicalVolume* logicPlane1 = new G4LogicalVolume(solidPlane1, G4Material::GetMaterial("AluminumAlloySeries5"), "Plane1"); + // new G4PVPlacement(0, G4ThreeVector(), logicPlane1, "Plane1", logicWorld, false, 0, true); + + // G4Box* solidPlane2 = new G4Box("Plane2", 10 * m, 10 * m, 10 * mm); + // G4LogicalVolume* logicPlane2 = new G4LogicalVolume(solidPlane2, G4Material::GetMaterial("Taparan"), "Plane2"); + // new G4PVPlacement(0, G4ThreeVector(0, 0, -6), logicPlane2, "Plane2", logicWorld, false, 0, true); + + // G4Box* solidPlane3 = new G4Box("Plane3", 10 * m, 10 * m, 5 * mm); + // G4LogicalVolume* logicPlane3 = new G4LogicalVolume(solidPlane3, G4Material::GetMaterial("AluminumAlloySeries5"), "Plane3"); + // new G4PVPlacement(0, G4ThreeVector(0, 0, -14.5), logicPlane3, "Plane3", logicWorld, false, 0, true); + // 节点舱 ConstructSectionSphere(logicWorld, (5.18 + 0.33) * m); // 过渡段 1