飛騰處理器上向量三角函數的設計實現與優化

沈 潔 龍 標 姜 浩 黃 春

(國防科技大學計算機學院 長沙 410073)

三角函數在導航、測繪、工程、物理學等科學計算類應用中使用廣泛，由于三角函數算法實現的復雜性，這些函數往往是應用中最耗時的浮點計算函數，因此如何提高三角函數性能成為科研人員設計和優化浮點計算函數的一個重點.

與此同時，單指令多數據(single instruction multiple data, SIMD)向量化對發揮當今處理器浮點計算性能至關重要[1-2].從x86處理器上的SSE(streaming SIMD extensions)，AVX(advanced vector extensions)到最新的AVX512F(AVX-512 founda-tion)向量指令集，從ARM處理器的ADVSIMD(advanced SIMD)到SVE(scalable vector extension)向量指令集[3]，可以看出隨著體系結構的發展，處理器向量寬度正在逐步加寬、并變得更加靈活可用.然而，目前主流處理器平臺上的SIMD向量指令集僅包括對基礎的加、減、乘、除等算術運算以及邏輯運算、訪存操作、位運算等基本操作的支持，對于復雜的浮點三角函數(以及其他超越函數)沒有提供直接的支持.

科研人員通常采用軟件方式設計數學庫以實現三角函數等超越函數.當前，應用廣泛的數學庫有GNU libm庫(GNU math library)[4]、Intel libimf庫(Intel math library)[5]、Intel SVML庫(Intel short vector math library)[6]等.其中GNU libm庫和Intel libimf庫是標量數學庫，未能利用SIMD向量部件提高函數性能，而Intel SVML庫雖然是向量數學庫，但屬于Intel設計和開發的商用數學庫，僅面向x86處理器平臺實現向量化.近來，科研人員提出了Sleef開源向量數學庫[7-8]，該數學庫具有良好的模塊化和跨平臺特性，可支持多種SIMD向量指令集，在多種處理器平臺上實現超越函數向量化.然而，Sleef為了實現向量化，以增加冗余計算為代價，換取函數中分支語句的減少.這樣，雖然增加了SIMD向量部件的利用率，但冗余計算在一定程度上還是導致了性能損失.

因此，本文基于Sleef開源向量數學庫，面向飛騰處理器平臺，進一步優化向量三角函數的性能.首先，我們通過完備深入的測試，發現Sleef向量三角函數實現中的性能瓶頸.然后，我們借鑒GNU libm庫實現方法中的部分分支算法，針對性地提出了一種支持分支處理的向量三角函數優化方法.基于此，我們還使用了分支預測、精度削減、查表優化、條件傳送優化、精度修正、性能修復等優化技術以盡可能減少引入分支后帶來的向量化開銷.本文所提優化方法結合Sleef向量數學庫和GNU libm標量數學庫各自的優勢，有效提高了向量三角函數的性能.最后，我們將所提優化方法實現為Sleef的一個模塊，移植入Sleef數學庫中.我們在飛騰處理器上進行實驗驗證，結果表明我們提出的優化方法在保證向量三角函數精度的前提下，有效提高了向量三角函數性能，與原始Sleef向量三角函數相比，加速比可達1.12～3.97倍.

1 相關背景

1.1 三角函數及其軟件實現

三角函數一般通過規約、近似和重建3個步驟獲得近似值[9-10].1)規約.根據三角函數的對稱性和周期性，將任意區間映射到[0,π/4]區間.2)近似.通過泰勒級數展開計算[0,π/4]區間的三角函數近似值.例如，sinx可以計算為

由于sinx在該區間內是收斂的，因此通過展開若干次可以獲得sinx的近似值.3)重建.任意區間三角函數值可以通過[0,π/4]區間三角函數值表示.此外，對于一些特殊值的三角函數，可以繞過上述常規過程進行特殊處理(例如0、無窮大等，可以不經過計算直接返回結果).

