周琳，陳憲，黃江濤，*，鄧俊，高正紅

1.中國空氣動(dòng)力研究與發(fā)展中心，綿陽 621000

2.西北工業(yè)大學(xué) 航空學(xué)院，西安 710072

雷達(dá)隱身技術(shù)是提高目標(biāo)生存能力的技術(shù)，低可探測特征是未來軍用飛機(jī)應(yīng)具備的重要特征之一［1-2］。隨著反隱身技術(shù)的發(fā)展，飛行器的氣動(dòng)/隱身設(shè)計(jì)受到越來越多重視?；诎殡S方程的雷達(dá)散射截面（RCS）梯度求解可以通過一次控制方程求解、一次伴隨方程求解實(shí)現(xiàn)目標(biāo)靈敏度的快速計(jì)算，是氣動(dòng)/隱身精細(xì)化設(shè)計(jì)的基礎(chǔ)。

近年來，基于電磁場積分方程伴隨方程的靈敏度計(jì)算方法及相應(yīng)的梯度優(yōu)化受到廣泛關(guān)注［3-7］?；诎殡S方程的電磁場優(yōu)化早期主要應(yīng)用在天線設(shè)計(jì)領(lǐng)域，Georgieva 等［8］在2002 年將伴隨方程引入電磁場頻域計(jì)算，研究了基于矩量法的伴隨方程求解，并基于伴隨方法對(duì)Yagi-Uda 天線進(jìn)行了優(yōu)化；Toivanen 等［9］將自動(dòng)微分法引入基于矩量法的電磁場伴隨方程的設(shè)計(jì)變量靈敏度求解中；Kataja 等［10-11］基于矩量法（MOM）［12］、多層快速多極子算法（MLFMA）［13-14］推導(dǎo)了電場積分方程的連續(xù)伴隨表達(dá)式。

基于伴隨方法的RCS 梯度計(jì)算在飛行器隱身優(yōu)化問題中研究較晚，公開文獻(xiàn)較少。Bondeson 等［15］推導(dǎo)了二維積分形式麥克斯韋方程的連續(xù)伴隨方程，并對(duì)簡單外形的RCS 進(jìn)行了優(yōu)化；周琳等［3-5］在MOM、MLFMA 的基礎(chǔ)上推導(dǎo)了混合場積分方程的離散伴隨方程，并對(duì)典型飛翼布局進(jìn)行了氣動(dòng)隱身優(yōu)化。Li 等［6］在MOM 的基礎(chǔ)上提出一種基于伴隨方程和自動(dòng)微分的高效設(shè)計(jì)變量梯度求解方法，并對(duì)機(jī)翼進(jìn)行了氣動(dòng)隱身優(yōu)化。

雖然近年來基于伴隨方程的氣動(dòng)隱身優(yōu)化取得了較大進(jìn)展，但尚無公開文獻(xiàn)探討RCS 的表面靈敏度求解方法。表面靈敏度可以反映表面網(wǎng)格位置的改變對(duì)目標(biāo)函數(shù)的影響，一方面對(duì)優(yōu)化設(shè)計(jì)有重要指導(dǎo)作用；另一方面設(shè)計(jì)變量的梯度可以根據(jù)表面靈敏度直接求出，從而克服經(jīng)典基于伴隨方程的設(shè)計(jì)變量RCS 梯度求解需反復(fù)填充阻抗矩陣的問題［3-4，8］，在經(jīng)典伴隨方程梯度求解的基礎(chǔ)上進(jìn)一步提升設(shè)計(jì)變量梯度的求解效率。飛行器氣動(dòng)優(yōu)化領(lǐng)域伴隨方法已較為成熟，已采用連續(xù)伴隨方程［16］、自動(dòng)微分方法［17］、網(wǎng)格伴隨方程等方法［18-19］對(duì)目標(biāo)的表面靈敏度的計(jì)算方法進(jìn)行了深入研究，相對(duì)于氣動(dòng)問題，電磁場求解問題的阻抗矩陣是稠密矩陣，對(duì)于表面網(wǎng)格量龐大的復(fù)雜問題，計(jì)算和存儲(chǔ)開銷較大，且電磁場問題中的阻抗積分項(xiàng)與幾何信息的依賴關(guān)系復(fù)雜，運(yùn)算多為復(fù)數(shù)運(yùn)算，難以手動(dòng)推導(dǎo)得到。以上2 點(diǎn)都使電磁場問題的表面靈敏度求解更具挑戰(zhàn)。

