Báo Cáo Nghiên Cứu: Phương Pháp Hybrid Monte Carlo Và Ứng Dụng Trong Việc Xác Định Hiệu Suất Của Detector NaI(Tl)

Tóm tắt nghiên cứu

Trong kỹ thuật phổ kế gamma và phân tích hoạt độ phóng xạ, việc xác định chính xác hiệu suất tổng (Total Efficiency) và hiệu suất nội (Intrinsic Efficiency) của đầu dò nhấp nháy $\text{NaI(Tl)}$ là yêu cầu tiên quyết. Nghiên cứu do ThS. Hoàng Đức Tâm (Khoa Vật lý, Trường Đại học Sư phạm TP.HCM) chủ nhiệm, phối hợp cùng Phòng Thí nghiệm Vật lý Hạt nhân (Trường Đại học Khoa học Tự nhiên - ĐHQG TP.HCM) thực hiện, tập trung giải quyết bài toán: Làm thế nào để tính toán nhanh chóng, chính xác hiệu suất của detector $\text{NaI(Tl)}\ 3'' \times 3''$ cho các cấu hình nguồn điểm (Point Source) và nguồn đĩa (Disk Source) mà không phải chịu chi phí tính toán nặng nề của phương pháp Monte Carlo truyền thống?

Nghiên cứu ứng dụng phương pháp Hybrid Monte Carlo (HMC) – sự kết hợp giữa kỹ thuật lấy mẫu ngẫu nhiên Monte Carlo để xác định hướng bay photon và phương pháp giải tích để tính chính xác quãng đường chuyển dời trong tinh thể. Nhóm tác giả đã phát triển chương trình tính toán CalcTotEff trên nền tảng Fortran và xây dựng phần mềm giao diện Visual Basic.

Kết quả mô phỏng được đối chứng với dữ liệu thực nghiệm đo từ nguồn chuẩn $^{137}\text{Cs}$ ($661.66\text{ keV}$) trên hệ phổ kế đa kênh 8192 kênh. Độ sai lệch giữa mô phỏng và thực nghiệm đạt mức dưới $10%$ ở hầu hết các khoảng cách. Đặc biệt, nghiên cứu đã làm sáng tỏ quy luật phụ thuộc phi tuyến của hiệu suất nội theo tỉ số hình học khoảng cách – bán kính detector ($d/R_d$), mở ra giải pháp tối ưu hóa cấu hình đo đạc trong nghiên cứu hạt nhân và công nghiệp.

Bối cảnh và tầm quan trọng

Thực trạng ứng dụng detector $\text{NaI(Tl)}$

Kể từ khi được Robert Hofstadter phát minh vào năm 1948, đầu dò nhấp nháy $\text{NaI(Tl)}$ (Sodium Iodide pha tạp Thallium) đã trở thành công cụ ghi nhận bức xạ gamma phổ biến nhất thế giới. Với ưu thế vượt trội về hiệu suất ghi nhận cao, khả năng vận hành ổn định ở nhiệt độ phòng và chi phí chế tạo hợp lý, $\text{NaI(Tl)}$ được ứng dụng rộng rãi từ quan trắc môi trường, y học hạt nhân đến kiểm tra không phá hủy (NDT) trong công nghiệp như dò tìm tắc nghẽn và khuyết tật đường ống dẫn dầu khí.

       [Nguồn phóng xạ] (Điểm / Đĩa)
              │
              │  (Góc phát ngẫu nhiên θ, φ) ──► [Lấy mẫu Monte Carlo]
              ▼
   ┌───────────────────────┐
   │  Detector NaI(Tl)     │
   │  (Kích thước 3" x 3") │ ──► [Tính quãng đường giải tích Λ]
   └───────────────────────┘
              │
              ▼
   [Hệ số suy giảm μ(E)] ──► S(E) = 1 - exp(-μΛ) ──► [Hiệu suất tổng & nội]

Khoảng trống nghiên cứu và tính cấp thiết