在通用處理器上，三角函數等超越函數通常以軟件方式集成在數學庫中實現.程序員調用數學庫中的三角函數接口，由數學庫中的相關代碼實現相應的函數運算.目前，常用的數學庫有Linux操作系統上的GNU libm庫、面向Intel處理器的Intel libimf庫、ARM發布的libmathlib庫(ARM optimized routines)[11]等，其中GNU libm和ARM libmathlib是開源數學庫，Intel libimf是Intel特別為其處理器設計和優化的商用數學庫.這些數學庫都屬于標量數學庫，一次調用計算一個浮點操作數.因此，為了充分利用處理器上的向量部件和向量指令集以優化程序中的三角函數計算，科研人員設計和開發了向量三角函數.

1.2 向量三角函數及其軟件實現

相比標量三角函數，向量三角函數內部使用SIMD指令實現向量化，根據處理器向量部件的寬度，一次調用同時計算多個浮點操作數，因而其理論性能是標量三角函數的數倍.在實現時，向量三角函數仍保留規約、近似和重建3個核心步驟，并使用SIMD指令對應的內部函數(intrinsics)或直接使用匯編指令進行向量化.

目前比較成熟的實現是Intel SVML數學庫，這是Intel開發的一套商用向量數學庫，雖然已經部分開源并納入GNU的標準C庫(GNU glibc)中[12]，但其應用依然局限于x86平臺.近來，科研人員相繼實現并發布了一些開源向量數學庫.其中，奈良先端科學技術大學院大學開發的Sleef開源向量數學庫，為多種硬件平臺上的向量超越函數實現提供了支持.Lauter[13]設計的向量libm庫(vector libm)實現了雙精度三角函數，但其計算精度較低.歐洲核子研究組織發布了VDT向量數學庫(vectorised math)，其實現依賴于第三方庫進行多項式近似，因而其可向量化的程度受限[14].本文是基于Sleef向量數學庫設計和優化向量三角函數的，主要針對12個三角函數進行設計實現與優化，這12個函數具體如表1所示：

Table 1 Vector Trigonometric Functions to be Optimized表1 待優化向量三角函數

此外，硬件實現超越函數雖然能很快返回計算結果，但僅有少量非通用體系結構提供超越函數的硬件實現[15]，因而應用受限，不屬于本文的研究范疇.

1.3 Sleef向量數學庫

Fig. 1 The code structure of Sleef vector math library圖1 Sleef向量數學庫代碼結構

Sleef是一款基于C語言編寫的多平臺向量數學庫，提供了單精度和雙精度2個版本，并實現了C99標準下所有實浮點數學函數的向量化.Sleef支持多種處理器架構和SIMD指令集，目前Sleef支持x86(SSE2，SSE4.1，AVX，FMA4，AVX2+FMA3，AVX512F指令集)，AArch64(ADVSIMD，SVE指令集)和PowerPC64(VSX指令集).

由于各SIMD向量指令集的向量類型和intrinsics函數存在較多差異，為了支持多種處理器平臺，Sleef采用模塊化封裝和抽象，面向不同指令集抽象出一套通用的Sleef向量類型(包括vdouble，vfloat，vint，vint2，vopmask，vmask共6種向量類型)和Sleef intrinsics函數(包括向量浮點和整型算術運算、邏輯運算、位運算、類型轉換、條件傳送等函數).Sleef代碼結構如圖1所示，底層在頭文件中實現Sleef向量類型和Sleef intrinsics函數到各指令集向量類型和intrinsics函數的映射(一個指令集對應一個頭文件)，上層使用Sleef向量類型和Sleef intrinsics函數實現各超越函數算法，編譯器在編譯時會根據所在處理器平臺自動將Sleef向量類型和Sleef intrinsics函數替換為指定指令集的向量類型和intrinsics函數，從而達到編寫一份Sleef超越函數算法、在不同處理器平臺上編譯運行的目的.正是這種封裝和模塊化設計使得Sleef具有良好的可擴展性和可移植性，只要符合Sleef代碼規則，就可以方便地增加和修改向量數學函數或增加對新的處理器平臺和指令集的支持.本文基于Sleef向量數學庫在飛騰處理器平臺上優化向量三角函數.

