張耀鋒，尹 遠，張春雷，蔣文麗，蘇 俊

(北京師范大學 核科學與技術學院，北京 100875)

時間投影電離室(TPC)是典型的氣體探測器，依據電離信號收集與放大部件的不同，TPC主要可分為絲狀電極TPC探測器和氣體電子倍增膜(GEM)類型的TPC探測器[1-3]。在TPC探測器設計及后續(xù)實驗階段，對探測器進行模擬研究，在加快設計進度、獲得優(yōu)化實驗參數、提高實驗效率等方面具有重要意義。GARFIELD程序是由歐洲核子研究中心(CERN)開發(fā)的適用于氣體探測器的模擬程序，為當前國內外主流氣體探測器模擬程序之一[4]。該程序主要采用有限邊界法，實現對具有細絲電極的氣體探測器進行電場計算以及后續(xù)的模擬優(yōu)化。對于細絲電極的電場計算，GARFIELD程序僅適用于具有無限長細絲結構的二維探測器模型。此外，GARFIELD程序也可讀入外部有限元程序計算得到的電場數據，以進行后續(xù)的探測器模擬計算。讀入外部程序計算數據時往往受到許多限制，如對于ANSYS有限元電場計算數據，二維模型只能采用plane 121元素單元，三維模型只能采用solid 123元素單元。在GARFIELD程序讀入有限元電場數據后，還需進行插值計算，以便進行數據顯示與后續(xù)粒子漂移模擬。現有GARFIELD程序中插值部分存在數據完備性檢驗過于嚴格致使插值計算失敗、電場數據經插值后連續(xù)性不好、部分位置無法輸出有效電場數據等問題，導致程序對探測器的模擬無法進行。

有限元程序ANSYS常用于包含細絲電極結構的探測器的電場模擬。采用solid 122元素單元以及掃掠方式進行探測器主體結構模型構建，可有效避免細絲結構由于大的長寬比(長度與直徑比)帶來的因元素數目過多導致的計算失敗的問題。鑒于GARFIELD模擬程序不能直接讀入由solid 122元素單元計算獲得的電場數據，本文擬對GARFIELD程序進行修改編譯，并引入Quadratic Shepard插值法，以獲得連續(xù)性很好的電場分布數據，用于后續(xù)的探測器模擬。

1 細絲結構探測器的3D電場計算

常見的TPC探測器采用數百根甚至上千根不同直徑的細絲電極作為探測器的信號收集裝置，可實現電離信號的雪崩放大以及收集測量[5]。細絲電極的主要特點是長度大、直徑小，具有大的長寬比。典型的電極細絲長約1 m、直徑為25 μm，其長寬比大于104。大的長寬比給細絲電極的3D電場計算帶來了巨大挑戰(zhàn)。以有限元程序ANSYS為例，用于靜電場計算的3D元素單元主要有solid 122和solid 123兩類[6]，如圖1所示?？梢妔olid 123元素單元為基本的三菱錐體類型，solid 122元素單元為類四棱柱類型。

圖1 ANSYS用于3D電場計算的元素單元結構Fig.1 Structure of element type for 3D electric field calculation in ANSYS

圖2 用于絲狀電極的兩種元素單元網格劃分方法Fig.2 Two meshing ways of element type for wire electrode

這兩種元素單元常通過不同的網格劃分方法得到。如圖2a所示，通過自由網格劃分方法可得到solid 123元素單元，而通過掃掠式網格劃分方法則得到solid 122元素單元。在細絲結構的網格劃分過程中，由于三棱錐體表面方向的不規(guī)則性，致使細絲結構整體劃分網格的元素基本尺寸較小，約為細絲直徑的1/10。因此，若采用solid 123元素單元進行細絲結構的網格模型構建，會造成整體元素數目過多。對于包含數千根細絲的探測器模型，總的元素數目更會急劇增加，導致網格劃分失敗或電場計算失敗。而采用掃掠方式進行細絲結構劃分，細絲縱向solid 122元素單元的基本尺寸與截面位置的元素基本尺寸無關，可大幅減小細絲縱向上的元素網格密度，從而將探測器整體模型的元素數目控制在可接受的范圍，以便實現探測器的整體3D電場計算。