Để định lượng hoạt độ phóng xạ tuyệt đối của mẫu, bắt buộc phải biết hàm truyền hiệu suất của đầu dò. Tuy nhiên, hai cách tiếp cận truyền thống đều bộc lộ hạn chế:

  1. Phương pháp giải tích thuần túy: Chỉ giải được các bài toán hình học đối xứng cực kỳ đơn giản; khi áp dụng cho nguồn diện tích (nguồn đĩa) hoặc hình học khối phức tạp, các tích phân góc khối trở nên cực kỳ nan giải.
  2. Phương pháp Monte Carlo đầy đủ (Full MC như MCNP, GEANT4): Mô phỏng từng lịch sử tương tác của hạt, đòi hỏi tài nguyên phần cứng lớn, thời gian tính toán kéo dài và quy trình thiết lập thông số phức tạp.

Phương pháp Hybrid Monte Carlo (HMC) do S. Yalcin đề xuất năm 2007 xuất hiện như một hướng tiếp cận đột phá: dung hòa tính linh hoạt của Monte Carlo và tốc độ vượt trội của giải tích. Tại Việt Nam, việc làm chủ thuật toán, phát triển phần mềm tính toán độc lập và kiểm chuẩn thực nghiệm trên detector kích thước tiêu chuẩn $3'' \times 3''$ mang tính thời sự cao, giúp giảm thiểu sự phụ thuộc vào các gói phần mềm thương mại đắt đỏ.

Methodology và Phương pháp tiếp cận

1. Mô hình toán học Hybrid Monte Carlo

Nguyên lý của phương pháp HMC phân tách quá trình vận chuyển photon thành hai giai đoạn độc lập:

  • Giai đoạn ngẫu nhiên (Monte Carlo): Vị trí phát photon trên bề mặt nguồn đĩa bán kính $R_s$ được chọn ngẫu nhiên qua biến đổi $r_i = R_s \sqrt{\xi}$ (với $\xi$ là số ngẫu nhiên phân bố đều $[0, 1]$ được tạo bởi hàm RANDOM_NUMBER() trong Fortran). Hướng phát photon được đặc trưng bởi góc cực $\theta \in [0, \pi/2]$ và góc phương vị $\phi \in [0, 2\pi]$.
  • Giai đoạn giải tích: Dựa vào góc phát $\theta$ và hình học tinh thể trụ bán kính $R_d$, chiều cao $H$, thuật toán xác định quỹ đạo photon đi vào tinh thể và thoát ra ở mặt đáy hay mặt bên. Từ đó, quãng đường chuyển dời $\Lambda(\theta, \phi)$ được tính chính xác bằng hình học giải tích:

$$\Lambda = \begin{cases} \dfrac{H}{\cos\theta} & \text{khi photon thoát ra từ mặt đáy} \ \dfrac{R_d - r_i \cos\phi}{\sin\theta} - \dfrac{d}{\cos\theta} & \text{khi photon thoát ra từ mặt bên} \end{cases}$$

Xác suất hấp thụ photon mang năng lượng $E$ trong tinh thể được tính qua định luật suy giảm:

$$S(E) = 1 - e^{-\mu(E) \cdot \Lambda}$$

Trong đó, hệ số suy giảm tuyến tính toàn phần $\mu(E)$ được trích xuất từ cơ sở dữ liệu chuẩn quốc tế XCOM (NIST).