1.4 數值結果精度和誤差

為了更加精確地分析向量三角函數實現的精度，本文采用ulp(unit in the last place)衡量計算結果誤差.根據Goldberg定義[16-17]，ulp是最接近實數真實值的2個浮點數之間的距離.假設函數真值為實數x，其對應的浮點數計算結果X可表示為d0.d1d2d3…dp-1βe，x∈[βe,βe+1)，則

ulp(x)=βmax(e,emin)-p+1，

其中,β是基，e是階數，d0.d1d2d3…dp-1是尾數，p是精度即浮點數尾數的有效位數(在IEEE754標準下單精度p=24，雙精度p=53).

本文實驗中，X由待測三角函數計算得出，x由任意精度算法庫MPRF庫[18]中對應的三角函數計算同樣的操作數得出，通過計算x與X的絕對誤差|x-X|的單位ulp數，從而評價待測三角函數的精度.一般情況下，三角函數等超越函數的精度需控制在0.5ulp～1.0ulp以內.

Fig. 2 Vector trigonometric function optimization圖2 向量三角函數優化方法

由于科學計算類應用中存在大量可并行的三角函數數值計算，這些計算十分耗時，因此本文將面向飛騰處理器，基于Sleef向量數學庫，在保證一定精度的前提下優化向量三角函數性能.

2 向量三角函數設計與優化

2.1 優化方法

我們調研了向量數學庫中三角函數的實現方法，發現雖然在算法上向量三角函數仍然采用規約、近似和重建3個步驟，但會盡量避免使用標量三角函數實現中廣泛采用的分支處理，這是因為如果同一向量中的不同元素進入不同分支，那么向量化將受到影響，不同分支中的元素不能同時進行運算.Sleef向量數學庫正是如此，為了充分利用SIMD指令的并行處理能力，減少分支操作，Sleef向量數學庫在實現向量三角函數算法時不區分操作數的特征，采用統一算法處理所有操作數.該實現方法雖然盡可能避免了算法中的分支語句，提升了處理器SIMD向量單元的利用率，但也導致了向量三角函數在處理一些操作數時的冗余計算(例如cos函數在計算絕對值小于2-27的操作數時，可以不需要經過計算而直接將1.0作為計算結果返回).相反，libm數學庫作為標量數學庫，為了盡可能減少函數串行計算的計算量，對三角函數的操作數進行了分段，利用三角函數的性質，對不同的操作數范圍采用不同的算法(分支)實現運算，從而提高函數運算的速度.

本文將結合Sleef向量化和libm分段計算的優勢，以libm的分段計算為指導，將分段計算引入Sleef向量數學庫中，從而進一步提高向量三角函數的性能.換言之，我們提出一種在向量數學庫中使用分支以優化函數性能的方法.本文的基本優化方法如圖2所示.

針對每個三角函數，首先我們對libm三角函數中的操作數范圍分支進行篩選，當某個或某些操作數范圍分支的計算量較小、執行速度較快時(也就是說該分支值得被向量化用以替換原始Sleef函數)，則將其提取出來并使用Sleef intrinsics函數進行向量化.在實際應用中，優化后的三角函數在運行時會判斷向量操作數opt中各元素是否均在可優化操作數范圍內，如果是則進入向量化的libm分支，否則進入原始Sleef分支，最后計算出向量結果ret返回.對libm函數分支進行篩選時，既要保證函數運算的時間優于原始Sleef函數運算時間，又要保證函數運算的結果滿足精度要求(最大誤差控制在1.0ulp以內).在實現時，我們通過大量實驗確定每個三角函數可移植的操作數范圍分支，同時在向量化libm函數分支的過程中，我們還使用了分支預測、精度削減、查表優化、條件傳送優化、精度修正和性能修復一系列優化技術以進一步加快函數運算速度.