由于探測器包含數目較多的細絲電極，常規(guī)的有限元電場計算方法在模型構建、網格劃分、邊界條件施加等過程中需要大量重復性的微區(qū)操作，極大地降低了整體電場的計算效率。利用有限元程序的腳本代碼如ANSYS程序的APDL語言，可實現模擬程序的自動運行，從而提高模擬計算效率。但腳本語言的運行需預先對各電極絲進行準確的幾何邊界定位并進行特定單元的模擬運行核驗，并在此基礎上對細絲電極進行逐一建模以及網格劃分。因此腳本語言的引入并不能提高太多的計算效率，只是方便了具有相似配置結構的探測器模型的重復計算。為進一步提高TPC多絲探測器的3D電場計算效率，采用節(jié)點組合法進行電場計算更為有效。該方法首先構建不包含任何細絲單元的宏觀探測器結構模型，并進行網格劃分。然后確定包含所有細絲電極的微小區(qū)域，對區(qū)域內的網格進行多次加密，直至元素節(jié)點分布能完整表征細絲截面形狀。最后依據細絲電極的位置和半徑信息，采用APDL腳本語言自動選取邊界對應位置的節(jié)點，進行邊界條件施加。對單絲模型采用節(jié)點組合法處理并進行邊界條件施加后模型截面節(jié)點分布結果和采用常規(guī)掃掠式網格劃分及邊界條件施加后截面節(jié)點分布結果示于圖3。由圖3可知，在適當進行節(jié)點加密后，可得到與常規(guī)方法幾乎完全一致的細絲單元的截面節(jié)點分布。這種逐步加密的方法，使得節(jié)點分布沿徑向向外存在一定規(guī)律的節(jié)點密度衰減，這種不同的節(jié)點密度分布與細絲電極徑向電場分布的指數衰減規(guī)律基本保持一致，從而有助于使用有限節(jié)點數獲得較準確的計算結果。

圖3 節(jié)點組合法(a)和常規(guī)掃掠式網格劃分法(b)在邊界施加后的細絲截面節(jié)點分布Fig.3 Distribution of node for wire section after boundary condition applied by node-construction method (a) and regular sweeping method (b)

2 Quadratic Shepard插值法優(yōu)化應用

2.1 有限元數據讀入

在采用solid 122元素單元對包含細絲結構的TPC探測器進行3D電場計算后，可獲得節(jié)點位置信息以及電勢分布數據。為使GARFIELD程序能讀入電場數據，對GARFIELD程序進行了修改編譯，以便能使用該程序通過電場計算函數得到任意位置處電場矢量及電勢的值，具體過程如下：1) 在程序即將進行電場計算時，將程序代碼中的ICTYPE變量直接賦值為0，使得程序能轉向處理有限元數據的代碼；2) 針對處理有限元數據的代碼，增加能讀入有限元電場計算數據的相應代碼，并進行插值計算，最終得到任意位置的電場和電勢數據，返回指定函數。所有新增加代碼采用Fortran語言編寫，并與GARFIELD程序進行交叉編譯。編譯完成后，程序可正常運行。

2.2 Quadratic Shepard插值法應用

GARFIELD程序需給出任意位置的電場數據，以進行探測器模擬。由于ANSYS有限元程序計算得到的結果為所有節(jié)點處的數據，因此在數據讀入GARFIELD程序后需采用插值法進行數據處理。GARFIELD程序原始的插值方法為雙線性插值法，選取當前位置所處的體元素的所有節(jié)點為插值數據來源。如對于solid 123元素單元，則插值數據最多選取10個節(jié)點。如此小數目的插值數據源會給插值結果帶來影響。如數據源點空間分布較均勻，則插值結果較為準確。反之，則計算結果誤差較大。而ANSYS對于細絲電極電場的計算采用節(jié)點加密方法實現，因此節(jié)點數目密度分布非常不均勻。如，在接近細絲表面處，節(jié)點數目密度較大；而遠離細絲表面的區(qū)域，節(jié)點數目密度較小。這兩處區(qū)域內的節(jié)點密度相差數十倍到上百倍。相應的電場計算數據如采用GARFIELD程序原始的插值代碼會導致完備性檢驗失敗而無法輸出結果或輸出的結果連續(xù)性差。因此，針對細絲電極的電場分布特點，需引入新的插值方法，以獲得連續(xù)性較好的電場分布結果。

圖4 solid 122元素單元對應的插值處理方法Fig.4 Interpolation method for data from solid 122 element type

Quadratic Shepard插值法采用距離倒數加權的最小二乘法實現插值計算[7]，該方法一般不限定插值所需的數據源點的數目，可使用多個參數進行插值的圓滑性調節(jié)，最終獲得連續(xù)性分布較好的插值結果。在實際插值應用過程中，常選取插值點臨近區(qū)域的數十至上百個節(jié)點進行插值。依據電場計算時采用的元素單元類型的不同，分別采用不同的插值數據點選取方式。對于solid 122元素單元構建的模型，由于在縱向上采用掃掠方式劃分網格，在縱向上分布的節(jié)點數密度相對于截面方向小得多?？刹捎脠D4所示的處理方法，插值點位于(x0,y0,z0)處，先在與插值節(jié)點相鄰的兩個橫向平面z=z1和z=z2上進行二維插值，得到相應的兩組插值結果f(x0,y0,z1)和f(x0,y0,z2)，然后在縱向上進行一維線性插值，則最終結果為：

f(x0,y0,z0)=f(x0,y0,z1)+

(1)

對于solid 123元素單元構建的模型，可直接采用Quadratic Shepard的3D插值法，選取所插值節(jié)點附近適當數目的原始數據點進行插值計算，最終得到準確的插值結果。