+-------------------------------------------------------------------------+
|                    LƯU ĐỒ THUẬT TOÁN HYBRID MONTE CARLO                 |
+-------------------------------------------------------------------------+
| 1. Khởi tạo thông số: Năng lượng E, hình học detector (Rd, H), nguồn (d)|
| 2. Vòng lặp N hạt (N ≈ 2.10^6 - 2.10^8 photon):                         |
|    a. Sinh số ngẫu nhiên chọn tọa độ phát r_i, góc cực θ, góc phương vị φ|
|    b. Kiểm tra điều kiện hình học: Photon có đi vào detector không?     |
|       - Không -> Bỏ qua.                                                |
|       - Có    -> Xác định thoát đáy hay thoát bên -> Tính quãng đường Λ |
|    c. Tính xác suất tương tác: S_i(E) = 1 - exp(-μ(E) * Λ)               |
| 3. Tính hiệu suất tổng: ε_t = (1 / 2N) * Σ S_i(E) (nguồn điểm)          |
| 4. Tính hiệu suất nội: ε_int = ε_t / Hình_học_khối                     |
+-------------------------------------------------------------------------+

2. Thiết kế chương trình và Bố trí thực nghiệm

  • Phần mềm: Chương trình CalcTotEff viết bằng ngôn ngữ Fortran 95/2003 chạy trên nền Plato IDE (Silverfrost), sau đó được đóng gói thành phần mềm ứng dụng giao diện trực quan bằng Visual Basic.
  • Thực nghiệm kiểm chuẩn:
    • Đầu dò: Tinh thể nhấp nháy $\text{NaI(Tl)}$ kích thước $3'' \times 3''$ (đường kính $7.62\text{ cm}$, chiều cao $7.62\text{ cm}$).
    • Nguồn chuẩn: Nguồn điểm $^{137}\text{Cs}$ có hoạt độ danh định $1\ \mu\text{Ci}$, phát tia gamma đơn năng $661.66\text{ keV}$ (xác suất phát $I_\gamma = 85.1%$).
    • Hệ ghi đo: Phổ kế đa kênh MCA 8192 kênh điều khiển bởi phần mềm ADMCA.
    • Xử lý số liệu: Phổ thực nghiệm được trừ phông môi trường, loại bỏ đỉnh tia X tán xạ $32\text{ keV}$ của $^{137m}\text{Ba}$, tính diện tích đỉnh toàn phần (Net Peak Area) và hiệu chỉnh sai số theo định luật lan truyền sai số thống kê.

Phát hiện chính của nghiên cứu

      Hiệu suất nội (ε_int)
         ▲
         │ *                                                    * * * (Bão hòa)
         │  *                                                * *
         │   *                                            * *
         │    *                                        * *
         │     *                                    * *
         │      *                                * *
         │        * *                        * *
         │            * * * *        * * * *
         │                    * * * (Cực tiểu)
         └───────────────────────┼───────────────────────────────► Tỉ số d/Rd
                                d/Rd ≈ 0.5 - 1.0

1. Độ chính xác vượt trội của mã nguồn CalcTotEff

Kết quả tính toán hiệu suất tổng bằng chương trình HMC được so sánh đối chiếu trực tiếp với các công bố quốc tế kinh điển (Nakamura, Belluscio, Vegors, Selim, Yalcin):

Năng lượng (keV) Khoảng cách $d$ HMC (CalcTotEff) Nakamura et al. Selim et al. Sai khác (%)
661.66 $0.001\text{ cm}$ (tiếp xúc) 0.4612 0.4608 0.4610 $< 0.1%$
661.66 $10.0\text{ cm}$ 0.0140 0.0139 0.0141 $< 0.7%$
1332.0 $10.0\text{ cm}$ 0.0118 0.0117 0.0118 $< 0.5%$

Sự trùng khớp gần như tuyệt đối giữa kết quả của nhóm nghiên cứu với các tác giả tiền bối đã khẳng định tính đúng đắn và độ tin cậy của thuật toán.

2. Sự tương thích giữa mô phỏng và đo đạc thực nghiệm $^{137}\text{Cs}$