1) 分支預測.如圖2所示，我們的優化方法是以分支的形式添加至原始Sleef代碼中的，這就在代碼中引入了執行分支語句、進行分支預測所帶來的開銷.為了盡可能減少分支語句對函數性能的影響，降低分支預測失敗的開銷，我們使用各個編譯器支持的分支預測函數(如gcc編譯器的_builtin_expect函數)對2個分支執行的概率進行調整，以大概率執行原始Sleef分支，小概率執行向量化的libm分支，這是因為我們通過實驗確定的每個函數的可優化操作數范圍相對于整個實數域仍是一個小范圍，因此我們認為設定以大概率執行原始Sleef分支是合理的.只有一個例外是atanf函數，因為經過實驗我們發現其可優化操作數范圍是整個實數域(詳見2.2節)，因此實現時是將原始Sleef分支整個替換為向量化的libm分支，因此不會受到分支預測優化的影響.Sleef庫中已經將各個編譯器的分支預測函數進行了打包，函數接口分別是LIKELY()和UNLIKELY().

2) 精度削減.libm的部分三角函數在一些操作數范圍內的計算精度很高(最大誤差在0.5ulp以內)，代碼中將對計算結果的精度進行判斷，如果精度不夠則進入下一層分支進行細化計算，從而達到更高的精度要求.由于Sleef的精度控制相比libm寬松一些，因此在移植libm分支時我們可以刪除一些高精度分支，以降低計算開銷.

3) 查表優化.標量三角函數實現時通常將函數計算所用到的常系數寫在表中，并使用查表操作以減少函數計算量.原始Sleef向量數學庫在實現時盡量規避了查表操作，這是因為進行向量查表會涉及到使用聚合(gather)訪存指令，將內存中不同位置的常系數值加載到一個向量寄存器中，供向量中的不同元素使用，這種聚合訪存因訪存不連續往往需要多次訪存操作才能完成，從而影響三角函數向量化的性能.我們發現當函數中需要進行查表操作，而查表所需表的大小小于向量寬度時，仍然可以使用查表操作加速函數計算.首先，使用向量load指令先將表中的數據從內存一次加載至向量寄存器中，然后使用特定的permute，blend，shuffle向量操作指令將這些數據存放至向量中的指定位置，從而代替向量gather操作，這樣既保留了查表操作并減少了函數計算量，同時又降低了向量查表可能引入的多次內存讀寫的開銷.目前支持特定permute，blend，shuffle向量操作的指令集有AVX2，AVX512F，SVE.

4) 條件傳送優化.向量代碼中可以使用類似三目運算符x=a?b:c功能的簡單條件運算，編譯器編譯時會將其翻譯成向量條件傳送指令.例如ADVSIMD指令集中的BSL，BIT，BIF指令.本質上這些指令可以等價為若干按位邏輯操作的組合，如圖3所示：

Fig. 3 Conditional move operation is implemented by bit operation圖3 條件傳送操作實現為按位邏輯操作

當優化代碼中使用這種條件傳送運算的操作數是一些特殊值時(比如1個源操作數為0或者2個源操作數互為相反數)，我們可以使用更加簡單的按位邏輯操作組合而不是直接使用向量條件傳送語句，從而達到提升性能的目的.向量條件傳送優化如圖4所示.

Fig. 4 Conditional move operation optimization圖4 條件傳送優化

5) 精度修正.當計算某幾個特定操作數結果有誤時(如atanf函數在計算|x|≥234的操作數時)，如圖5我們將在代碼最后使用向量條件傳送語句來修正結果向量中的錯誤元素.

6) 性能修復.我們發現優化代碼在處理可優化操作數范圍內的一些特殊操作數時，在運行時可能會進入非規格化數處理的情況中而變得十分緩慢.為了解決這個問題，我們使用了類似于精度修正的方法來規避這種潛在的性能問題.如圖6所示，我們的優化代碼在開始計算之前會先判斷操作數是否落入特殊操作數范圍中，并使用條件傳送語句將操作數中可能導致性能問題的元素替換為普通操作數.計算結束后再使用同樣的條件傳送語句將計算結果(即錯誤的結果)替換為正確值.由于導致被優化函數陷入非規格化數的緩慢分支對應的特殊操作數范圍可以確定出來，其計算結果往往可以被近似為一個常量，因此可以提前計算得出正確的計算結果，并在代碼中直接使用.

