電漿火炬電源非線性且隨機控制的研究

(1)

行政院原子能委員會

委託研究計畫研究報告

電漿火炬電源非線性且隨機控制的研究

A Study of Nonlinear and Stochastic Behavior in Plasma Torch

計畫編號：972001INER019 受委託機關(構)：國立交通大學電機與控制工程學系計畫主持人：廖德誠教授核研所參與人員：李恆毅、黃世文聯絡電話：(03) 5712121 ext. 54363 E-mail address：[email protected] 報告日期：中華民國 97 年 12 月 23 日

(2)

目錄 目錄... I 中文摘要... 1 英文摘要... 2 壹、計畫緣起與目的... 3 一、計畫緣起... 3 二、計畫目的... 5 三、研究項目與進度... 7 貳、研究方法與過程... 8 參、主要發現與結論... 13 一、直流電漿火炬動態模型之建立...13 (一) 雙對流理論 ...14 (二) 直流電漿火炬之數學模型 ...18 二、模擬及分析電漿火炬非線性之動態特性...33 三、協助建立井式電漿火炬混沌參數量測儀器之實驗平台...68 四、井式電漿火炬混沌參數之初步量測與資料分析...75 五、電漿火炬混沌現象成因之分析與探討...91 六、混沌信號判斷方法及控制對策...92 七、結論...93 肆、參考文獻... 96

(3)

中文摘要 本計畫擬配合原能會核研所之研究，探討井式電漿火炬之非線性特性，並驗證理論分析成果之成效。在近年來已發表文獻中指出，大氣電弧電漿裝置之顫動現象可能為混沌動態行為。此研究結果指出，電漿火炬之顫動現象可能為混沌行為而非隨機行為，而這樣的動態行為是有可能被控制的。由過去的研究發現，混沌動態行為有可能是因為系統不穩定周期性行為所衍生的，如果這些不穩定周期性行為發生的成因能夠被探討並加以控制且使其穩定化後，將可消除電漿裝置的顫動現象以改善系統之性能。本計畫應用分叉理論及非線性控制理論針對井式電漿火炬之非線性特性加以探討，找出其造成系統混沌現象的成因，並探討混沌現象發生的可能條件。此外，本計劃也協助核研所之井式電漿火炬實驗平台的建立與量測資料分析以初步確認理論分析之成果。

(4)

Abstract

It is known that a DC power-driven plasma torch might exhibit random behavior which might affect the performance of the torch. A recent study shows that such a random behavior phenomenon might be a deterministic type chaotic behavior. One of the major goals of this project is to verify such a finding. In the project, we have cooperated with the researchers from the Institute of Nuclear Energy to work on the dynamical study of the well-type DC power plasma torch. First, an analytical model was derived from theoretical point of view. It was followed by the dynamical analysis of the derived system model via the help from general nonlinear system theory. Numerical studies of the plasma torch’s behavior were also carried on by using the code Matlab and AUTO. Moreover, we have worked with the researchers from the Institute of Nuclear Energy to build up a test-bed for well-type DC power plasma torch for experimental study. A preliminary of sensor measurement, data collection and numerical analysis of the experimental data were also covered in the proposed tasks.

(5)

壹、計畫緣起與目的 隨著科技的日新月異，人類的生活因科技的進步而愈來愈便 利。然而，科技產生的污染與廢棄物未因時代的進步而減少，反而造成許多無法處理的污染源與垃圾。目前，一般的廢棄物處理最常使用高溫燃燒與垃圾掩埋的方式，利用高溫燃燒將垃圾體積縮小以便於提高垃圾掩埋場的使用年限。但是高溫燃燒會產生有毒氣體與灰渣，且垃圾掩埋場會造成人民的抗爭與場地的限制，甚至可能造成土地二次污染。高溫電漿火炬的開發即是解決廢棄物處理的有效方法之一。因此，為了解決廢棄物處理問題，本計畫配合國家「電漿焚化熔融處理有害廢棄物產業化應用與發展」施政目標，探討高功率直流電漿火炬的動態行為且分析高溫電漿火炬目前尚須改善的問題，以提升國家環境保護技術，朝零廢棄物的目標發展。一、計畫緣起在地球上，物質大部分是呈現固態、液態和氣態。當物質為氣態時，電子在電場束縛下圍繞原子核旋轉，若氣體被加熱，其電子的熱運動動能就會增加。持續加熱至電子的熱運動動能超過原子核對它的束縛，電子就成為自由電子，這種過程稱之為電離。若氣體之電離度為 100％時，即氣體被完全電離，則氣體就會游離成電漿。電漿依溫度與用途可分為熱電漿（Thermal Plasma）和冷電漿（Cold Plasma），本計畫以熱電漿系統為研究重點。電漿火炬則為熱電漿源的一種；電漿火炬中心的溫度高達攝氏一萬度以上，其熱輻射可使熱傳效率優於傳統火焰，大

(6)

部分所有的物質，在這種高溫之下會裂解成簡單的原子，甚至難以破壞的碳氫化合物，包括多氯聯苯，也會被分解成簡單無害的氣體。因此，處理灰渣熔融，高濃度廢溶劑、低放射性廢棄物及低濃度有害廢氣時，電漿火炬也都可以有效的處理且安全固化焚化爐所產生之高污染灰渣。近年來，工業上已廣泛的使用電漿裝置，例如：處理廢棄物、電漿噴塗、半導體製成等等。然而，不同用途的電漿裝置內部顫動（Fluctuation）現象會影響系統的使用效能與使用週期。雖然已有許多相關的研究文獻探討電漿火炬的顫動現象，但仍侷限於實驗觀察與分析，對於不同電漿裝置內部顫動現象的分析與了解仍有待深入探討。在過去的研究，分析與探討電漿火炬之動態行為仍以實驗量測為主，改善顫動現象的效果有限且必須花費相當多的設備成本才能量測到狀態之數據。因此，學者開始朝數學理論發展，建立電漿火炬之數學模型，希望藉由數學理論分析判斷與估測影響電漿火炬使用效率的原因。研究學者依據電漿火炬之物理特性開始建立系統之數學動

態模型，透過能量守恆理論推導出 PDE (Partial Differential

Equation) 形式之數學模型。然而，PDE 形式之數學模型並不易求解與分析，利用數值方法求解也必須仰賴高速電腦進行平行

運算才能縮短計算時間與誤差。在 2000-2004 年，印度學者

Ghorui、Sahasrabudhe、Murthy, Das 和 Venkatramani 則引用三對流(Triple Convection)理論推導電漿火炬之數學模型，其目的是將在空間變化的系統狀態投影至時間變化的非線性三階振幅方程式。如此一來，即可引用發展成熟的系統分析理論探討電漿

(7)

火炬之動態行為且擬定電漿火炬之性能改善對策。另外，學者 Ghorui, Sahasrabudhe, Murthy, Das 和 Venkatramani 透過三階振幅方程式與實驗量測發現電漿裝置內部顫動現象可能是一種混沌

動態行為（Chaotic dynamical behavior）。學者提出顫動現象並非

以往認知的隨機行為而是會呈現週期性變化的混沌現象，電漿裝置內部的顫動行為就有可能可被控制以改善其特性。但從非線性理論分析三階非線性振幅方程式之結果發現，仍有其他因素可能使電漿火炬進入混沌行為。因此，顫動現象之成因須再深入探討，並分析使電漿火炬發生顫動現象的各種可能發生之原因。這些研究結果將可提供更完善的理論依據以避免系統發生不可預期的結果。二、計畫目的近年來，電漿火炬已被廣泛的應用於各種不同的工業用途，例如電漿噴塗、線電弧噴塗、電漿化學蒸汽沈積（CVD）、電漿粉末合成、有毒廢棄物分解、珠化處理、電漿熔融、電漿精鍊及奈米製成等。然而，受限於電漿火炬系統本身之硬體架構，系統仍有許多問題尚待解決且電漿火炬的動態分析仍以實驗量測為主，實驗流程往往需花費相當多的時間與成本。另外，高溫電漿火炬的內部動態為三度空間中的流體變化，其物理與化學反應相當複雜，且火炬的內部與外部溫度相當高(約可產生 5,000～20,000℃)，實驗數據不易量測。因此，透過實驗量測往往無法有效分析其行為特性。基於上述考量，近年來之相關研究文獻開始建立數學分析模型本以探討電漿火炬之動態特性。在數學理論方面，研究文獻建立偏微分方程式(Partial Differential

(8)

Equation, PDE)形式之數學分析模型以探討電弧電漿之動態行為。然而，PDE 之動態模型不易求解，數值分析往往需要透過高速電腦進行平行運算且邊界條件、程式撰寫及相關參數設定會影響準確度。在2000-2004 年，學者 Ghorui、Sahasrabudhe、Murthy、Das 和Venkatramani 則是將 PDE 形式之動態方程式投影至常微分方