Khi kiểm chứng với thực nghiệm nguồn chuẩn $^{137}\text{Cs}$, hiệu suất tổng thu được từ phần mềm phản ánh rất sát kết quả đo đạc:

  • Tại các khoảng cách $d \ge 10\text{cm}$ ($10\text{cm}, 15\text{cm}, 20\text{cm}, 30\text{cm}$), sai số tương đối giữa mô phỏng HMC và thực nghiệm đều dưới $4%$, nằm trọn vẹn trong dải sai số thống kê của phép đo.
  • Hiện tượng bất thường tại cự ly cực gần ($d = 5\text{cm}$): Độ chênh lệch giữa lý thuyết và thực nghiệm tăng lên khoảng $9.2%$. Nguyên nhân vật lý được nhóm nghiên cứu chỉ ra là do ở cự ly quá gần, hiệu ứng trùng phùng tổng (True Coincidence Summing) và hiện tượng tán xạ ngược từ vỏ bọc detector trở nên rõ rệt, trong khi mô hình HMC lý thuyết đã lược bỏ sự suy giảm qua vỏ nhôm và không khí.

3. Quy luật phụ thuộc phi tuyến của hiệu suất nội vào tỉ số $d/R_d$

Một trong những phát hiện có giá trị khoa học sâu sắc nhất của đề tài là việc mô tả chi tiết đường cong biến thiên của hiệu suất nội $\epsilon_{int}$ theo tỉ số $d/R_d$ (khoảng cách nguồn chia cho bán kính detector):

  1. Vùng cự ly cực gần ($d/R_d < 0.01$): Hiệu suất nội đạt giá trị rất cao do hầu hết các photon đi vào mặt trước tinh thể với quãng đường xuyên thấu lớn.
  2. Vùng chuyển tiếp giảm dần: Khi tăng khoảng cách từ $d/R_d = 0$ đến $1.0$, hiệu suất nội giảm mạnh và đạt giá trị cực tiểu tại vùng $d/R_d \approx 0.5 - 1.0$. Cơ chế: góc mở photon nghiêng khiến nhiều hạt đâm xuyên qua góc rìa và thoát ra từ mặt bên trước khi kịp tương tác hoàn toàn.
  3. Vùng phục hồi và bão hòa ($d/R_d > 10$): Khi nguồn đặt cách xa ($d/R_d \ge 10$), chùm tia photon tới detector gần như song song và đồng trục. Quãng đường photon đi trong tinh thể đạt cực đại xấp xỉ chiều cao $H$, dẫn tới hiệu suất nội tăng trở lại và bão hòa hoàn toàn khi $d/R_d \ge 100$.

Kết quả này hoàn toàn thống nhất với các quan sát của Jehouani (2000) và Ogundare (2008), mang lại cơ sở lý thuyết định hình vị trí đặt mẫu tối ưu.

4. Quy luật suy giảm theo năng lượng

Đối với cả nguồn điểm và nguồn đĩa ($R_s \le R_d$), hiệu suất nội giảm đơn điệu khi năng lượng bức xạ gamma tăng từ $100\text{ keV}$ lên $2000\text{ keV}$. Điều này hoàn toàn phù hợp với tiết diện tương tác hiệu ứng quang điện và tán xạ Compton của vật liệu $\text{NaI}$.

Đóng góp khoa học và Giá trị thực tiễn

                  CÁC TRỤ CỘT ĐÓNG GÓP CỦA ĐỀ TÀI
  ┌─────────────────────────────────────────────────────────────┐
  │ 1. ĐÓNG GÓP PHƯƠNG PHÁP & HỆ THỐNG MÃ NGUỒN                │
  │    - Tối ưu hóa thuật toán Hybrid Monte Carlo               │
  │    - Xây dựng mã nguồn Fortran (CalcTotEff) & GUI VB        │
  ├─────────────────────────────────────────────────────────────┤
  │ 2. ĐÓNG GÓP LÝ THUYẾT & VẬT LÝ THỰC NGHIỆM                  │
  │    - Khảo sát biến thiên phi tuyến ε_int theo d/Rd          │
  │    - Xác lập ngưỡng bão hòa d/Rd > 10 và cực tiểu d/Rd ≈ 1  │
  ├─────────────────────────────────────────────────────────────┤
  │ 3. ỨNG DỤNG THỰC TIỄN & GIẢNG DẠY                           │
  │    - Tối ưu hóa cự ly đo phổ gamma trong phòng thí nghiệm   │
  │    - Ứng dụng trong kỹ thuật NDT, kiểm tra khuyết tật ống   │
  │    - Phục vụ đào tạo chuyên ngành Vật lý Hạt nhân           │
  └─────────────────────────────────────────────────────────────┘