需要指出的是，原則上我們盡可能對每個向量三角函數使用上述所有優化方法，但由于各個函數內在算法的差異性，上述優化方法不一定適用于每一個函數.其中，分支預測優化應用于所有函數中，精度削減應用于所有雙精度函數中，條件傳送優化和精度修正只使用在asinf，acosf，atanf函數中，性能修復只使用在atanf函數中，而查表優化則因為ADVSIMD指令集不支持permute操作而只使用在了SVE指令集的atanf函數中.

Fig. 5 Precision correction圖5 精度修正

Fig. 6 Performancefixing圖6 性能修復

2.2 具體實現

結合Sleef向量數學庫和libm標量數學庫，我們一共對12個三角函數(包括單精度和雙精度)的向量版本進行了優化.我們通過大量實驗確定了各個函數可優化的操作數范圍，具體如表2所示.我們發現如果libm某個函數的某個分支代碼中包含大量的分支跳轉語句，該分支將不值得移植和向量化，除非這些條件跳轉比較簡單可以轉換為條件傳送.另外，atanf函數的可優化操作數范圍是整個實數域，因為我們發現向量化后的libm函數在整個范圍內的性能均優于原始Sleef函數.

在代碼移植和向量化的過程中，我們在原Sleef向量數學庫的基礎上增加了7個新的Sleef intrinsics函數，如表3所示.其中，vint2，vdouble，vfloat分別代表整型向量類型、雙精度浮點數向量類型和單精度浮點數向量類型，vgatherr_vi2_p_vi2()與vgatherr_vf_p_vi2()為前文所提的為查表優化打包的函數.

Table 2 The Optimization Range of Vector Trigonometric Functions表2 各向量三角函數可優化操作數范圍

Table 3 New Sleef Intrinsics表3 優化過程中新增加的Sleef Intrinsics函數

優化的代碼是以模塊化方式移植入Sleef向量數學庫中，從而保證這部分優化代碼的獨立性，也方便后續Sleef更新后，能將優化代碼快速合并到新的Sleef中.具體實現時，我們將每個函數的優化代碼放置于一個單獨xxxx_opt.c文件中，例如單精度函數放置于opt_include/single/xxxx_opt.c文件中，而雙精度函數放置于opt_include/double/xxxx_opt.c文件中，所有優化函數共同需要的常量和函數則統一放置于頭文件opt_include/opt_vecf_helper.h中，為函數優化新增的intrinsics函數則放置于opt_include/sleef_arch_patch.h中，通過宏定義和include語句將這些文件導入Sleef源碼中，在編譯時通過控制編譯選項-DCMAKE_C_FLAGS=“-DOPT_LB”來選擇是否開啟優化，圖7顯示了模塊化優化方式的偽代碼.

最后需要說明的是，雖然我們面向飛騰處理器進行實驗，但本文所述優化方法和基于Sleef的模塊化實現也可以擴展至其他處理器平臺和向量指令集.

Fig. 7 Porting optimized code into Sleef圖7 以模塊化方式植入優化代碼

3 實驗結果與分析

3.1 實驗設計與環境

本文結合Sleef和libm對向量三角函數進行優化，并在飛騰處理器上對優化后的向量三角函數進行功能驗證、精度實驗和性能實驗.其中，功能驗證是驗證優化后的向量三角函數功能是否正確，精度實驗是具體分析優化后向量三角函數計算精度的變化，性能實驗是對比我們提出的優化方法是否能進一步提升向量三角函數的性能.

本文實驗主要在飛騰處理器上進行(FT-1500A/16，兼容ADVSIMD指令集，向量寬度為128 b).此外，為了進行全面的功能驗證，我們還使用了x86處理器平臺以驗證SSE2，AVX，AVX2指令集，以及ARM Instruction Emulator模擬器以驗證SVE指令集.我們使用的Sleef版本為3.3.1版本，libm版本為GNU Glibc Linaro-2.20版本.

3.2 功能驗證

我們對所有優化后的向量三角函數進行嚴格的功能性驗證，具體包括6項驗證內容.