2.3 插值結果分析

圖5 電勢沿徑向分布的插值計算結果與理論計算結果Fig.5 Interpolation and theory result of potential distribution along radial direction

為檢驗Quadratic Shepard插值法的準確性，選取單絲電極結構模型進行電場有限元計算及插值結果分析。細絲電極直徑為20 μm，長度為60 cm，施加電勢為1 000 V。在距中心10 cm處同心圓柱外表面設為0 V電勢區(qū)域，該區(qū)域縱向(z方向)長度為100 cm。對該模型采用節(jié)點組合方法進行3D電場有限元計算，之后將數據讀入GARFIELD程序實現Quadratic Shepard插值計算以及電場數據輸出。在模型中心平面處，電勢沿徑向分布的插值計算結果與理論計算結果的對比示于圖5，可見兩者符合很好，最大相對偏差在0.5%以內。此外，還對該平面處電場強度沿徑向分布進行分析比較，以便對電勢分布的連續(xù)性進行評價，結果示于圖6。由圖6可見，由插值計算得到的電場結果與理論計算結果符合較好，兩者之間的相對偏差較圖5有所增大，最大相對偏差在3%以內。其中偏差較大的位置均對應連續(xù)加密的分界區(qū)域，在該區(qū)域內由于節(jié)點加密導致節(jié)點空間分布呈現嚴重偏心表象。由于Quadratic Shepard插值法在計算過程中需選取計算目標點周圍一定數目的節(jié)點數據，而這些選取的節(jié)點數據分布的空間不均勻性是影響插值結果準確程度的重要因素。上述插值計算過程中，實際選取其周圍距離最近的120個數據點進行插值。對于某一數據點，采用最小距離法所選取的原始節(jié)點空間分布如圖7a所示，數據分布的偏心程度異常嚴重。為改善插值結果，對最小距離法進行修正，加入修正因子，如式(2)所示：

圖6 電場強度沿徑向分布的插值計算結果與理論計算結果(a)及比較(b)Fig.6 Interpolation and theory calculation result of electric field strength radial distribution (a) and comparison (b)

十字標線中心為所需計算的數據位置圖7 最小距離法(a)和最小距離+修正因子法(b)所得節(jié)點分布Fig.7 Node distribution for minimum distance method (a) and modified minimum distance method (b)

(2)

式中：(x0,y0)為所需要插值計算的數據點坐標；ri為每個節(jié)點與數據點(x0,y0)之間的距離；mi為每個節(jié)點的修正因子。mi的計算方法為：選取距離該節(jié)點最近的周圍n個節(jié)點，取其中的最大距離。將最終獲得的fi最小的120個數據點作為插值點陣列。其中n是優(yōu)化的參數。如果n=30，則圖7a的節(jié)點選取結果改變?yōu)閳D7b。此時插值計算結果的相對偏差由2.9%減小為0.17%。在實際插值計算中，可針對不同位置的節(jié)點選取不同的優(yōu)化參數n，來獲得較準確的插值結果。

采用上述電場數據插值計算方法，對包含多絲電極結構的TPC探測器進行3D電場的有限元計算。在中心平面處分別進行GARFIELD程序原始插值方法的電場計算輸出和Quadratic Shepard插值法的電場計算輸出，結果示于圖8。TPC探測器門電極細絲位于y=-1.4 cm處，相鄰細絲x方向間距0.1 cm。細絲電極電勢施加按照計數分為奇、偶兩組，分別施加Vg+V0和Vg-V0電勢。其中Vg=-115 V、V0=30 V。Vg為門電極細絲處于測量打開狀態(tài)下施加的基礎電勢，而±V0為探測器向測量關閉狀態(tài)下轉換時在電極細絲上按照奇偶兩組分別施加的偏壓電勢。兩種插值方法所得結果的差別主要表現在細絲電極臨近區(qū)域。GARFIELD程序原始插值方法計算數據分布連續(xù)性較差，部分區(qū)域發(fā)生錯誤輸出，會導致后續(xù)模擬程序中斷；而Quadratic Shepard插值法對于電極臨近區(qū)域能獲得連續(xù)性較好的數據分布，有助于獲得準確的探測器模擬結果。圖8結果表明，Quadratic Shepard插值法是一種適用于細絲電極結構電場有限元程序計算的有效的數據后處理方法。

圖8 TPC探測器門電極細絲附近區(qū)域等電勢分布輸出結果Fig.8 Output of equal potential distribution near gating wire of TPC detector

3 結論

利用有限元程序ANSYS并結合節(jié)點組合方法實現了包含細絲電極結構的TPC探測器的3D電場計算。在GARFIELD模擬程序讀入電場數據后，采用Quadratic Shepard插值法對電場數據進行插值計算，在周圍插值點分布極端不均勻情況下獲得了較為準確的結果。Quadratic Shepard插值法有助于實現GARFIELD程序對包含細絲結構的TPC探測器的3維尺度的準確模擬。