Đóng góp về mặt phương pháp luận

  • Rút ngắn thời gian tính toán: Phương pháp HMC giúp giảm thời gian chạy mô phỏng từ hàng giờ (ở các mã hạt nhân lớn) xuống chỉ còn vài giây đến vài chục giây trên máy tính cá nhân thông thường mà vẫn đảm bảo độ chính xác tương đương.
  • Tự chủ công nghệ phần mềm: Việc phát triển thành công công cụ giao diện Visual Basic giúp các kỹ thuật viên và sinh viên có thể thao tác nhập thông số hình học, năng lượng và nhận kết quả tức thì mà không cần cài đặt các môi trường mô phỏng phức tạp.

Đóng góp thực tiễn và Khuyến nghị quy trình đo

Từ phát hiện về tỉ số $d/R_d$, nghiên cứu đưa ra khuyến nghị thực nghiệm quan trọng:

  • Để tối đa hóa độ chính xác: Khi tiến hành ghi nhận phổ kế gamma với đầu dò $\text{NaI(Tl)}\ 3'' \times 3''$, nên bố trí nguồn ở vị trí rất gần ($d/R_d < 0.01$) để lấy số đếm cao hoặc ở vị trí đủ xa ($d/R_d > 10$) để hiệu suất nội ổn định ở vùng bão hòa, hạn chế sai số do dịch chuyển vị trí nguồn.
  • Tránh đặt nguồn trong vùng $0.5 \le d/R_d \le 1.0$ vì đây là vùng hiệu suất nhạy cảm nhất với sự thay đổi vị trí, dễ gây sai số lớn trong phân tích hoạt độ.

Đối tượng quan tâm và Lợi ích ứng dụng

Nghiên cứu mang lại giá trị thiết thực cho nhiều nhóm đối tượng:

  • Nhà nghiên cứu & Kỹ thuật viên Vật lý hạt nhân: Sở hữu một công cụ tính toán hiệu suất chuẩn hóa, tin cậy để hiệu chỉnh đầu dò và phân tích phổ gamma từ các mẫu môi trường, địa chất, sinh học.
  • Kỹ sư kiểm tra không phá hủy (NDT): Ứng dụng mô hình truyền qua của chùm tia gamma và độ nhạy của detector $\text{NaI(Tl)}$ để nâng cao độ chính xác khi định vị khuyết tật đường ống, đo bề dày kim loại hoặc giám sát mức chất lỏng trong bồn chứa công nghiệp.
  • Cơ quan quản lý & An toàn bức xạ: Nắm bắt phương pháp đánh giá nhanh hoạt độ các nguồn phóng xạ vô chủ hoặc sự cố bức xạ dựa trên dữ liệu phổ thực địa.
  • Giảng viên và Sinh viên chuyên ngành Vật lý/Kỹ thuật hạt nhân: Tài liệu tham khảo mẫu mực về việc ứng dụng kết hợp giữa mô hình giải tích, kỹ thuật mô phỏng Monte Carlo và đối chứng thực nghiệm.

Các câu hỏi thường gặp (FAQ)

1. Phương pháp Hybrid Monte Carlo vượt trội hơn Monte Carlo thuần túy ở điểm nào?