1) 優化后函數在可優化操作數范圍內的計算精度是否滿足設計時要求.該項驗證要求待測試向量操作數中的各個元素不一致，且最終每個元素對應的計算結果的最大誤差控制在1.0ulp以內.

2) 優化后函數在可優化操作數范圍邊界及其附近的計算精度是否滿足設計時要求.同樣要求待測試向量操作數中的各個元素不一致以及最終結果的最大誤差控制在1.0ulp以內.這是用以檢驗被優化函數在算法邊界處的處理是否正確.

3) 優化后函數在可優化操作數范圍外的計算精度是否仍然滿足設計時要求，這是用以檢驗本文提出的優化方法不會對可優化操作數范圍外的計算產生精度方面的負面影響.

4) 優化后函數對浮點數域內的特殊浮點常量的計算結果是否正確.包括驗證對浮點數±0，±∞，NaN的計算結果是否等于對應的C99標準下規定的浮點值，驗證±1.0、正負最大最小非規格化數、正負最大最小規格化數的計算結果是否正確，以及驗證被優化函數對大浮點數的處理是否正確(所有雙精度函數檢驗操作數為±1E+10的計算結果，所有單精函數檢驗操作數為±1E+7的計算結果).

5) 優化后函數在計算與算法相關的特殊浮點常量時的計算結果是否正確.其中，雙精度sin，cos，tan函數和單精度sinf，cosf，tanf函數需要檢查操作數為π/4的倍數時的計算結果，雙精度asin，acos函數和單精度asinf，acosf函數需要檢驗操作數趨于±1.0時的計算結果，這些計算結果的誤差要控制在1.0ulp以內.

6) 多種情況組合下的功能驗證，將前5項中列出的各種情況下的浮點操作數組合到同一向量中，測試每個向量元素的計算結果是否達到精度要求(同樣要滿足最大誤差在1.0ulp以內)，這是用以檢驗優化函數中的分支處理與特殊值處理等是否正確.

我們在飛騰處理器平臺、x86處理器平臺和ARM Instruction Emulator模擬器上對不同SIMD向量指令集進行功能驗證.實驗證明對于飛騰處理器平臺的ADVSIMD指令集、x86平臺的SSE2，AVX，AVX2指令集、以及ARM Instruction Emulator模擬器上的SVE指令集，優化后的向量三角函數均通過了功能驗證，這表明我們對向量三角函數的優化達到了Sleef庫的標準，可以移植入Sleef庫中使用.

3.3 精度實驗

對每個優化后的向量三角函數，精度實驗以指定步階從小到大測試可優化操作數范圍內的浮點數，我們使用任意精度MPFR庫作為參考標準，計算優化后三角函數的計算結果與MPFR庫中相同函數的計算結果的具體誤差(以ulp形式表示).實驗將循環迭代次數控制在107次，并統計出最大誤差和平均誤差的ulp值.我們在飛騰處理器上進行精度實驗，結果如表4所示：

Table 4 Precision Experiment Results表4 精度實驗結果

可以看出，相比優化前、優化后的向量三角函數在可優化操作數范圍內的精度變化很小，最大誤差的變化在0.37ulp以內，平均誤差的變化在0.17ulp以內，此外，除sinf，tanf，acosf函數外，其余函數在優化后的最大誤差還有所降低(最大降低了0.37ulp).這說明本文所提出的優化技術對函數精度的影響較小，函數精度變化是符合要求的.

3.4 性能實驗

我們從2個方面分析評估本文所提出的向量三角函數優化方法：1)通過與libm標量三角函數以及原始Sleef向量三角函數進行對比，以分析優化方法對向量三角函數的性能提升；2)分析所提優化方法所引入的開銷.所有性能實驗在飛騰處理上進行.

3.4.1 優化方法帶來的性能提升

對于每個向量三角函數，我們在可優化操作數范圍內產生單精度或雙精度浮點隨機數以測試函數性能.實驗中，我們一共產生108個浮點隨機數(存放在數組中)，依次使用數組中的隨機數(以向量寬度為單位，雙精度為2個隨機數，單精度為4個隨機數)做三角函數計算，并取函數計算的平均耗時作為衡量函數性能的指標.