程式(Ordinary Differential Equation, ODE)。研究學者將空間變化之系統狀態投影至時間變化之振幅方程式，動態行為之特性探討就可引用系統分析理論解釋電漿系統之顫動行為。其振幅方程式表示如下： 3 0 1 2F F F F

F&&&+Ω &&+Ω & +Ω =−

透過文獻之模擬分析與實驗驗證發現電弧電漿裝置弧根的顫動現象可能是一種混沌動態行為（Chaotic dynamical behavior）而非以往認知的隨機行為，則電漿之顫動行為可能可以經由控制來改善其特性。在模擬分析方面，學者Ghorui、Sahasrabudhe、 Murthy、Das 和 Venkatramani 則是藉由改變系統參數Ω₀發現當系統參數進入某一個範圍內會發生混沌現象。然而，從文獻建立的三階非線性振幅方程式發現混沌現象的發生可能不只Ω₀改變，即Ω₁與Ω₂的改變也可能使系統發生混沌行為。因此，造成混沌現象的參數範圍仍需再深入分析以找出使電漿系統發生顫動現象的原因。本計畫之主要目的是建立可量測混沌參數之非線性混沌分析模型，並且分析其熱電漿火炬之動態特性。另外，協助架設井式電漿火炬參數量測之實驗平台，並且透過實驗初步量測火

(9)

炬之數據。引用系統分析理論初步分析量測之數據及分析電漿火炬之顫動現象。最後，提出混沌行為之判斷方法與控制對策以提供日後相關研究者之參考依據。三、研究項目與進度本計畫之研究工作主要在於分析直流電漿火炬之顫動現象，其研究項目如下所示： 1. 建立電漿火炬可量測混沌參數﹙包括電壓、電流、聲音、光線強度等﹚之非線性混沌分析模型。 2. 根據上述建立之模型，模擬及分析電漿火炬非線性之動態特性。 3. 協助準備及建立井式電漿火炬混沌參數量測儀器之實驗平台。 4. 井式電漿火炬混沌參數之初步量測與資料分析。 5. 電漿火炬混沌現象及成因之分析與探討。 6. 提出混沌信號判斷方法及控制策略。

(10)

貳、研究方法與過程 本計畫的主要目的是藉由數學理論分析直流電漿火炬顫動行為的成因，也探討直流電漿火炬可能發生的動態特性。在過去，電漿系統內部的顫動現象被認定為隨機的動態行為而無法估測其系統狀態。然而，近年有研究文獻提出電弧電漿裝置之顫動現象可能為混沌行為而非以往認知的隨機行為，則系統之顫動現象可能可以經由控制改善火炬的性能。本計畫利用非線性理論與分叉理論探討電漿火炬之顫動現象，分析可能使系統發生混沌現象之分叉參數的變化以及分叉參數的變化範圍。以下說明此計畫之研究方法與過程： 1. 建立電漿火炬可量測混沌參數之非線性混沌分析模型為了能分析電漿火炬之動態行為，必須先建構電漿火炬之數學模型，以方便對直流電漿火炬系統之模擬及分析。然而，電漿火炬之動態會隨著空間和時間變動，經由電磁理論與守恆定理推導之數學模型為PDE 形式。PDE 形式之動態方程式不易求解，且電漿火炬之 PDE 數學模型可能無法直接引用系統分析理論與控制理論。若透過數值方法求PDE 之解，則必須依賴高速平行電腦計算且可能會造成計算時間過長與誤差過大，將 PDE 形式之數學方程式投影至 ODE 形式數學動態模型是必要的。因此，本計畫依據P. H. Coullet 和 E. A. Spiegel 提出之 normal form 理論、雙對流理論以及系統物理特性，透過數學推導以建

構系統分析模式。經由normal form 理論與雙對流理論可以建構

(11)

2. 依據建立之模型，模擬及分析電漿火炬非線性之動態特性此部分主要是針對建立之數學分析模型探討直流電漿火炬之動態行為，探討系統之動態特性。實際系統在運作時，往往受限於本身的物理特性與硬體架構而無法完全理想化操作。當系統操作在極限邊緣時，系統之動態可能會因參數稍微變動而使系統由穩定狀態改變為不穩定或發生振盪週期。探討系統參數之變化與系統動態行為之間的關聯性可以避免系統改變穩定性。因此，本計畫應用非線性理論及分叉理論以推導系統之靜態與動態工作解及其成立之系統參數解析條件。另外，我們應用非線性系統模擬軟體，例如：Matlab 及 Auto，以驗證理論推導之成果及 non-local 系統靜態與動態特性。對系統的分叉現象分析，我們採用Auto 這套數學軟體來繪製系統之分叉圖，分析系統可能發生的動態特性。另外，研究文獻提出電弧電漿裝置之顫動現象可能為混沌行為，但文獻中僅考慮單一參數變化之動態分析；本計畫則考慮系統可能發生混沌現象之原因，且找出改變系統穩定性的參數範圍。我們也探討學者 A. K. Das 沒有分析的系統參數以及模擬可能之動態特性；從分叉圖可以知道系統動態發生改變的範圍。接著，此部分選取分叉圖上每個區域之分叉參數繪製時間響應圖和相位圖，模擬系統之動態特性以及震盪週期之變化。另外，為了探討發生混沌現象之成因，我們透過雙參數分叉理論探討二個參數同時發生變化之系統動態特性。我們也選取雙參數分叉圖上每個區域之分叉參數以及選取系統初始值繪製時間響應圖和相位圖，驗證系統之動態特性以及振盪週期之

(12)

變化。 3. 協助準備及建立井式電漿火炬混沌參數量測儀器之實驗平台除了理論分析工作之外，在本計畫中，我們透過與原子能委員會核能研究所物理組同仁之合作，協助架構相關的實驗設備，並以實際系統之實驗量測數據來評估理論分析之成果。茲簡述如下： (1) 本計畫透過與原能會核研所物理組討論直流電漿火炬實驗平台之操作範圍與量測儀器之使用規格和限制，規劃直流電漿火炬的相關實驗平台。 (2) 協助建立一適當的量測系統在此部分之工作重點在於建立實驗平台，實驗設備包含直流電漿火炬平台以及量測系統，例如：架設紅外線測溫儀觀察火炬的溫度變化、架設電壓電流探棒、架構流量計量測直流電漿火炬內部的流量變化、架設麥克風接收直流電漿火炬的音波以及設計ADC將擷取資料存入電腦以提供日後分析。量測儀器之架構方面，本計畫主要是針對電漿火炬之顫動現象探討，架設之重點放在電壓、電流、光和聲音之訊號量測，藉由實驗量測之數據以利於系統之模擬驗證。 (3) 評估量測相關的周邊系統之操作特性以及確認各項儀器之使用規格，例如：最大氣體流量、磁場、音波量測、溫度量測之範圍、外加電流或外加電壓。確認量測平台之操作範圍以及儀器規格可以提高實驗數據之準確性，避免擷取到錯誤資料。

(13)