Hybrid Monte Carlo chỉ sử dụng kỹ thuật ngẫu nhiên để chọn góc phát và điểm phát ban đầu, sau đó áp dụng hình học giải tích để tính trực tiếp quãng đường photon đi xuyên tinh thể. Nhờ đó, phương pháp loại bỏ việc mô phỏng từng bước tán xạ vi mô lặp đi lặp lại, giúp tốc độ tính toán nhanh hơn hàng trăm lần nhưng vẫn giữ độ chính xác cao đối với bài toán xác định hiệu suất tổng.

2. Tại sao hiệu suất nội lại đạt cực tiểu khi tỉ số $d/R_d \approx 0.5 - 1.0$?

Tại dải khoảng cách này ($d$ xấp xỉ bán kính tinh thể $R_d$), góc mở chùm tia từ nguồn tới mép detector tạo ra một tỷ lệ lớn các photon đi chéo xuyên qua các góc mép của tinh thể. Các photon này có quãng đường chuyển dời $\Lambda$ ngắn và dễ thoát ra ngoài qua mặt bên trước khi tương tác truyền năng lượng, khiến xác suất ghi nhận nội tại giảm xuống mức thấp nhất.

3. Nguyên nhân gây sai số gần $10%$ tại khoảng cách đo $5\text{cm}$ là gì?

Ở cự ly $5\text{cm}$, góc khối ghi nhận rất lớn, dẫn tới hai hiện tượng vật lý chưa được tính đến trong mô hình HMC lý thuyết: hiệu ứng trùng phùng tổng (True Coincidence Summing) khi nhiều photon đến detector cùng lúc và sự suy giảm/tán xạ qua nắp vỏ nhôm bảo vệ tinh thể.

4. Thuật toán của đề tài có áp dụng được cho các loại detector khác như HPGe hay LaBr3 không?

Hoàn toàn có thể. Thuật toán cốt lõi của HMC chỉ yêu cầu thay đổi thông số hình học (kích thước bán kính, chiều cao) và cập nhật bộ hệ số suy giảm tuyến tính toàn phần $\mu(E)$ tương ứng với vật liệu bán dẫn HPGe hoặc tinh thể nhấp nháy $\text{LaBr}_3\text{(Ce)}$.

5. Ứng dụng thực tế rõ nét nhất của phần mềm phát triển trong đề tài là gì?

Phần mềm cho phép các phòng thí nghiệm xác định ngay hiệu suất detector ứng với nguồn điểm hoặc nguồn đĩa ở bất kỳ khoảng cách nào mà không bắt buộc phải mua sắm đầy đủ bộ nguồn chuẩn đắt tiền để hiệu chuẩn thực nghiệm cho từng vị trí.

Kết luận

Đề tài nghiên cứu khoa học cấp cơ sở mã số CS.60 của ThS. Hoàng Đức Tâm đã chứng minh tính ưu việt, độ chính xác và tính thực tiễn cao của phương pháp Hybrid Monte Carlo trong việc xác định hiệu suất tổng và hiệu suất nội của detector nhấp nháy $\text{NaI(Tl)}\ 3'' \times 3''$. Việc hoàn thiện mã nguồn CalcTotEff cùng giao diện trực quan Visual Basic và xác lập quy luật phụ thuộc hình học $d/R_d$ là đóng góp ý nghĩa cho cộng đồng nghiên cứu vật lý bức xạ tại Việt Nam.

Trong tương lai, hướng nghiên cứu có thể mở rộng mô hình HMC để tính toán thêm thành phần tán xạ nhiều lần, hiệu ứng tán xạ trong vật liệu bao bọc (vỏ nhôm, lớp phản xạ $\text{MgO}$, đệm quang học) và mở rộng tính toán cho các hình học nguồn phức tạp như nguồn trụ thể tích hay mẫu môi trường dạng hộp Marinelli.

[!TIP] Khám phá thêm: Bạn có thể tham khảo các công cụ tự động hóa tính toán hoặc lập lịch đo lường thực nghiệm bằng cách sử dụng lệnh /schedule hoặc lập kế hoạch chi tiết cho các nghiên cứu mô phỏng bức xạ tiếp theo với /plan.