Fig. 8 Performance comparison of original Sleef functions and libm functions圖8 Sleef函數優化前與libm函數性能對比

首先，我們對比分析Sleef向量三角函數優化前與libm標量三角函數在可優化操作數范圍內的性能差異，結果如圖8所示.可以看出，盡管使用ADVSIMD向量指令，除異常的asinf和acosf函數外，只有一部分Sleef向量三角函數(sinf，cos，cosf，tanf)比libm標量三角函數快(最高只有1.64倍加速比)，另外一部分Sleef向量三角函數與libm性能相當(sin和atanf)，其余函數的性能比libm差(tan，asin，acos，atan)，這是因為原始Sleef向量三角函數在計算可優化操作數范圍時的算法冗余而libm采用分段計算在可優化操作數范圍內降低了函數計算量.對于asinf和acosf函數，libm性能低下的原因是我們使用的libm版本內部實現時使用了軟件開方(sqrt)，導致函數異常耗時.

其次，我們對比分析Sleef向量三角函數優化后與libm標量三角函數的性能差異，結果如圖9所示.我們發現向量三角函數優化后與libm標量三角函數相比，雙精度下平均加速比為1.74倍，單精度下平均加速比為3.06倍.這一性能表現一定程度上符合飛騰處理器128 b向量寬度下向量函數與標量函數的理想加速比(雙精度下2倍加速比，單精度下4倍加速比).

Fig. 9 Performance comparison of optimized Sleef functions and libm functions圖9 Sleef函數優化后與libm函數性能對比

最后，我們對比分析Sleef向量三角函數優化前與優化后的性能差異，結果如圖10所示.可以看出使用本文所提出的優化技術，在可優化操作數范圍內，所有向量三角函數的性能均得到了提升，平均加速比為2.04倍(最低加速比為1.12倍，最高可達3.97倍).

Fig. 10 Performance comparison of Sleef functions before and after optimization圖10 Sleef函數優化前與優化后性能對比

3.4.2 優化方法引入的開銷

為了評估所提優化方法引入的開銷，我們進行如下實驗,對比分析Sleef向量三角函數優化前與優化后在可優化操作數范圍之外的性能，以檢查我們的優化方法在非可優化操作數范圍內對原始Sleef向量三角函數的性能影響，實驗結果如圖11所示.可以看出使用本文所提優化方法后，在非可優化操作數范圍內，所有三角函數的性能變化均很小，其中雙精度三角函數的性能波動維持在4.3%以內，單精度三角函數的性能波動維持在3.8%以內(由于atanf函數的可優化操作數范圍為整個實數域，因此不存在對比).換言之，我們實驗確定的可優化操作數范圍是正確的，并且我們的優化方法開銷很低，在可優化操作數范圍之外對原始Sleef向量三角函數的性能影響很小.

Fig. 11 Performance comparison of Sleef functions before and after optimization with operands outside the optimization range圖11 在可優化操作數范圍外，Sleef函數優化前與優化后性能對比

綜上所述，功能驗證、精度實驗和性能實驗證明，本文提出的向量三角函數優化方法，在功能上滿足設計要求，能夠在保證計算精度和較小開銷的前提下有效提高向量三角函數的性能.

4 結束語

三角函數等超越函數是科學計算類應用中最耗時的浮點計算函數，利用SIMD向量化技術提高函數性能是設計和優化數學庫的一個重要研究方向.本文基于Sleef開源向量數學庫，以libm分段計算為指導，設計和提出了一種支持分支處理的向量三角函數優化方法，該方法降低了原始Sleef向量三角函數中冗余計算的開銷，同時保證了函數計算精度.在飛騰處理器平臺上的實驗表明，本文所提方法有效提高了向量三角函數性能.未來我們將進一步研究向量指數函數、對數函數、冪函數的性能優化以及代碼的自動向量化工作，使得更多的科學計算應用可以得益于超越函數向量化而提升性能.