本文在矩量法的基礎(chǔ)上提出基于自動(dòng)微分的雷達(dá)散射截面表面靈敏度計(jì)算方法，以及適用于表面靈敏度存儲(chǔ)的稀疏矩陣存儲(chǔ)技術(shù)。該方法將基于自動(dòng)微分法的阻抗矩陣偏導(dǎo)數(shù)求解嵌入阻抗矩陣并行填充循環(huán)中。對(duì)于任意數(shù)量設(shè)計(jì)變量和入射角度，將需進(jìn)行一次阻抗矩陣、阻抗矩陣偏導(dǎo)數(shù)求解，從而使RCS 梯度計(jì)算時(shí)間與設(shè)計(jì)變量數(shù)無關(guān)；當(dāng)入射角度改變時(shí)，計(jì)算所有設(shè)計(jì)變量梯度約為8 次矩陣乘，顯著提升基于伴隨的梯度求解效率。

1 基于伴隨方程的RCS 梯度計(jì)算

典型雷達(dá)散射截面σ的減縮問題可表示為

式中：I為表面感應(yīng)電流的基函數(shù)系數(shù)；D為設(shè)計(jì)變量向量；V-ZI=0 為離散形式的麥克斯韋方程組。

引入拉格朗日乘子Λ可以構(gòu)造目標(biāo)函數(shù)L

式（2）對(duì)設(shè)計(jì)變量D求導(dǎo)，有

令式（3）右端第1 項(xiàng)為0 即可得到伴隨方程，又由

伴隨方程可以表示為

求式（5）得到伴隨變量Λ后，設(shè)計(jì)變量的RCS 梯度可以表示為

設(shè)計(jì)變量D的梯度也可以根據(jù)表面靈敏度間接求解。將式（6）中的D用r替換得到基于伴隨方程的表面靈敏度計(jì)算方法

設(shè)計(jì)變量D與表面網(wǎng)格r的關(guān)系可以表示為

式中：Np為表面網(wǎng)格點(diǎn)數(shù)量；為RCS關(guān)于表面網(wǎng)格點(diǎn)rj=［x，y，z］在3 個(gè)坐標(biāo)方向的靈敏度；為表面網(wǎng)格點(diǎn)rj與設(shè)計(jì)變量的關(guān)系，由網(wǎng)格變形算法決定。由于設(shè)計(jì)變量位置、入射角度的改變不影響?Z/?r的結(jié)果，因此梯度計(jì)算過程中只需要一次?Z/?r求解。而且基于表面靈敏度的設(shè)計(jì)變量梯度求解方法不需要反復(fù)填充阻抗矩陣Z，矩陣-矢量乘積Z的計(jì)算次數(shù)少，從而提升梯度求解計(jì)算效率。

2 RCS 表面靈敏度存儲(chǔ)方法

采用表面靈敏度計(jì)算設(shè)計(jì)變量梯度的關(guān)鍵在于表面靈敏度的計(jì)算、存儲(chǔ)及表面靈敏度與感應(yīng)電流的矩陣-矢量乘。本節(jié)首先針對(duì)表面靈敏度存儲(chǔ)量大的問題提出了表面靈敏度矩陣的稀疏存儲(chǔ)技術(shù)。