4. 井式電漿火炬混沌參數之初步量測與資料分析此部分主要是經由實驗量測數據初步分析直流電漿火炬發生的顫動現象，找出影響直流電漿火炬產生顫動的主因，以便於日後研擬控制對策與設計。在實驗量測方面，本計畫擷取電壓、電流、光和聲音之實驗數據，透過時間響應和頻譜分析系統之動態特性。從實驗結果可以觀察到電壓、電流、光和聲音均有顫動行為的發生，頻譜分析也顯示系統會發生週期變化。 5. 電漿火炬混沌現象及成因之分析與探討此部分之工作重點是經由數學理論分析可以評估直流電漿火炬可能發生混沌現象的成因。本計畫希望透過非線性理論和分叉理論探討系統可能發生混沌行為的分叉參數變化以及參數範圍。因此，我們建構系統分叉圖以估測造成混沌行為的可能原因，並模擬分析電漿火炬之時間響應圖和相位圖觀察系統的動態特性。另外，除了電漿火炬之混沌行為，本計畫探討系統可能產生的其他動態特性，例如：電漿火炬之動態行為從穩定狀態改變至振盪週期的特性分析，最後進入混沌現象到系統狀態發散之參 數變化。 6. 提出混沌信號判斷方法及控制策略已知高溫電漿火炬之物理特性相當複雜，系統狀態並不易量測且數據擷取受限於量測儀器之使用規格而無法精確量測實驗數據。然而，學者A. K. Das 提出顫動現象可能為一混沌行為，那電漿火炬可能可以控制改善其性能。從數學理論分析結果可以預估混沌現象的系統參數範圍；當分叉參數值開始減少且選

(14)

取特定區域的系統初始值會導致系統由穩定狀態轉變為一倍週期之振盪行為，在由一倍週期依序轉換成二倍、四倍和八倍週期，最後進入混沌現象到系統狀態發散；由數學理論估測分叉參數之變化範圍，監控系統動態之變化以避免進入混沌行為。另外，混沌現象是一個奇特的動態行為，系統進入混沌行為時，狀態會呈現有界且不規格之週期變化；這種動態特性會造成無法預估系統之行為，即無法控制。因此，此部分之研究方法必須依據上述之數學理論分析發生混沌行為的參數範圍，設計監控系統與控制器補償系統狀態以避免產生不規則之週期變化。

(15)

參、主要發現與結論 本計畫之研究方法與過程主要是針對直流電漿火炬建立數學分析模型，並透過非線性理論與分叉理論分析直流電漿火炬之顫動現象。另外，本計畫也考慮可能使系統發生混沌現象的其他因素，探討發生混沌行為的各種可能性，並與學者 A. K. Das 之研究結果做比較。依據非線性系統分析之方法，搭配分叉理論可以發現其他參數也發生變化時，其他的系統參數仍會因變化而產生分叉現象，即系統可能發生混沌的動態行為。因此，本計畫之研究成果提供其他可能發生混沌行為之系統參數範圍與其他可能因素。在實驗方面，我們依據核研所現有的直流電漿火炬系統協助建立實驗之量測硬體架構，且針對系統之數據量測提供相關意見與建議。另外，透過核研所物理組安排實驗並初步量測相關數據，本計畫引用系統分析理論初步分析量測數據與探討火炬之動態行為。最後，藉由理論探討與實驗初步量測分析之結果提出直流電漿火炬發生顫動現象的可能原因與初步研擬控制對策，以提供未來火炬效能改善之依據。以下將逐一說明本計畫之研究成果。一、直流電漿火炬動態模型之建立直流電漿火炬系統動態模型是經由能量守恆、金屬守恆、動量守恆、質量守恆以及 Maxwell 理論所建立，並搭配符合實際物理意義設定邊界條件。建立之數學模型是考慮在空間中變化，即數學分析模型是以偏微分的形式呈現。PDE 形式之動態模型不易求解且無法有效擬定改善對策。因此，本計畫依據 P. H.

(16)