式中：Ne為基函數(shù)個(gè)數(shù)。由于σ和V的求解表達(dá)式較為簡潔，不涉及復(fù)雜的積分計(jì)算，因此?σ/?r、?V/?r可以通過自動(dòng)微分直接求解，通過向量-向量相乘計(jì)算ΛT?V/?r，計(jì)算量較小。式（9）的難點(diǎn)在于?Z/?r的求解和的計(jì)算。理論上?Z/?r是維度為Ne×Ne×Np的大規(guī)模三維矩陣，無法直接存儲(chǔ)，本節(jié)通過分析?Z?r的特點(diǎn)，提出一種點(diǎn)索引稀疏矩陣存儲(chǔ)技術(shù)，將的計(jì)算量和存儲(chǔ)量降到優(yōu)化可承受的范圍。

以電磁場積分方程求解最常使用的RWG（Rao-Wilton-Glisson）基函數(shù)［20］為例進(jìn)行推導(dǎo)，RWG 基函數(shù)定義在一對(duì)相鄰三角形的公共邊上（見圖1），第n條邊對(duì)應(yīng)的基函數(shù)fn(r)為

圖1 基函數(shù)RWG 示意圖［21］Fig.1 RWG basis function diagram［21］

由RWG 基函數(shù)的特點(diǎn)可知，每一阻抗矩陣元素Zmn至多只與8 個(gè)表面網(wǎng)格點(diǎn)的位置有關(guān)（圖2）。因此，雖然?Z/?r理論上是維度為Ne×Ne×Np的大型三維矩陣，但?Z/?r中非零元素?cái)?shù)量只有Ne×Ne×8，因此可以建立以網(wǎng)格點(diǎn)為索引的稀疏矩陣存儲(chǔ)技術(shù)，從而顯著降低?Z/?r的存儲(chǔ)需求。

圖2 Zmn存儲(chǔ)示意圖Fig.2 Zmn storage diagram

圖3 雙錐體外形及設(shè)計(jì)變量位置Fig.3 Design variables of double-ogive test case

在程序處理中，開辟整型索引數(shù)組Zdiff_index（Ne，Ne，8），并按圖2 中Index 數(shù)組中點(diǎn)的記錄順序存儲(chǔ)與阻抗元素Zmn有關(guān)的8 個(gè)網(wǎng)格點(diǎn)的全局編號(hào)。開辟復(fù)數(shù)數(shù)組Zdiff（Ne，Ne，8）存儲(chǔ)Zmn關(guān)于每個(gè)相關(guān)網(wǎng)格點(diǎn)的導(dǎo)數(shù)?Zmn/?r。若與基函數(shù)m、n相關(guān)的4 個(gè)三角形間存在重合頂點(diǎn)，則與基函數(shù)Zmn相關(guān)的網(wǎng)格點(diǎn)數(shù)量將小于8，此時(shí)在索引數(shù)組Zdiff_index 和Zdiff 數(shù)組的相應(yīng)位置補(bǔ)0。?V/?r的稀疏存儲(chǔ)方法與此非常相似，不再贅述。

3 RCS 表面靈敏度計(jì)算方法

以電場積分方程為例詳細(xì)分析表面靈敏度的計(jì)算技術(shù)。為將?Z/?x的計(jì)算嵌入阻抗矩陣Z的填充循環(huán)［22］，從而利用Z計(jì)算的中間結(jié)果，實(shí)現(xiàn)?Z/?r的高效求解，首先將Zmn的求解拆分成乘積項(xiàng)M1、M2和M2的組合，并將乘積項(xiàng)進(jìn)一步拆分成積分項(xiàng)G 和H的組合。其中乘積項(xiàng)的表達(dá)式為

積分項(xiàng)G和H的表達(dá)式為

當(dāng)場三角形與源三角形重合時(shí)積分項(xiàng)需采用奇異積分計(jì)算方法，不重合時(shí)可以采用高斯積分法。則的線性組合。