Coullet 和 E. A. Spiegel 以及 A. K. Das 提出之理論推導以時間變化為基礎之非線性振幅方程式。另外，A. K. Das 提出的研究成果僅討論單參數變化造成混沌現在之參數範圍，而本計畫利用非線性理論與分叉理論分析直流電漿火炬數學模型，探討其他可能發生變化之系統參數與混沌動態行為之間的關聯性。最後，分析模擬之結果與A. K. Das 提出之成果作比較。 (一) 雙對流理論雙對流理論是P. H. Coullet 和 E. A. Spiegel 於 1983 年所提出，其構想是建立以時間變化為基礎之動態方程式，以利於引用系統分析理論探討其動態特性。雙對流理論是將數學模型正規化，經由數學投影簡化成正規形式。因此，本計畫利用正規化理論將PDE 形式之直流電漿火炬數學模型簡化為 ODE 形式之非線性振幅方程式。首先，考慮一非線性方程式 ) (U N U M U t = λ + λ ∂ (1) 其中 t 表示時間、λ=(λ1,λ2,K,λp)表示 p 個參數的集合、 ) , , (U₁ Up U = _K 則表示定義在實函數空間的系統狀態以及Mλ 和N_λ分別表示線性和非線性操作子(Operator)且N_λ(U)至少為二次式。當M_λ和Nλ包含空間微分，我們可以將邊界條件指定為 0 ) ( = + λ λU C U B on ∂V (2) 其中 V 表示物理空間(Physical Space)；B_λ和Cλ分別表示線性與非線性操作子。假設M_λ的頻譜為離散的。為了簡化表示，

(17)

我們僅考慮在非退化之頻譜。首先，假設此問題為同次的（Homogeneous），則可以得到下列形式的解 k t e U = ηΦ (3) 其中 ) , ( ) , ( ) , ( 2 2 1 1 k x W k x W k x W N k N k k k Ω Ω Ω = Φ M (4) k L Ω 為常數且WL為位置函數。當(1)為 ODE 時，則WL =1， N L=1,2,_K 。向量 k 表示用來描繪空間架構函數WL的參數集 合。對於每一個 k，(1)和(2)可以簡化為特徵值問題 Mλ|k Φk =ηΦk (5) 和 Bλ|k Φk =0 on ∂V (6) 特徵方程式則可以為下列之形式 ( ; ) ( ) 1 ₀( ) 0 1 λ η + + λ = + η ≡ λ η − − k N k N N N k a a P _L (7) 對於(7)的每一個根，我們可以找到振幅 k i Ω 以及相符合的正規形式Φk。因此，將(7)整理後可以得到 (η;λ)≡

∏

(η;λ)=0 k N k P P (8) 假設(8)有 l 個零根(η=0)和 m 對虛根(Reη=0和Imη≠0)，則(8) 有 d =l+2m (9) 個Reη=0的根。這些 d 個根則稱為 critical。假設在λ₀附近的Δ 區域內，系統存在有很小的|Reη|及系統的其他根(Reη≤η₀ <0) 也在Δ區域內，其中η₀為實常數；則ODE 形式的(1)之解則可

(18)

能使η₀ →−∞。在(8)中我們可以在Δ區域找到二種可能性的多項式：一種是degree 為 d 且根為 critical 的多項式(λ=λ₀)，此多項式稱為critical polynomial；另外一個多項式則是根在Δ區域從虛軸遠離。因此，若考慮最簡單的例子：l =1和m=0， critical polynomial 則為η+μ₀k(λ)_{，在此情況會發生靜態分叉} (Stationary Bifurcation)；另外，當l =0和m=1(即d =2)，則會

發生 Hopf bifurcation ，其中 critical polynomial 為

) ( ) ( 2 1 2 ₊_μ _λ _η₊_ω _λ η k 。考慮d =2的情況下，定義 codimension 2 表面之條件 μ1(λ)=0 k 和 0 ) ( 0 λ = μk (10) 考慮一般式，存在有n=l+m個critical 條件定義 codimension n 的表面，則將此表面稱為polycritical surface。當我們可以將問題簡化為ODE 時，求得的解(隨時間變化)會接近 polycritical surface。接著，我們將依據 polycritical surface 簡化問題。假設(1) 存在有polycritical surface [μi(λ);ωj(λ);η]=0 C P , i=0,1,_K,n−1, j =0,1,_K,m−1 (11) (11)有 l 個η=0的根和 m 對η=±iω_j 的根。因此，系統的 polycritical surface 有 d =l+2m 個 critical modes 且有

codimension n=l+m。假設有任意點λ=λ₀存在於 polycritical

surface ，從 (3) 可得到 generalized critical modes φi ,

d

i=1 ,2,_K, ；modes 的集合構成了 generalized critical surface 的

基底(basis)。這些集合從(3)提供了完整的 stable modes fj。因 此，U 可以分成二部分：一個是 generalized critical space，另

(19)

外一個為stable space。我們將 U 改寫成以下之式子：

∑

= β + φ α = d i j j j i i t x t f x t x U 1 ) ( ) ( ) ( ) ( ) , ( (12) (1)則可以寫成 α& =Jα+F(α,β) (13) β&=Lβ+G(α,β) (14) 其中 J 和 L 為矩陣且 F 和 G 為非線性函數。非奇異矩陣 (nonsingular matrix)L 的所有特徵值均在負實部，而 J 有零實部 且 J 為 Jordan form。接著，引進非線性轉換 β=b+B(α), ∂_αB|_α₌₀=0 (15) 將(15)帶入(14)，則可以得到 L b LB (J ) B F( ,b B) B G( ,b B) dt d ₌ ₋ _α _∂ ₋ _α ₊ _∂ ₊ _α ₊ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ₋ α α (16) 若可以選擇 B 使得 (Jα∂α−L)B=−F(α,B)∂αB+G(α,B) (17) 則b=0定義為(13)和(14)的不變流形(invariant manifold)，且對於很小的|α|，此不變流形是線性穩定的。將 β B= (α) (18) 帶入(13)，則可以得到 α& =Jα+F(α,B(α)) (19) 因此，對於μ≡(μ₀,μ₁,K,μn₋₁)=0，我們可以知道系統有 n 個振幅方程式(amplitude equations)存在使得(17)可以求解的，即可以利用perturbation 理論求解。當根離開 polycritical surface 時，(19)可以變形為 α& =K_μα+h_μ(α) (20)

(20)

其中K_μ表示當 J 的所有 elements 在擾動以及h_μ則表示為非線性向量(vector-valued)函數，即K₀ =J ，K_μ 的特徵方程式為 critical polynomial C P 。 C P 有 d 個參數；一般來說，K_μ有d2個參數。我們要將K_μ放入正規形式且包含最少的參數。因此，可以利用以下非奇異之座標轉換得到上述結果： A=C−1αC (21) 另外，J 的變形則稱為 Jordan-Arnold form，符號定義為J_μ。將 μ K 放入Jordan-Arnold form，則(20)可以改寫為 A& = JμA+gμ(A) (22) 其中 J J₀ = (23) 表示為另外一個非線性向量函數。因此，(22)即為振幅方程式的正規化形式。接下來，本計畫將依據此雙對流理論推導直流電漿火炬之數學動態模型。 (二) 直流電漿火炬之數學模型直流電漿火炬之建立主要為了探討弧根之動態行為，建立能解釋顫動行為之數學動態模型。直流電漿火炬之動態行為是由熱守恆(conservation of heat) 、金屬守恆 (conservation of metal) 、動量守恆 (conservation of momentum) 、質量守恆 (conservation of mass)及 Maxwell’s equation 所決定且依據電漿

火炬之物裡特性找出電漿系統之boundary conditions。以下將

依據P. H. Coullet 和 E. A. Spiegel 與 A. K. Das 之研究理論為基礎，建立直流電漿火炬之非線性方程式。建立電漿火炬，我們

(21)

有以下六個假設: 1. 弧根可近似為對稱的圓柱體 2. (u, v, w)表示在r、θ和 z 方向之速度 3. T(r,θ,z,t)和S(r,θ,z,t)分別表示產生顫動現象時的溫度和金屬揮發濃度及分別存在有靜態組合T₀(r,θ,z)和S₀(r,θ,z) 4. z 方向且高度為 d 的區域內存在溫度梯度ΔT d和金屬濃度梯度ΔS d 5. 在 r 方向存在均勻溫度 6. 對於圓柱形的電弧圓柱體的半徑 R 而言，系統操作在固定 電流的電源供應且電流密度J₀趨近於常數描述電漿火炬顫動現象的動態方程式分別可由下述之方程式共同決定： (1) 熱守恆 T T v d T w t T ₋ Δ ₊ _⋅_∇ ₌_κ_∇2 ∂ ∂ _r r (24) 其中

T：fluctuating temperature component

s k ρ = κ ：熱擴散率 z w v r u vr= ˆ+ θˆ+ ˆ：速度 k：熱傳導率 ρ：密度 s：比熱 (2) 金屬守恆 S S v d S w t S S 2 ∇ κ = ∇ ⋅ + Δ − ∂ ∂ _r r (25)

(22)

其中 S：金屬濃度 κS：材料的傳導率 (3) 動量守恆 v v p v v g T g S J B t vr _r _r _r _r _r r_× r ρ + β − α + ∇ + ∇ ρ − = ∇ ⋅ + ∂ ∂ 0 2 0 1 1 ₍₂₆₎ 其中 壓力(p)造成的垂直應力 黏性(v)造成的切線方向應力 流體和材料加入的擴張所造成的引力電流密度(J)伴隨著磁場(B)的相互影響所造成的 Lorenz force gr：重力加速度 (4) 質量守恆 ∇ vrr ⋅ =0 (27) (5) Maxwell’s equation 1 ( ) 0 B Jr ∇r × r μ = , t B E ∂ ∂ − = × ∇ r , ∇ Br⋅ r=0, ε ρ = ⋅ ∇ C Er (28) 其中 σC：導電率 μ₀：介電係數 (6) 電流密度 J C(E v B) r r r r × + σ = (29) (7) 狀態方程式

(23)

Δρ=ρ₀(−αT +βS) (30) 其中 α和β則定義於(30) E：電場 k：熱傳導率 ρ_C：電荷密度 考慮 B 為圓柱型對稱。因此(28)可以寫成 v B B B v t B B r r r r r r r r ) ( 2 ₊ _⋅_∇ ∇ η = ∇ ⋅ + ∂ ∂ 若考慮接近弧根處為一微小圓柱，則可近似為Br =B_θ(r)θˆ和 z J Jv= ₀ˆ。方程式則可以改寫為 B u J v B t B B 2 0 0 ) , ( −μ +η ∇ φ = ∂ ∂ (31) 其中 φ(B,v) =(B⋅∇)v+u(B r)， C B σ μ = η 0 1 _， θ = B B r：圓柱的半徑 為了在簡化上述之動態方程式，取其旋度移除壓力項。引進 stream function ψ使得 z dt dr u ∂ ψ ∂ = = ， r dt dz w ∂ ψ ∂ − = = ，Γ(ψ,ϕ)=−Vr⋅∇rϕ ψ ∇ = ∂ ∂ − ∂ ∂ = η 2 r w z u ， 2 0 I d J = ，B₀ =μ₀J₀R 2 接著，利用表一 dimensionless quantities 的形式和表二的 nondimensional numbers 帶入上式，可以將方程式簡化。

(24)

表一、Form of dimensionless quantities Dimensionless quantities Form

Temperature T T T′= _Δ Concentration gradient S S S′= _Δ Spatial coordinates d x x′= , d y y′= , d z z′= r d r′= ⋅ , θ′ d= ⋅θ, z′=d⋅z Velocity components _u u_d κ = ′ , d v v′= κ , d w w′= κ Time t′= _d2t_κ Magnetic field B0 B B′= 表二、Nondimensional numbers Rayleigh number Rl R_l g _kvTd 3 Δ α =

Absolute Raylleigh number Rs R_s g _kvSd

3 Δ α = Prandtl Number σ k v = σ

Magnetic Prandtle number ξ ξ= η_kB

Lewis number τ _τ₌ κ_κS Chandrasekhar number Q Bd_v Q= μρ0η 2 2 0 其中定義r=dr′、θ d= θ′、z=dz′。在nondimensionalization 之後，governing equation 的集合可以改寫如下：由(26)知： v v p v v g T g S J B t v r r r r r r_× r ρ + β − α + ∇ + ∇ ρ − = ∇ ⋅ + ∂ ∂ 0 2 0 1 1

(25)

(A1) t v ∂ ∂ _⇒_∵ ₌_η ∂ ∂ − ∂ ∂ = r w z u v curl ⇒ ⎟=η ⎠ ⎞ ⎜ ⎝ ⎛ ′ ∂ ′ ∂ − ′ ∂ ′ ∂ ⋅ κ r w z u d2 ⇒ η= κ₂ η′ d ∴ t d t v ′ ∂ η′ ∂ κ = ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ∂ 4 2 curl (A2) vr ∇⋅ vr ⇒ 取旋度，則 curl( vr⋅∇vr)=vr⋅∇η ⇒ ₄ ( , ) 2 4 2 η′ ψ Γ κ − = η′ ∇ ⋅ ′ κ = ∇ ⋅ d v d v vr r r (A3) ν⋅∇2vr ⇒ 取旋度，則 ⇒ ∇ = κ ν⋅∇ η′= κ ⋅κσ⋅∇ η′= κ σ⋅∇2η′ 4 2 2 4 2 4 3 2 ₎ (ν curl d d d vr (A4) Jr×Br ρ₀ 1 (A5) grαT =ΔTgαT′ ⇒ Tg d T d T Tg T Tg r ⋅∂r ′ κ ν κν α Δ ⋅ κ − = ′ −∂ α Δ = ′ α Δ ′ 3 4 2 ) ( ) ( curl R T d lσ⋅∂r ′ κ − = 2₄ _′ (A6) gβr S，同(A5)可得 ⇒ Sg d S d S Sg S Sg r S S r ∂ ′ κ κ ⋅ κ ν ν κ β Δ ⋅ κ = ′ ∂ β Δ = ′ β Δ ′ 3 4 2 ) ( curl R T d sστ∂r ′ κ = ₄2 _′ 由(A1)-(A6)可知 R T d d d t d lσ⋅∂r ′ κ − η′ ∇ ⋅ σ κ + η′ ψ Γ κ = ′ ∂ η′ ∂ κ ′ 4 2 2 4 2 4 2 4 2 ) , ( Q B Rd T R d s r′ σξ ∂z′ κ − ′ ∂ στ κ + ₄2 2 2₄

(26)

⇒ Q B R S R T R t′ =Γ ψ η′ +σ∇ η′− lσ∂r ′+ Sστ∂r ′− σξ ∂z ′ ∂ η′ ∂ ′ ′ ′ 2 ) , ( 2 ₍₃₂₎ 由(24)知： v T T d T w t T ₋ Δ ₊ _⋅_∇ ₌_κ_∇2 ∂ ∂ _r r (B1) t T k d T t T ′ ∂ ′ ∂ ⋅ Δ = ∂ ∂ 2 (B2) ₂ ⎟= κΔ₂ (−∂ ψ) ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ψ ∂ − Δ κ = Δ ⋅ ′ ⋅ κ = Δ ⋅ r d T r d T d T w d d T w (B3) ₂ ₂ ( ,T ) d T T v d T T vr⋅∇ = κΔ r′⋅∇ ′=−κΔ Γ ψ ′ ∵vr=urˆ+vθˆ +wzˆ ⇒ v d z w v r u d vr = κ( ′ˆ+ ′θˆ + ′ˆ)= κ r′ ∵ z T y T x T T ∂ ∂ + ∂ ∂ + ∂ ∂ = ∇ ⇒ T d T z T y T x T d T T _⎟⎟= Δ ∇ ′ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ′ ∂ ′ ∂ + ′ ∂ ′ ∂ + ′ ∂ ′ ∂ Δ = ∇ ∴ ₂ ₂ ( ,T ) d T T v d T T vr⋅∇ = κΔ r′⋅∇ ′=−κΔ Γ ψ ′ (B4) T d T T = κΔ ⋅∇ ′ ∇ κ 2 2 2 ∵ 2 2₂ ₂2 ₂2 z T y T x T T ∂ ∂ + ∂ ∂ + ∂ ∂ = ∇ ⇒ T d T z T y T x T d T T _⎟⎟= Δ ∇ ′ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ′ ∂ ′ ∂ + ′ ∂ ′ ∂ + ′ ∂ ′ ∂ Δ = ∇ 2 2 2 2 2 2 2 2 2 2 因此，從(B1)-(B4)可知 T d T d T T d T t T d T r ⋅∇ ′ Δ κ + ψ ∂ ⋅ Δ κ − ′ ψ Γ ⋅ Δ κ = ′ ∂ ′ ∂ ⋅ Δ κ 2 2 2 2 2 ( , ) ⇒ T T t T rψ+∇ ′ ∂ − ′ ψ Γ = ′ ∂ ′ ∂ ₍ _, ₎ 2 ₍₃₃₎ 由(25)知，同理可得

(27)

S S t S rψ+τ∇ ′ ∂ − ′ ψ Γ = ′ ∂ ′ ∂ ₍ _, ₎ 2 ₍₃₄₎ 其中κS =τ⋅κ。由(31)知： B v J u B t B B 2 0 0 ) , ( −μ +η ∇ φ = ∂ ∂ (C1) t B d B t B ′ ∂ ′ ∂ κ = ∂ ∂ 2 0 (C2) ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ′ κ + ′ ∇ ⋅ ′ κ = ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + ∇ ⋅ = φ r B d B v B d B r B u v B v B, ) ( ) 0₂ ( ) 0 ( (C3) u dR B u J = κ ′ μ 0 0 0 2 ∵ 2 0 0 0 R J B =μ ⇒ R B J₀ 0 0 2 = μ (C4) B d kB B B ∇ ′ ξ = ∇ η 2 2 0 2 從(C1)-(C4)可知： B d kB u dR B r B d B v B d B t B d B ₋ κ _′₊ ξ_∇ _′ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ′ κ + ′ ∇ ⋅′ κ = ∂ ′ ∂ κ 2 2 0 0 0 2 0 2 0 ₍ ₎ 2 假設r′=d⋅r且ξ=κ ηB ，則 B dR B v t B zψ+ξ∇ ′ ∂ − ′ ′ φ = ′ ∂ ′ ∂ ′ −1 2 2 ) , ( (35) 接下來，我們將利用正規化形式推導直流電漿火炬之振幅方程式。首先，將(32)-(35)改寫為下列形式 ∂tLY =MλY+N(Y) (36) 其中 L 為非奇異非線性操作子，M 表示為相依操作參數之線 性操作子 N 則表示完全非線性操作子。因此，我們可以得到 下列之關係式

(28)

⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ψ = B S T Y ， ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛∇ = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 2 L ， ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ′ ′ φ ′ ψ Γ ′ ψ Γ η′ ψ Γ = ) , ( ) , ( ) , ( ) , ( ) ( B v S T Y N ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ∇ ξ ∂ − ∇ τ ∂ − ∇ ∂ − ∂ σξ − ∂ στ ∂ σ − ∇ σ = − − λ 2 1 2 2 1 4 0 0 2 0 0 0 0 2 z r r z r s r l dR Q R R R M (37) 其中動態的控制參數為λ=(Rl,Rs,Q,σ,τ,ξ)。取(36)之線性部分，則可以得到 Y M LY t = λ ∂ (38) 依據雙對流之正規化形式之理論，(38)解的形式假設如下： st mn mn e Y Y = *Λ (39) 其中Ymn為4 個常數組成的向量。接著，定義運算符號“*”： _⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ = ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎝ ⎛ bd ac d c b a * 且 ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ π π π π = Λ ) . . cos( ) . . sin( ) . . sin( ) . . cos( ) . . sin( ) . . cos( ) . . sin( ) . . sin( z n r a m z n r a m z n r a m z n r a m mn 其中 m=0,1,2,_K 和n=0,1,2,K (40) 接著，將(40)帶入(39)，則 ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ψ′ = ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ′ ′ ′ ψ′ st mn st mn st mn st mn e z n r a m B e z n r a m S e z n r a m T e z n r a m B S T ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) sin( ) sin( (41) 將上式帶入(38)的左式，則

(29)

⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ψ′ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛∇ ∂ = ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ′ ′ ′ ψ′ ∂ st mn st mn st mn st mn t t e z n r a m B e z n r a m S e z n r a m T e z n r a m B S T L ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) sin( ) sin( 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 2 ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ∇ ψ′ ∂ = st mn st mn st mn st mn t e z n r a m B e z n r a m S e z n r a m T e z n r a m ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) cos( ) sin( 2 ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ψ′ π − − ∂ = st mn st mn st mn st mn t e z n r a m B e z n r a m S e z n r a m T e z n r a m n a m ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) cos( ) sin( ) ( 2 2 2 2 ⎪ ⎪ ⎭ ⎪ ⎪ ⎬ ⎫ ⎪ ⎪ ⎩ ⎪ ⎪ ⎨ ⎧ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ′ ′ ′ ψ′ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ⋅ π ⋅ ⋅ ⋅ ⋅ π ⋅ ⋅ ⋅ ⋅ π ⋅ ⋅ ⋅ ⋅ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛− − π ∂ = st mn mn mn mn t e B S T z n r a m z n r a m z n r a m z n r a m n a m * ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) cos( ) sin( 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 2 2 2 2 =∂tLmnY 將(41)帶入(38)的右式 MλY ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ′ ⋅ π ⋅ ⋅ ⋅ ψ′ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ∇ ξ ∂ − ∇ τ ∂ − ∇ ∂ − ∂ σξ − ∂ στ ∂ σ − ∇ σ = − − st mn st mn st mn st mn z r r z r s r l e z n r a m B e z n r a m S e z n r a m T e z n r a m dR Q R R R ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) sin( ) sin( 0 0 2 0 0 0 0 2 2 1 2 2 1 4 ⎪ ⎪ ⎩ ⎪ ⎪ ⎨ ⎧ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ π + ξ − π − π + τ − − π + − − σξ π στ − σ π + σ = − − ) ( 0 0 2 0 ) ( 0 0 0 ) ( 2 ) ( 2 2 2 2 1 2 2 2 2 2 2 2 2 1 2 2 2 2 2 n a m n dR n a m ma n a m ma Q R n maR maR n a m _l _s

(30)

st mn mn mn mn e B S T z n r a m z n r a m z n r a m z n r a m ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ′ ′ ′ ψ′ ⋅ ⎪ ⎪ ⎭ ⎪ ⎪ ⎬ ⎫ ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ⋅ π ⋅ ⋅ ⋅ ⋅ π ⋅ ⋅ ⋅ ⋅ π ⋅ ⋅ ⋅ ⋅ π ⋅ ⋅ ⋅ ) cos( ) sin( ) sin( ) cos( ) sin( ) cos( ) sin( ) sin( * mn st mn mn e Y M * ) ( Α = 其中 ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ ξ − π − τ − − − − πξσ στ − σ σ = − − 2 1 2 2 1 4 0 0 2 0 0 0 0 2 mn mn mn s l mn mn q dR n q ma q ma QR n R ma R ma q M (42) ⎟⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎜ ⎜ ⎝ ⎛− = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 2 mn mn q L (43) 2 2 2 2 2 ₌ ₊ _π n a m q_mn (44) 對應給定的特徵值，我們可以得到特徵方程式 det(Mmn −Lmns)=0 (45) 將(42)和(43)帶入(45)，則 ] 0 0 0 2 0 0 0 0 2 det[ 2 1 2 2 1 2 4 = ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ ⎛ − ξ − π − − τ − − − − − πξσ στ − σ + σ − − s q dR n s q ma s q ma QR n R ma R ma sq q mn mn mn s l mn mn ⇒q_mn2 s4 +(1+σ+τ+ξ)q_mn4 s3 +[(σ+τ+ξ+στ+σξ+τξ)q_mn6 −m2a2σR_l 2 2 ₄ 2 2 2_] 2 _{( ₎ 8 mn s n dQR s q R a m στ − π ξσ + στ+σξ+τσξ+τξ + − [ (m a Rl)( ) m a Rs(1 ) 4n QdR (1 )]qmn}s 2 2 2 2 2 2 2 2 _σ _ξ₊_τ ₊ _στ ₊_ξ ₋ _π _ξσ ₊_τ − + − ₊_στξ 10 ₊ 2 2_σξτ(₋ ₊ ) 4 ₋4 2_π2_ξστ −2 4 ₌0 mn mn s l mn m a R R q n dQR q q 因此，我們可以得到特徵方程式

(31)

2 ₁ ₀ 0 2 3 3 4 ₊_Π _s ₊_Π _s ₊_Π _s₊_Π ₌ s (46) 其中 6 2 2 2 2 2 2 0 =[qmn +m a (Rs −Rl)−4n π dQR ]σξτqmn Π − Π = στ+σξ+τξ+τξσ 6 + τ +ξ − ξ+τ 2 2σ 1 ( )qmn [ Rs(1 ) Rl( )]m a ₊_τ (1₊_ξ)₋4 2_π2_ξσ −2(1₊_τ) QdR n Rs 4 2 =(σ+τ+ξ+στ+σξ+τξ)qmn Π ₊_[(_τ ₋ ₎ 2 2 ₋₄ 2_π2_ξ −2_]_σ −2 mn l s R m a n dQR q R 2 3 =(1+σ+τ+ξ)qmn Π 假設s= iω且ω為實數，則代入(46)可以得到 0 0 1 2 2 3 3 4 ₋ _ω _Π ₋_ω _Π ₊ _ω_Π ₊_Π ₌ ω i i (47) 因此，從(47)可以觀察到可能的解 (i) Π₀ =0和ω=0 (ii) _ω2 ₌_Π₁_/_Π₃_和 2 ₀ 3 0 3 2 1 2 1 −Π Π Π +Π Π = Π ，其中Π₁≥0 從雙對流理論可以知道滿足 polycritical surface 的條件為 Π₀ =Π₁ =Π₂ =0 (48) 將條件(48)代入(46)則可以得到 critical 多項式： 2 ₁ ₀ 2 3 ) (s s c s cs c PC ≡ + + + (49) 其中 3 0 0 _Π Π = c ， ₂ 3 0 3 1 1 _Π Π − Π Π = c ， ₃ 3 0 2 3 0 3 2 2 _Π Π + Π Π − Π Π = c (50) ploycritical surface 的條件則為c₀ =c₁ =c₂ =0。因此，若n= m=1 時，可以得到 ) 1 )( 1 ( ) ( 2 6 − ξ − τ σ ξ + τ + σ = a q Rlo ， ) )( 1 ( ) 1 ( 2 6 τ − ξ − τ σ ξ + σ + τ = a q Rso

(32)

) )( 1 ( ) 1 ( 2 6 0 _σ _ξ₋ _ξ₋_τ ξ τ + σ + = c q Q (51) 在 critical 表面之任意點，正規化模式在s=0表面上滿足 0 *Λ = φ _mn C M 其中φ表示四個常數組成之向量。由於在 critical

表面時，det(MC)=0。因此可以得到null space 之基底

T R q d q a q a ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ξ π − τ − − = φ 1 ₂ ₂ 2₂ (52) 其餘的二個向量則可以由以下式子求得 M_Cψ=L_Cφ 和 MCχ=LCψ (53) 最後可以得到線性的動態方程式(38)。因此，由雙對流理論可

知系統可以分成二部分：一是generalized critical space，另一

則是stable modes。我們可以將 Y 寫成下面之方程式

Y =[A(t)φ+B(t)ψ+C(t)χ]*Λmn +S~(r,z,t) (54)

在上式中，S~表示 stable modes；定義符號Φ=[φ ψ χ]為 null

space 的基底且A=[A(t) B(t) C(t)]。從(52)和(53)可以得到

∑

= Φ Ω = Φ 3 1 i j C ji i C L M (55) 其中 ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ = Ω 0 0 0 1 0 0 0 1 0 (56) 因此，我們可以得到 A& =ΩA (57) 由於(57)是描繪系統在 critical 表面上之動態行為，我們在接下來的部分將探討系統在 critical 表面附近時的動態行為特性。定義符號：λ表示為接近critical 表面的參數值；Φλ =[φλ ψλ χλ]

(33)

為 critical 表面附近的基底向量空間；振幅方程式則為 A K A& = _λ ；Y =[A(t)φ_λ +B(t)ψ_λ +C(t)χ_λ]*Λ11+S~_λ。同理，依據雙對流理論，我們可以得到 ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ − − − = λ 2 1 0 1 0 0 0 1 0 c c c K (58) 因此，靠近 polycritical 的線性振幅方程式A& =K_λA可以得到以下之數學模型： ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ − − − = ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎣ ⎡ C B A c c c C B A 2 1 0 1 0 0 0 1 0 & & & (59) 將(59)整理後可得

&A&&+c₂A&&+c₁A& +c₀A=0 (60)

另外，我們接下來要考慮系統之非線性部分，將(60)擴充為非線性振幅方程式。考慮非線性則(60)可以改寫為

&A&&+c₂A&&+c₁A& +c₀A=g(A) (61)

其中g(A)表示非線性項。最後，依據(32)-(35)可得到 g(−A)=−g(A) (62) 從(62)可知，系統存在有最低階數之非線性項 ( 3) A O 。另外， ) (A g 可依據perturbation 理論表示為系統的狀態(A, B, C)之線 性組合，即 g(A)=t₁A3+t₂A2B+t₃AB2 +t₄B3 +t₅A2C+t₆ABC+t₇B2C 3 10 2 9 2 8AC t BC t C t + + + (63) 其中ti's為不同的常數係數。另外， 3 1A t 為(63)的 dominant term，則我們可以將三階非線性振幅方程式近似為

(34)

3 1 0 1 2A c A c A t A c

A&&+ &&+ & + =

& (64) 對於直流電漿火炬之行為特性可知，t₁等於1 和-1，則(64)為 3 0 1 2A cA c A A c

A&&+ &&+ &+ =±

& ₍₆₅₎ 為了簡化分析，引進新的時間比率s=εt和振幅A=ε3/2F且選擇 2 c = ε ，則(65)可以改寫為 3 0 1 2F F F F

F&&&+Ω &&+Ω & +Ω =± (66)

其中 ₂ 2 0 0 c c = Ω ， ₃ 2 1 1 c c = Ω ，Ω2 =1 且定義新的符號以簡化非線性振幅方程式的系統參數 _{δ q}ˆ ₌ 2(1₊_τ₊_σ₊_ξ) δˆ = 2₍ 6 − 2 2 + 2 2 −₄ π2 −2 2 ₎στξ 1 q q a m Rl a m Rs d R n Q a m Rl a m Rs d R n Q 2 2 2 2 2 2 2 2 ( ) (1 ) 4 (1 ) ˆ ₌₋ _σ _τ₊_ξ ₊ _στ ₊_ξ ₋ _π _σξ ₊_τ δ − ₍ ₎ 6 q στξ + τξ + σξ + στ + δˆ =− 2₍ 6 − 2 + 2 −₄ π2 −2 2 ₎στξ 3 q q a Rl a Rs d R n Q a q m Rl a q m Rs d R q n Q 2 2 2 2 2 2 2 2 2 2 4 4 ˆ ₌₋ _σ − ₊ _στ − ₋ _π − _σξ − δ ₍ ₎ 4 q τξ + σξ + στ + ξ + τ + σ + a m Rl a m Rs d R n Q 2 2 2 2 2 2 2 5 ( ) (1 ) 4 (1 ) ˆ ₌₋ _σ _τ₊_ξ ₊ _στ ₊_ξ ₋ _π _σξ ₊_τ δ − ₍ ₎ 6 q στξ + τξ + σξ + στ + δˆ = 2₍ 6 − 2 2 + 2 2 −₄ π2 −2 2 ₎στξ 6 q q a m Rl a m Rs d R n Q 則得到以下的系統參數 δ δ = ˆ ˆ 1 0 c 2 3 2 1 _ˆ ˆ ˆ ˆ δ δ − δ δ = c

(35)

₃6 2 5 4 2 _ˆ ˆ ˆ ˆ ˆ ˆ δ δ + δ δ − δ δ = c 最後，我們可以推導出系統參數(Ω₀、Ω₁和Ω₂) ₂ 6 5 2 4 5 1 0 ) ˆ ˆ ˆ ˆ ˆ ( ˆ ˆ δ + δ ⋅ δ − δ ⋅ δ δ ⋅ δ = Ω 3 6 5 2 4 7 3 8 2 1 ) ˆ ˆ ˆ ˆ ˆ ( ˆ ˆ ˆ ˆ δ + δ ⋅ δ − δ ⋅ δ δ ⋅⋅ δ − δ ⋅ δ = Ω Ω₂ =1 因此，最後我們可以利用非線性理論與分叉理論分析三階非線性振幅方程式(66)以及探討系統的動態行為。二、模擬及分析電漿火炬非線性之動態特性此部分之重點為分析電漿火炬之分叉特性，探討其顫動現象的發生原因。從學者A. K. Das 之研究發現電漿系統之顫動現象為一混沌行為，並非以往認知的隨機行為。然而，A. K. Das 的研究僅考慮單一系統參數發生變化時之系統動態，文獻討論此系統參數發生變化時的混沌現象且探討系統從穩定到發生混沌行為的參數範圍。經由分叉理論分析電漿系統之分叉特性發現影響電漿的因素可能不僅A. K. Das 考慮之單一參數，其他參數也可能影響系統使系統進入混沌現象甚至狀態發散。因此，本計畫之模擬分析將以前述建立之電漿火炬動態模型為基礎，探討電漿系統的所有動態範圍，分析系統的所有參數之範圍。另外，電漿系統之系統參數可能同時發生變化，系統可能會進入固定振盪週期之參數範圍。在接下來的部分，我們將分析系統參數之間的關聯性，探討可能發生顫動現象之因素。透過非

(36)

線性理論和分叉理論探討分叉參數變化之範圍，找出混沌現象的範圍。另外，我們也探討電漿火炬之其他可能發生的動態特性，找出系統可能發生的動態行為。在 A. K. Das 的研究中，學者僅將(66)之Ω₀設定為系統參數，其他則設定為常數，即Ω₁ =50及Ω₂ =1。在模擬分析時，本計畫則是將Ω₀、Ω₁和Ω₂設定為系統參數，模擬分析Ω₀、Ω₁和Ω₂ 發生變化時的關聯性以及探討系統參數影響顫動現象之成因。本計畫採用數值分析軟體 AUTO 模擬系統之分叉特性且利用 Matlab 模擬系統的時間響應以及相位圖。我們分 3 種可能性探討系統之動態行為：(i)僅Ω₀改變時之動態特性；(ii)考慮Ω₀或Ω₁ 改變時之動態特性；(iii)考慮Ω₀或Ω₂發生改變時之動態特性。以下將深入探討這三種情況之分叉特性。 Case 1: 僅改變Ω₀，Ω₁和Ω₂固定在此例子中，固定Ω₁=50和Ω₂ =1，將Ω₀設定為系統參數。因此，設x₁ =F、x₂ =F& 和x₃ =F&&，(66)可以改寫成狀態空間形式 x&₁ =x₂ x&₂ =x₃ (67) 3 1 3 2 1 3 x 50x x x x& =−μ − − − 其中Ω₀ =μ為分叉參數。在本計畫中之模擬圖，實線代表穩定的平衡點軌跡或穩定的limit cycle，虛線則表示不穩定的平衡點，模擬圖符號表示為 HB：Hopf bifurcation

(37)

SB：Stationary bifurcation PD：Period doubling 圖 1 為直流電漿火炬之分叉圖，表三列出系統之分叉點。從圖1 可以知道電漿火炬在μ=0會發生pitchfork bifurcation 且在0<μ<50時系統是漸進穩定的。當μ>50，則系統是不穩定的；當μ持續減小且通過原點時，系統由一個平衡點分叉為三個平衡點。另外，在μ=−25以及μ=50出現 Hopf bifurcation，即表示系統之Jacobian matrix 的特徵值會落在虛軸上，系統會出現週期解並產生振盪之動態行為。接著，當μ從μ>50持續減小且穿過μ =50時，系統之平衡點從 unstable focus 變成

stable focus，且產生一個 limit cycle。當μ減少並跨過μ=0時，系統由一個穩定的平衡點軌跡分叉為二個穩定的平衡點軌跡和一個不穩定的平衡點軌跡，即發生 supercritical pitchfork bifurcation。因此，系統的穩定性將取決於系統狀態之初始值。當μ<0且初始值在原點附近時，系統會被不穩定的平衡點吸引而造成不穩定。當μ減少至μ =−25，系統會發生二個 Hopf bifurcation(即 HB2 和 HB3)；這二個 HB 會分叉為一條穩定的 limit cycle 和一條不穩定的 limit cycle。當μ減少至μ =−108.823

和μ=−108.777時，穩定的limit cycle 軌跡則變成不穩定的 limit cycle 軌跡，此分叉現象表示系統會發生週期變化；因此，從模擬結果可知在μ=−108.823 和μ=−108.777 會發生 period doubling。從圖 1 可以知道，選擇 HB2 到 PD1(或 HB3 到 PD2)

之limit cycle 軌跡附近的初始值與分叉參數，系統之動態為一

(38)

則系統會發生二倍週期之振盪。當μ持續減少時，系統將發生週期性變化且進入混沌現象之區域。系統之行為開始出現週期性變化，則無法預估其動態行為；如此一來，我們就無法預期系統可能會發生之行為特性，導致無法應用系統分析理論與控制設計改善系統之穩定與使用性能；甚至μ仍持續變化減小，系統之狀態可能會發散而使系統發生不可預期的危險。從圖1 可以了解分叉參數與系統狀態行為之間的關係；接下來，我們透過時間響應與相位圖探討直流電漿火炬隨時間之動態行為。圖1. 直流電漿火炬之分叉圖 (μ v.s. x₁)

表三、Location of bifurcation points

) , ,

(x₁ x₂ x₃ μ

(39)

HB1 (0 ,0 ,0) 50 HB2 (5 ,0 ,0) -25 HB3 (−5 ,0 ,0) -25 PD1 (13.15 ,43.9 ,330.8) -108.823 PD2 (0.25 ,49.52 ,343.67) -108.777 我們在分叉圖上選取分叉點附近的分叉參數與初始值，從時間響應觀察系統之時域動態以及其相位圖。圖 2(a)為選取 3 = μ 且初始值設為(x₁₀,x₂₀,x₃₀)=(1,0,0)；圖 2(b)為選擇μ=3和初始值設為(x₁₀,x₂₀,x₃₀)=(−1,0,0)。圖 3(a)(b)則是選取μ=−5且初始值分別選取(x10,x20,x30)=(4,0,0)和(x10,x20,x30)=(−4,0,0)；圖 4(a)(b)則是選取μ=−5且初始值分別選取(x₁₀,x₂₀,x₃₀)=(2,0,0)和 ) 0 , 0 , 2 ( ) , , (x₁₀ x₂₀ x₃₀ = − 。從圖2 和圖 3 之時間響應圖與相位圖可以發現當0<μ<50及初始值選取在穩定平衡點附近時，系統狀態被穩定平衡點吸引而使系統狀態趨近於0。另外，圖 3 選取 0 25<μ< − 及初始值介於穩定平衡點軌跡與不穩定的 limit cycle 之間的值，模擬圖顯示此區域之系統狀態被穩定的平衡點吸引。同理，圖4 則在−25<μ<0的範圍選取介於二條穩定的平衡點軌跡之間的初始值，從時間響應圖與相位圖可以知道這個區域之初始值會被穩定的平衡點吸引，最後系統狀態趨近於穩定。另外，我們將分析系統之動態分叉，即當μ持續變化至發生Hopf bifurcation 時之行為特性，探討系統之振盪行為。接下來，我們分析分叉圖搭配時間響應圖和相位圖，解釋高溫直流電漿火炬之振盪週期以及週期變化。

(40)

(a) (b)

圖2. 時間響應(x₁v.s. t)和相位圖(x₁v.s.x₂)

(a) (b)

(41)

(a) (b)

圖4. 時間響應(x₁v.s. t)和相位圖(x₁v.s.x₂)

(a) (b)

(42)

(a) (b)

圖 6. 時間響應(x₁v.s. t)和相位圖(x₁v.s.x₂)

(a) (b)

(43)

(a) (b)

(44)

(a) (b)

(45)

(a) (b) 圖 12. 時間響應(x₁v.s. t)和相位圖(x₁v.s.x₂) 圖 5(a)(b)選取μ=60且初始值(x₁₀,x₂₀,x₃₀)分別設為(5,0,0) 和(−5,0,0)，搭配圖 1 之分叉圖可以發現在此區域內之平衡點為不穩定的。圖 6(a)(b)選取初始值(x₁₀,x₂₀,x₃₀)=(10,0,0)且μ分別為5 和-3；圖 6 之時間響應圖和相位圖表示在此區域內之系統狀態被穩定的平衡點吸引，系統將趨近於穩定。圖7(a)(b)~ 圖 9(a)(b)則選取μ=−50及初始值(x₁₀,x₂₀,x₃₀)分別選取(15,0,0) 和(10,0,0)、(7,0,0)和(2,0,0)及(−5,0,0)和(−8,0,0)。從圖 7~圖 9 之時間響應與相位圖的模擬結果說明了這些區域會發生 Hopf bifurcation，系統會產生振盪之週期解，其模擬之行為特性符合圖1 分叉圖之分析結果。另外，從圖 1 可以知道當μ<−108.8 時，系統的振盪週期會發生變化甚至進入混沌行為之區域。我們將逐步探討從 Hopf bifurcation 延伸軌跡的動態特性，解釋

(46)

系統從穩定狀態開始產生固定週期的振盪行為之參數變化，最後進入非固定週期之混沌行為到系統狀態發散。我們在圖 7~ 圖 9 之時間響應圖與相位圖說明系統之振盪週期變化以及混沌現象。圖10(a)(b)~圖 12(a)(b)則選取相同的μ=−50和不同的初始值(x₁₀,x₂₀,x₃₀)=(20,0,0)和(17,0,0)、(13,0,0)和(5,0,0)以及 (1,0,0)和(-1,0,0)。圖 10~圖 12 表示這些區域的分叉參數和初始值會使系統發散，系統會因初始值的不同而造成動態行為由振盪行為轉變成發散。圖 13(a)(b)~圖 16(a)(b)使用相同的初始值(x₁₀,x₂₀,x₃₀)= ) 0 , 0 , 5 ( 且分別選取不同的μ值，即μ=-20、-80、-115、-120、 -122.15、-130、-144 和-146。從圖 13~圖 16 可以看出，μ=−20 為穩定的；當μ=−80時，則發生一倍週期的振盪行為；當 115 − = μ 時，則發生二倍週期的振盪行為；而當μ=−120和 15 . 122 − 時，則分別發生四倍週期和八倍週期的振盪行為；當 130 − = μ 和−144時，系統之振盪行為產生不規則之週期變化，即表示系統之狀態發生混沌現象；最後，當μ=−146時，系統之動態行為是不穩定的。因此，透過模擬圖之分析結果發現，系統會隨著μ值變化而使系統發生週期變化甚至是混沌行為。在此 case 1 中，我們僅考慮Ω₀發生變化時之系統可能發生的動態行為。此分析結果包含學者A. K. Das 探討之參數範圍；另外，case 1 也探討比文獻更廣的參數變化範圍，說明系統狀態發生不穩定的參數區域，此結果將可以避免系統發生不可預期的危險。

(47)

(a) (b)

(48)

(a) (b)