在阻抗矩陣Zmn填充循環(huán)的基礎(chǔ)上采用自動(dòng)微分法求解?Z/?r，采用Tapenade 自動(dòng)微分工具［23］。自動(dòng)微分方法通過分析程序源碼，將復(fù)雜計(jì)算分解為簡單運(yùn)算的組合，根據(jù)鏈?zhǔn)角髮?dǎo)法則得到復(fù)雜目標(biāo)的梯度，自動(dòng)微分法不產(chǎn)生截?cái)嗾`差，避免了差分法依賴步長的問題。?Z/?r的具體程序求解方法見算法1。其中最多與一條公共邊相連的2 個(gè)三角形有關(guān)（6 個(gè)網(wǎng)格點(diǎn)），當(dāng)考慮積分項(xiàng)在r=［x，y，z］3 個(gè)方向的導(dǎo)數(shù)時(shí)，分別為維度為1×3、3×3 的矩陣。

?M1/?r、?M2/?r、?M3/?r最多與一條公共邊相連的2 個(gè)三角形有關(guān)（6 個(gè)網(wǎng)格點(diǎn)），根據(jù)鏈?zhǔn)角髮?dǎo)法則，結(jié)合?G/?r、?H/?r的求解結(jié)果，可以得到

采用式（20）～式（22）結(jié)合算法1 可以得到?Z/?r在r=［x，y，z］3 個(gè)方向上的導(dǎo)數(shù)。由于實(shí)際優(yōu)化中通常重點(diǎn)關(guān)注目標(biāo)函數(shù)在z方向的梯度，因此程序?qū)崿F(xiàn)時(shí)可以僅對(duì)z方向的梯度進(jìn)行計(jì)算和存儲(chǔ)，使向量-向量乘和矩陣-向量乘為標(biāo)量乘，從而進(jìn)一步降低計(jì)算量和存儲(chǔ)量。同時(shí)，采用自動(dòng)微分工具計(jì)算積分項(xiàng)和乘積項(xiàng)的偏導(dǎo)數(shù)時(shí)，需重復(fù)計(jì)算積分項(xiàng)和乘積項(xiàng)求解的中間量，可以通過修改自動(dòng)微分得到的源碼，降低重復(fù)計(jì)算量，從而進(jìn)一步降低計(jì)算量。

完成?Z/?r的計(jì)算和存儲(chǔ)后，需根據(jù)式（9）求解目標(biāo)函數(shù)對(duì)網(wǎng)格位置的梯度dσ/dr，其核心是基于稀疏矩陣存儲(chǔ)格式計(jì)算對(duì)于每個(gè)網(wǎng)格點(diǎn)ri，可以寫為

從算法1 可以看到，提出的表面靈敏度計(jì)算方法中矩陣元素填充方法與矩量法基本一致，因此可以直接采用矩量法的并行策略［22］，較易實(shí)現(xiàn)并行求解，從而進(jìn)一步提高計(jì)算效率。激勵(lì)向量的自動(dòng)微分梯度與阻抗矩陣的梯度計(jì)算非常類似，這里不再贅述。

4 表面靈敏度計(jì)算量分析

本節(jié)重點(diǎn)關(guān)注基于經(jīng)典伴隨方法和基于表面靈敏度方法求解設(shè)計(jì)變量梯度的計(jì)算量。工程問題中的隱身優(yōu)化目標(biāo)通常為給定角域的RCS 均值，記優(yōu)化問題考慮的設(shè)計(jì)變量數(shù)為N，考慮的入射角度數(shù)為Na。2 種方法的計(jì)算量對(duì)比見表 1。

表1 設(shè)計(jì)變量梯度計(jì)算量對(duì)比Table 1 Gradient calculation amount comparison of design variables

梯度計(jì)算過程中需要擾動(dòng)設(shè)計(jì)變量，生成新外形表面網(wǎng)格并計(jì)算新外形的阻抗矩陣Z(D+ΔD)，并進(jìn)行矩陣-向量乘、向量-向量乘得到ΛTZ(D+ΔD)。由于每個(gè)入射角度的感應(yīng)電流不同，因此計(jì)算N個(gè)設(shè)計(jì)變量在Na個(gè)入射角度的梯度需要進(jìn)行N次表面網(wǎng)格變形計(jì)算、N次矩陣元素填充、N×Na次矩陣-矢量相乘運(yùn)算、N×Na次矢量-矢量相乘運(yùn)算。

對(duì)于基于表面靈敏度的梯度計(jì)算方法，表面靈敏度計(jì)算過程時(shí)需采用自動(dòng)微分法計(jì)算?Z/?r，計(jì)算量大于阻抗矩陣填充，通過修改自動(dòng)微分代碼，復(fù)用阻抗矩陣填充的中間量可以提升?Z/?r的計(jì)算效率。由于?Z/?r與設(shè)計(jì)變量數(shù)量、入射角度無關(guān)，在梯度計(jì)算過程中僅需一次?Z/?r求解，因此雖然?Z/?r的計(jì)算量較大，但對(duì)梯度計(jì)算效率影響不大。?Z/?r所需的存儲(chǔ)量約為阻抗矩陣存儲(chǔ)量的8 倍。得到?Z/?r后計(jì)算該過程需要8Na次矩陣-向量相乘和8Na次向量-向量相乘運(yùn)算。得到表面靈敏度后，根據(jù)式（8）計(jì)算設(shè)計(jì)變量梯度，計(jì)算時(shí)需依次擾動(dòng)設(shè)計(jì)變量，由網(wǎng)格變形算法（如自由曲面造型法）進(jìn)行表面網(wǎng)格變形計(jì)算得到dr/dD，與網(wǎng)格靈敏度dσ/dr相乘得到設(shè)計(jì)變量的梯度dσ/dD。在基于伴隨方程的氣動(dòng)優(yōu)化中，通常引入網(wǎng)格伴隨方程［18-19］實(shí)現(xiàn)空間網(wǎng)格的dr/dD高效求解。由于基于麥克斯韋積分方程的RCS 求解僅需要對(duì)物體的表面進(jìn)行網(wǎng)格劃分，計(jì)算所需的網(wǎng)格量遠(yuǎn)小于氣動(dòng)求解，采用有限差分法計(jì)算dr/dD計(jì)算量較小。

5 算例驗(yàn)證

本節(jié)采用雙錐體、球體、F22 外形和飛翼模型對(duì)提出的表面靈敏度計(jì)算方法進(jìn)行驗(yàn)證，并探討典型飛行器的表面敏感性特征。下文算例均以x-y平面為工作平面，計(jì)算RCS 在z方向的表面靈敏度。

1）雙錐體模型

雙椎體模型長0.190 5 m，由2 個(gè)錐角分別為22.62°、46.4°的半錐拼接而成（見圖 3），是EMCC（Electromagnetic Code Consortium）［24-25］提供用于校驗(yàn)計(jì)算電磁學(xué)代碼的標(biāo)準(zhǔn)算例之一。采用雙錐體算例對(duì)采用的電磁場求解器和提出的表面靈敏度計(jì)算方法進(jìn)行驗(yàn)證。雙錐體算例計(jì)算采用的雷達(dá)波頻率f=9 GHz，單站散射垂直極化（VV），定義沿x軸正向入射時(shí)入射角φ=0°，沿x軸負(fù)向入射時(shí)φ=180°。計(jì)算采用的網(wǎng)格平均尺寸約為λ/10，未知量總數(shù)N=11 094。

圖4 為采用矩量法求解器的計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果［25］的對(duì)比，可見計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果吻合較好，證明采用的求解器具有較高精度，滿足RCS評(píng)估和靈敏度計(jì)算的需求。

圖4 雙錐體外形RCS（f=9 GHz，VV 極化）Fig 4 RCS of double-ogive test case（f=9 GHz，VV）

采用域元法進(jìn)行表面網(wǎng)格變形，采用表面靈敏度計(jì)算得到的設(shè)計(jì)變量梯度與有限差分計(jì)算結(jié)果對(duì)比見圖5?？梢钥吹剑? 種梯度計(jì)算方法得到的結(jié)果具有良好的一致性，證明提出的表面靈敏度計(jì)算方法能夠通過表面靈敏度的計(jì)算結(jié)果獲取任意設(shè)計(jì)變量的梯度。

圖5 表面靈敏度計(jì)算梯度與有限差分結(jié)果對(duì)比Fig.5 Comparisons of design variable gradient calculation by surface sensitivity and finite difference

圖6 雙錐體外形表面靈敏度（f=9 GHz，VV 極化）Fig.6 Surface sensitivity of double-ogive（f=9 GHz，VV）

圖7 表面網(wǎng)格點(diǎn)梯度對(duì)比Fig.7 Gradients of surface mesh

圖8 球體算例的表面電流和表面靈敏度（f=1 GHz，VV）Fig.8 Surface current and surface sensitivity（f=1 GHz，VV）

圖9 表面靈敏度計(jì)算中間量分布（f=1 GHz，VV）Fig.9 Intermediate quantity results（f=1 GHz，VV）

圖10 類F22 戰(zhàn)斗機(jī)表面靈敏度（沿鼻錐方向入射）Fig.10 Surface sensitivity of F22（incident towards nose）

圖11 典型飛翼布局表面靈敏度（沿鼻錐方向入射）Fig.11 Surface sensitivity of flying wing（incident towards nose）

采用本文方法計(jì)算φ=0°，90°的表面靈敏度見圖 6。與表面電流的分布特點(diǎn)相似，單一入射頻點(diǎn)和入射角度的表面靈敏度具有正負(fù)交替的條紋，且正負(fù)條紋間隔約為λ4。

對(duì)比采用經(jīng)典伴隨方法和表面靈敏度方法求解設(shè)計(jì)變量梯度的計(jì)算效率和內(nèi)存需求。計(jì)算采用Intel i7-9700 處理器8 核心并行計(jì)算。在計(jì)算效率方面，雙錐體外形的阻抗矩陣Z填充用時(shí)11 s；僅計(jì)算z方向靈敏度時(shí)填充?Z/?r和Z矩陣用時(shí)35 s；除去求解控制方程和伴隨方程的耗時(shí)，采用經(jīng)典伴隨方法計(jì)算36 個(gè)設(shè)計(jì)變量在180 個(gè)入射角度的梯度用時(shí)413 s，采用基于表面靈敏度的梯度計(jì)算方法用時(shí)19 s。在內(nèi)存開銷方面，經(jīng)典矩量方法內(nèi)存占用1.78 G，表面靈敏度方法內(nèi)存占用16.72 G。可以看到，雖然偏導(dǎo)數(shù)矩陣?Z/?r的填充速度較慢，但由于基于表面靈敏度方法求解設(shè)計(jì)變量梯度時(shí)不需要填充擾動(dòng)外形的阻抗矩陣，因此當(dāng)設(shè)計(jì)變量較多時(shí)效率具有優(yōu)勢。

2）球體模型

球體模型半徑r=0.2 m，計(jì)算入射頻率f=1 GHz，單站散射垂直極化（VV），表面網(wǎng)格尺寸Δx=0.05 m；表面網(wǎng)格的總數(shù)N=197。

電磁波沿x軸負(fù)向入射時(shí)，采用經(jīng)典伴隨法（Adjoint_diff）和本文表面靈敏度計(jì)算方法（Adjoint_surface）求解模型在z方向的表面靈敏度dσ/dz，計(jì)算結(jié)果對(duì)比見圖 7?？梢钥闯? 種方法計(jì)算得到的梯度完全重合，證明本文表面靈敏度方法具有較高精度。表面感應(yīng)電流和表面靈敏度計(jì)算結(jié)果見圖 8，表面靈敏度空間分布特征和表面電流存在一定差異，對(duì)于表面感應(yīng)電流，入射波直接照射的區(qū)域明顯高于其他區(qū)域，而表面靈敏度的幅值與電流強(qiáng)度不直接相關(guān)，感應(yīng)電流強(qiáng)的區(qū)域不一定具有更大的表面靈敏度，因此在優(yōu)化過程中難以直接根據(jù)電流的分布特點(diǎn)判斷區(qū)域的RCS 減縮潛力和方向，需進(jìn)行表面靈敏度求解。

為進(jìn)一步探討表面靈敏度的特性，從式（9）可以看到，表面靈敏度dσ/dz的計(jì)算結(jié)果由?σ/?z、線性疊加得到，計(jì)算電磁波沿x軸負(fù)向入射、垂直極化時(shí)每一部分的空間分布特征，計(jì)算結(jié)果如圖 9。可以看到，3 部分計(jì)算結(jié)果均具有明顯的相位特征，空間分布與dσ/dz基本一致，與感應(yīng)電流有差異顯著。

3）F22 戰(zhàn)斗機(jī)外形

分析典型戰(zhàn)斗機(jī)的RCS 表面靈敏度特征，采用F22 戰(zhàn)斗機(jī)外形，計(jì)算模型機(jī)身總長約19 m，半展長約7 m?？紤]入射頻率f=200 MHz，計(jì)算網(wǎng)格未知量數(shù)N=27 595。

沿鼻錐方向入射時(shí)水平極化（HH）和垂直極化（VV）的表面靈敏度計(jì)算結(jié)果見圖 10?？梢钥闯?，與雙錐體和球體算例相似，2 種極化方式入射時(shí)表面靈敏度均具有正負(fù)相間的特征，但2 種極化的表面靈敏度空間分布存在較大差異。對(duì)于垂直極化，戰(zhàn)斗機(jī)機(jī)身的表面靈敏度大于機(jī)翼的表面靈敏度，且厚度較大的區(qū)域（如座艙）靈敏度更高，平尾邊緣的靈敏度較大。對(duì)于水平極化，靈敏度較大的區(qū)域主要為機(jī)身和機(jī)翼的邊緣，平尾的表面靈敏度較小，水平極化的靈敏在幅值上整體略小于垂直極化的表面靈敏度，且與幾何的厚度沒有明顯關(guān)聯(lián)。表面靈敏度的計(jì)算結(jié)果表明，靈敏度的空間分布與極化方式關(guān)聯(lián)明顯，極化方式的RCS 減縮需重點(diǎn)關(guān)注的區(qū)域有所區(qū)別。

4）典型飛翼外形

計(jì)算采用的飛翼布局根弦長約12 m，半展長約10 m，采用NACA65016 對(duì)稱翼型生成?？紤]入射頻率f=200 MHz，計(jì)算網(wǎng)格的未知量總數(shù)N=28 086，分別計(jì)算入射波為水平極化（HH）、垂直極化（VV）時(shí)的表面靈敏度。

沿鼻錐方向入射時(shí)的表面靈敏度計(jì)算結(jié)果見圖 11，其中水平極化的靈敏度明顯高于垂直極化的靈敏度，2 種極化方式的靈敏度幅值的差異主要受當(dāng)前狀態(tài)的RCS 影響，當(dāng)沿鼻錐入射時(shí)，兩種極化的RCS 分別為RCSHH=-4.9 dB、RCSVV=-21.0 dB，RCS 與靈敏度正相關(guān)。垂直極化靈敏度較大的區(qū)域主要為中央體上下表面，厚度越大的區(qū)域靈敏度越大，該現(xiàn)象與F22戰(zhàn)斗機(jī)的計(jì)算結(jié)果一致；而水平極化靈敏度較大的區(qū)域?yàn)闄C(jī)翼的前后緣，中央體前緣的靈敏度顯著高于其他區(qū)域。

圖12 分別為y=6 m 站位VV、HH 極化方式入射時(shí)表面靈敏度的計(jì)算結(jié)果，由于布局采用對(duì)稱翼型，2 種極化方式入射時(shí)上下表面的靈敏度符號(hào)相反、幅值幾乎一致，且都表現(xiàn)出了強(qiáng)烈的正負(fù)震蕩特性，正負(fù)波峰的間隔約為入射波長的1/4。

圖12 y=6 m 表面靈敏度分布Fig.12 Surface sensitivity at y=6 m

圖13 為入射波水平極化（HH）時(shí)前向威脅扇區(qū)±30°角域的平均感應(yīng)電流和表面靈敏度?？梢钥吹?，機(jī)翼前緣的平均感應(yīng)電流較強(qiáng)，從前緣向后緣電流幅度逐漸降低。角域平均靈敏度仍表現(xiàn)出正負(fù)相間的特征，但由于不同角度入射時(shí)表面靈敏度存在相位疊加，平均表面靈敏度的明暗條紋不如單一角度入射時(shí)顯著。圖 14 為Slice 1、Slice 2 剖面（見圖 11）的上下表面設(shè)計(jì)變量在前向威脅扇區(qū)的平均梯度，其中網(wǎng)格變形采用自由曲面造型技術(shù)（FFD），設(shè)計(jì)變量位置如圖 11。雖然表面靈敏度呈現(xiàn)出顯著正負(fù)交替的特征，設(shè)計(jì)變量梯度變化較為和緩。上下表面設(shè)計(jì)變量梯度不完全對(duì)稱是由于設(shè)計(jì)變量上下表面在z方向的位置不完全對(duì)稱造成。

圖13 角度平均感應(yīng)電流和平均表面靈敏度（f=200 MHz，HH）Fig.13 Angle-averaged current and sensitivity（f=200 MHz，HH）

圖14 Slice 1 和Slice 2 上下表面設(shè)計(jì)變量平均梯度Fig.14 Angle-averaged design variable gradients at Slice 1 and Slice 2

6 結(jié)論

提出了一種基于伴隨方程和自動(dòng)微分的雷達(dá)散射截面表面靈敏度計(jì)算方法，并結(jié)合基函數(shù)特點(diǎn)提出了適用于表面靈敏度的稀疏矩陣存儲(chǔ)和矩陣-矢量相乘方法，降低了表面靈敏度的計(jì)算、存儲(chǔ)開銷。該方法可以通過一次阻抗矩陣偏導(dǎo)數(shù)求解得到所有網(wǎng)格點(diǎn)的表面靈敏度，當(dāng)改變?nèi)肷浣嵌葧r(shí)，求解所有設(shè)計(jì)變量梯度的計(jì)算量約為8 次矩陣-矢量乘，極大提高了基于伴隨方程的梯度計(jì)算效率。

計(jì)算典型入射條件下雙錐體、球體、F22 戰(zhàn)斗機(jī)、典型飛翼布局的表面靈敏度，證明了提出的表面靈敏度計(jì)算方法的準(zhǔn)確性，該方法能夠通過表面靈敏度快速完成設(shè)計(jì)變量梯度計(jì)算。表面靈敏度的空間分布特征受極化方式影響顯著，具有明顯的相位特征，可以為隱身優(yōu)化設(shè)計(jì)提供指導(dǎo)。

本文僅研究了基于矩量法的表面靈敏度計(jì)算方法，但是矩量法的計(jì)算量和存儲(chǔ)量較大，限制了其在電大尺寸問題中的應(yīng)用，后續(xù)將進(jìn)一步研究將該方法擴(kuò)展到基于多層快速多極子算法的RCS 表面靈敏度求解中。