T.C. ULUDAĞ ÜNİVERSİTESİ FEN BİLİMLERİ ENSTİTÜSÜ DİELEKTRİK ISITMANIN FDTD METODUYLA MODELLENMESİ HATİCE BİÇER YÜKSEK LİSANS TEZİ ELEKTRONİK MÜHENDİSLİĞİ ANA BİLİM DALI BURSA 2006 T.C. ULUDAĞ ÜNİVERSİTESİ FEN BİLİMLERİ ENSTİTÜSÜ DİELEKTRİK ISITMANIN FDTD METODUYLA MODELLENMESİ HATİCE BİÇER YÜKSEK LİSANS TEZİ ELEKTRONİK MÜHENDİSLİĞİ ANA BİLİM DALI BURSA 2006 T.C. ULUDAĞ ÜNİVERSİTESİ FEN BİLİMLERİ ENSTİTÜSÜ DİELEKTRİK ISITMANIN FDTD METODUYLA MODELLENMESİ HATİCE BİÇER YÜKSEK LİSANS TEZİ ELEKTRONİK MÜHENDİSLİĞİ ANABİLİM DALI Bu Tez 28.07.2006 tarihinde aşağıdaki jüri tarafından oybirliği/oy çokluğu ile kabul edilmiştir. Prof.Dr.ALİ OKTAY Prof. Dr. A. Burak POLAT Prof. Dr. Naim DEREBAŞI (Danışman ) i ÖZET Endüstriyel mikrodalga ısıtma, kurutma ve dielektrik malzemelerin işlenmesi, mikrodalga enerjinin sağladığı avantajlar sayesinde hızla ve geniş bir sahada gelişmektedir. Mikrodalga ısıtma sürecinde, malzeme içindeki sıcaklık dağılımı güç kaynağının ürettiği elektrik alanın değeri ile Isı Transfer Denklemi’nin çözümü yapılarak elde edilmektedir. Isı transfer eşitliğinin çözümü için farklı sayısal teknikler bulunmaktadır. Bu tezde zamana bağlı Sonlu-Farklar-Zaman-Domeni (FDTD) yöntemi uygulanmıştır. Bir dalga kılavuzu içindeki dielektrik malzemenin mikrodalga ısıtması, deneysel ve sayısal olarak analiz edilmiştir. Bu amaçla malzeme olarak tuğla, motor yağı ve su seçilmiştir. Deneysel düzenek WR340 dalga kılavuzundan yayılan TE101 tek modu içermektedir. Isıtma koşullarında, malzemedeki sıcaklık değerleri sıcaklık ölçüm cihazı ile ölçülmüştür. Sayısal olarak; dalga kılavuzundaki elektromanyetizma analizi için, Maxwell Denklemleri’nin temeli Yee notasyonu olan 3 Boyutlu-FDTD yöntemi kullanılarak ayrıklaştırılması ile yapılmıştır. Dielektrik malzemenin üzerindeki sıcaklık dağılımını elde etmek için, Isı Transfer Eşitliği Crank-Nicholson şeması göz önüne alınarak 2Boyutlu-FDTD yöntemi yardımıyla ayrıklaştırılmıştır. Kılavuzdaki malzemelerin deneysel çalışması yapılmış ve sayısal sonuçlar ile karşılaştırılarak zamana bağlı eğrileri elde edilmiştir. Anahtar Kelimeler Sonlu Farklar Zaman Domeni yöntemi (FDTD), sayısal analiz, Maxwell Denklemi, dalga kılavuzu, elektromanyetik alanlar, ısı transfer eşitliği, Crank- Nicholson, mikrodalga ısıtma, dielektrik. ii ABSTRACT Industrial microwave heating, drying and processing of dielectiric materials are rapidly growing in number and scope because of the advantage of microwave energy. During microwave heating, the temperature distribution within the material is obtained by solving the heat transfer equation with the electric field as the power source. There are different techniques that have been used to solve the heat transfer equation. Due to its time- dependency the techniques are usually employed in the time domain. In this thessis the finite difference time domain (FDTD) method is used . Microwave heating of a dielectric in a waveguide is experimentally and numerically analyzed by the FDTD method. For numeric and experimantal purpose, engine oil, clinker and water are chosen. The experimental set up consists of a single mode TE101 cavity excited by a WR340 waveguide with a coupling iris. For heating conditions, the temperature distribution on water inside the cavity is measured with temperature measurement device. For computational analysis of electromagnetic fields inside the waveguide, The Maxwell’s Equations are considered using 3-Dimension FDTD method based on Yee’s notation. The heat transfer equation is also considered using 2- Dimension FDTD with Crank-Nicholson schema in order to obtain temperature distribution on the dielectric materials. The comparison of experimental heating data and computational results for materials are presented as a function of time. Keywords Finite-difference-time-domain-method (FDTD), numeric analysis,Maxwell Equation, waveguide applicator, electromagnetic fields, heat transfer equation, Crank- Nicholson, mikrowave heating, dielectric. iii İÇİNDEKİLER Özet ................................................................................................................................ i Abstract .......................................................................................................................... ii İçindekiler ...................................................................................................................... iii Simgeler Dizini .............................................................................................................. v Şekiller Dizini ................................................................................................................ vi 1 Giriş............................................................................................................................. 1 2 Kaynak Araştırması..................................................................................................... 4 3 Materyal ve Yöntem.................................................................................................... 7 3.1 Mikrodalga Isıtmanın Elektromanyetizma temelleri ........................................... 7 3.2 Isı Transfer Denklemi ........................................................................................... 10 3.2.1 İletim ile ısı taşınması ................................................................................. 13 3.2.2 Işıma ile ısı iletimi..................................................................................... 14 3.2.3 Taşınma ile ısı iletimi ................................................................................ 14 3.3 Mikrodalga Isıtmanın Sayısal Modellenmesi........................................................ 15 3.3.1. FDTD Yöntemi (Sonlu Farklar Zaman Domeni Yöntemi) ...................... 15 3.3.1.1 FDTD Yönteminde Türev İfadelerinin Taylor Serileri İle Analizi ................................................................................................... 16 3.3.1.2 Birinci Dereceden Türevin Sonlu Farklar Yaklaşımı .................................................................................... 17 3.3.1.3 İkinci Dereceden Türevin Sonlu Farklar Yaklaşımı ................................................................................... 19 3.4 Elektromanyetik Alanların FDTD ile Analizi........................................................ 20 3.4.1 FDTD Yönteminin Uygulanmasında Temel İlkeler ..................................... 20 3.4.2 FDTD Yönteminin Maxwell Denklemlerine Uygulanması .......................... 22 iv 3.4.3 FDTD Yönteminde Sınır Koşulları ............................................................... 27 3.4.3.1 Birinci Dereceden Mur Yutucu Sınır Koşulu Yönteminin Matematiksel İfadesi............................................... 29 3.4.3.2 Dielektrik Yüzey Üzerinde Sınır Koşullarının FDTD Yöntemi ile Analizi ............................................ 32 3.4.4 Kararlılık İçin Hücre Boyutunun Belirlenmesi ............................................. 35 3.4.5 Kararlılık İçin Zaman Adımı Boyutu ........................................................... 36 3.5 Isı Transfer Denkleminin FDTD Yöntemi ile Analizi ........................................... 37 3.5.1 Crank-Nicolson Yöntemi İle Isı Denkleminin Ayrıklaştırılması…………..38 3.5.2 Isı Transfer Denkleminde Sınır Koşulu ....................................................... 40 3.6 Gauss-Seidel Yöntemi............................................................................................ 44 3.7 Malzeme Parametreleri ......................................................................................... 47 3.8 Deneysel ve FDTD Sonuçlar ................................................................................ 49 4 Araştırma Sonuçları ve Tartışma .............................................................................. 54 5 Kaynaklar ................................................................................................................. 60 Teşekkür Özgeçmiş v SİMGELER DİZİNİ erx e r y e r z kartezyen koordinatlar birim vektörleri ε r bağıl dielektrik geçirgenlik ε r′ dielektrik sabiti (bağıl dielektrik geçirgenliğinin gerçel kısmı) ε r′′ dielektrik kayıp faktörü (bağıl dielektrik geçirgenliğinin sanal kısmı) ε 0 serbest uzay dielektrik geçirgenliği σ iletkenlik (S/m) µr bağıl manyetik geçirgenlik µ0 serbest uzay manyetik geçirgenlik ω açısal frekans (rad/sn) r E elektrik alan vektörü (V/m) r H manyetik alan vektörü (A/m) r D deplasman akı yoğunluğu r B manyetik akı yoğunluğu r J akım yoğunluğu (A/m2) P mikrodalga kayıp güç yoğunluğu (W/m3) nr Birim vektör k Isı iletkenlik katsayısı (W/mºC) Cp Özgül ısı (J/kgºC ) ρ yoğunluk (kg/m3) Hc konveksiyon katsayısı(W/m2ºC) ∆x konum adımı (mm) ∆t zaman adımı (sn) C0 ışık hızı vi ŞEKİLLER DİZİNİ Şekil 3.1: Elektrik ve Manyetik alanların herhangi bir t anında yayınımı..................... 8 Sekil 3.2: Elektromanyetik alanın malzemeye etkisi..................................................... 11 Şekil 3.3: Yee Birim Hücresi......................................................................................... 21 Şekil 3.4 Yutucu maddeye gelen dalganın durumu...................................................................................28 Şekil 3.5: MUR yönteminde açık uçlara yutucu yüzey yerleştirilmesi ......................... 31 Şekil 3.6 İki dielektrik ortam arasındaki sınır şartı........................................................ 33 Şekil3.7: Farklı dielektrik değerlere sahip hücrelerde E alanı....................................... 34 Şekil 3.8: Crank-Nicolson yöntemi için Sonlu Fark Molekülleri.................................. 39 Şekil 3.9: Malzeme yüzeyindeki sınırlar ve ısının dağılım yönü .................................. 40 Şekil 3.10: Tuğla için dielektrik kayıp sabitinin sıcaklıkla değişimi............................. 47 Şekil 3.11: Suyun Sıcaklığa Bağlı Dielektrik Kayıp Faktörü ε" ................................... 48 Şekil 3.12: Deney Düzeneği .......................................................................................... 49 Şekil 3.13: P=100Watt için suyun deneysel ve FDTD sonuçları .................................. 51 Şekil 3.14: P=180w ve P=72w için motor yağının deney sonuçları.............................. 51 Şekil 3.15: Kullanılmış motor yağının deneysel ve FDTD sonuçları............................ 52 Şekil 3.16: P=72W için motor yağının deneysel ve FDTD sonucu ............................. 52 Şekil 3.17: P=360W için tuğla (Al2O3SiO2) deneysel ve FDTD sonucu ...................... 53 Çizelge 3.1: Tuğlanın termo-fiziksel malzeme özellikleri ............................................ 47 Çizelge 3.2: Motor yağının termo-fiziksel malzeme özellikleri .................................... 48 Çizelge 3.3: Suyun termo-fiziksel malzeme özellikleri................................................. 48 1 1. GİRİŞ Mikrodalga enerjisinin dielektrik kayıplarla ısıya dönüşmesi uygulama olarak endüstride kurutma, pişirme, sinterleme olarak karşımıza çıkar ayrıca bilimsel çalışmalarda da en çok çalışılan konulardandır. Kısa etki zamanına sahip olması ve malzemeyi içten ısıtmaya başlaması sebebiyle endüstride, ısıl işlemlerde geleneksel ısıtma yöntemlerinin yerini mikrodalga enerjisi almaktadır. Mikrodalga enerjisinin endüstride ısıtma olarak kullanılması 1950’lerden beri geliştirilen bir çalışmadır. 1921 yılında A.W. Hull tarafından bulunan ve etkin olarak 2. Dünya Savaşında radarlarda kullanılmaya başlanan magnetron teknolojisindeki gelişmeler doğrultusunda endüstriyel mikrodalga çalışmaları hızlanmış ve yaygınlaşmıştır. Diğer taraftan, 1954 yılında malzeme özelliklerinin gerçek verilerini elde etmek için Von Hippel tarafından önemli araştırmalar yapılmıştır. 100Hz-10GHz frekans bölgesinde organik ve inorganik birçok malzemenin özellikleri araştırılmıştır. Bu ilk çalışmalar endüstriyel yüksek frekans ve mikrodalga enerjisinin kullanımında esas teşkil etmiştir. Günümüzde, ileri endüstrilerin bir çok uygulamalarında ve tıpta bazı hastalıkların (kanser ve diyaterapi gibi..) tedavisinde yüksek ve çok yüksek frekanslı elektromanyetik dalgalar geniş ölçüde kullanılmakta, bir çok düzenek ve cihaz bulunmaktadır. Yine 1950’li yıllarda, malzeme özelliklerine ait verilere ek olarak magnetronların güç beslemeleri ve yardımcı cihazların tasarımında önemli gelişmeler olmuştur. Mikrodalga ısıtma ile özellikle kimyasal reaksiyonlarda ortaya çıkan sıcaklıkla etkilenen, kimyasal ve biyokimyasal analiz aracı olarak yararlanmak, onunla kimyasal bileşiklerin dielektrik davranışlarını doğru biçimde ortaya koymak ve mikrodalga saçılma etkilerini karakterize etmek mümkündür. Malzemelerin dielektrik kayıplarını kullanarak, mikrodalga enerjisi ile ısıtılması, klasik ısıtma işlemlerine göre birçok avantaj sunmaktadır. Klasik ısıtmada, malzemenin yüzeyi ısıtılmakta ve iç bölgesinin ısıtılması ise ısının yüzeyden içeri doğru iletim (conduction) ile yayılması ile sağlanmaktadır. Bu, malzemede istenmeyen neticeleri doğurabilecek, büyük sıcaklık gradyanının oluşmasına sebep olacaktır. Mikrodalga ile malzemelerin ısıtılmasında, mikrodalgalar malzemenin iç kısımlarına kadar daldıkları için daha düzgün dağılımlı ve hacimsel bir ısıtma sağlayacaktır. Bu da sıcaklık gradyanının çok küçük olmasını sağlayacaktır. Isı transferinin hızlı ve ısıtma sürecinin kısa olmasını sağlayacaktır. Mikrodalgaların kurutma, pişirme ve yüksek frekans 2 radyasyonunun özel etkilerinden yararlanarak yapılan endüstriyel uygulama imkânları oldukça geniştir. Son zamanlarda telekomünikasyon ve deteksiyon amaçlı kullanımlarda da mikrodalga enerjisi kullanılmaktadır. Maxwell denklemlerinin belirli bir geometriye uygulanmasıyla elektromanyetik alanların elde edilmesi, o geometriye ait elekromanyetik analizin yapılması anlamını taşır. Eğer fiziksel geometri yeterince basit ise Maxwell denklemlerin analitik çözümü ile elektromanyetik analizin yapılması mümkündür. Fakat içerisinde, gelişi güzel geometride, homojen olmayan kayıplı dielektrik malzeme içeren ve bünyesinde karmaşık fiziksel olayları barındıran mikrodalga ısıtma düzenlerinin elektromanyetik analizinin analitik olarak yapılması güçtür ve dolayısıyla bilgisayar ortamında modelleme gerektirir. Bu model, problemin fiziksel geometrisi sınır koşulları da dahil olmak üzere, elektromanyetik gücü sağlayan kaynak ve çözümün doğruluğuna ve hassasiyetine zarar vermeyecek bazı yaklaşıklıklar ihtiva eder. Hesaplamalı Elektromanyetik (Computational Electromagnetics-CEM ) günümüzde, bilgisayar teknolojisindeki gelişmelerle birlikte büyük bir ilerleme gösteren, analitik çözümleri çok zor olan büyük ve karmaşık elektromanyetik problemlerin bilgisayar ortamında çözümüyle ilgilenen yeni bir daldır. Özellikle yüzyılımızın ikinci yarısından sonra, birbirlerine göre değişik avantajlar sunan pek çok hesaplamalı elektromanyetizma yöntemleri geliştirilmiştir. Sonlu elemanlar yöntemi (Finite Element Method-FEM), Sonlu Farklar Zaman Domeni Yöntemi ( Finite Difference Time Method-FDTD), İletim Hat Matrisi Yöntemi (Transmission Line Matrix Method-TLM) ve Moment Yöntemi (Method of Moment-MOM) bu hesaplamalı elekromanyetizma yöntemlerinde en çok kullanılan ve üzerinde çalışmalar yapılanlarıdır. Hesaplamalı Elektromanyetizma yöntemlerine bakıldığında, bir kısmı frekans domeninde (harmonik maxwell denklemleri) ve bir kısmı da zaman domeninde (zamana bağlı maxwell denklemleri) iş görürler. Bu tezde, bir dalga kılavuzu içinde bulunan dielektrik malzemenin, FDTD- (Sonlu Farklar Zaman Domeni Yöntemi) ile modellenmesi incelenmiştir. Tezin malzeme ve yöntem kısmında, mikrodalga ısıtma sisteminin modellenmesi hakkında bilgi verilmiştir. Bu kısımda bir dalga kılavuzunda bulunan dielektrik malzemenin ısıtılması olayını incelerken, malzeme ile dalga arasındaki etkileşim, elektromanyetik olaylarla oluştuğundan bu oluşumu inceleyebilmek için, Maxwell Denklemlerinden 3 kısaca bahsedilmiştir. Isı transfer denkleminin çıkarılması ve ısı iletim yolları üzerinde durulmuştur. Mikrodalga ısıtmanın sayısal olarak modellenmesi ve sayısal yöntem olarak FDTD yöntemi hakkında bilgi verilip elektromanyetik alan analizi ve ısı transfer denkleminin analizi bu FDTD yöntemi ile yapılmıştır. Mikrodalga ısıtma sisteminin malzeme üzerindeki etkisi deneysel olarak da incelenmiştir. Isıtılan malzemelerin parametreleri ve mikrodalga ısıtma sonucu elde edilen sayısal sonuçlar ve deneysel çalışma yer almaktadır. Deneysel çalışma için gerekli deney düzeneği verilmiş FDTD modellemesi sonucu alınan veriler ile deneyden alınan veriler aynı eğride karşılaştırılmıştır. Son kısımda araştırma sonuçları ve tartışma yer almaktadır. 4 2. KAYNAK ARAŞTIRMASI 1966’da Yee’nin zamana bağlı Maxwell denklemlerini sonlu farkları kullanarak çözmesiyle başlayan FDTD (finite difference time domain=sonlu farklar zaman domeni) yöntemi, bugün elektromanyetikte en çok kullanılan sayısal yöntemlerin başında gelir. Yöntemin aslı, Maxwell denklemlerinde alan büyüklüklerinin zamana ve konuma bağlı türevleri yerine, bu büyüklüklerin Taylor serisine açılımı ile elde edilen fark denklemlerinin kullanılmasıyla, sürekli olan elektromanyetik alanların sonlu hacimlerle ve eş zaman aralıklarıyla örneklenmesine dayanır. Neticede türev ifadelerine sahip sürekli denklem, ayrık fark denklemlerine dönüşür. Bu fark denkleminin çözümü ile elde edilen değerler, süreli alanı ilgili zaman ve konumdaki değerine karşılık gelir. Elektromanyetik teoremi için Maxwell denklemleri FDTD yöntemine, benzetilerek analizi yapılır. Bu yöntemin düşük maliyetli olması yanında, hızlı bilgisayarlara uygulanması sonucu, sonuçlara daha hızlı ve yaklaşık bir değerde ulaşılmış olması gibi avantajları vardır. Diğer yöntemler tarafından çözülmesi aşırı derecede zor olan geometrik sorunlar için Maxwell denklemlerinin çözümlerini kolayca hesaplayarak bulabiliriz. Basitlik ve gücün bu birleşimi FDTD’ yi böyle popüler bir yöntem yapıyor. Mikrodalga jenaratörlü bir mikrodalga ısıtma sistemi için dalga kılavuzu, boşluk ve dielektrik yükleme gibi birçok sayısal konuları dikkate almak gerekir. Süreksiz ara yüzler ve sınır koşullarının hesaplanmasındaki yanlış varsayımlar, yanıltıcı benzetim sonuçlarına yol açabildiğinden, çözümün niteliği için bu kabuller çok önemlidir. İlk mikrodalga fırınlar 1940’lı yıllarda daha çok gıda endüstrisinde ortaya çıkmıştır. 1960’ların sonlarında dielektrik ısıtma için mikrodalga kılavuzun ilk kullanımı görülmüştür. 1979’da Holme ve Metaxas deneysel olarak mikrodalga işleme ile konveksiyonel kurutma teknikleri çalışmışlardır. Nylon halıları kurutmak için konveksiyonel sistemin parçası olarak mikrodalga uygulayıcıyı deneysel olarak kullanmışlardır. Nylon halıların mikrodalga ile kurutulması basamak basamak artar. Halıdaki nem miktarı azalırken kılavuzun duvarındaki güç kaybı da yükselir. Ayrıklaştırmada ve FDTD yönteminin temelinin oluşmasında Rıchtmeyer ve Morion (1967), Smith (1978), Anderson, Tannehill ve Pleicher (1984), Berezin ve Zhidkov (1965) ve Roache (1976) gibi araştırmacılar vardır. Son yıllarda, mikrodalga 5 ısıtma sistemlerine elektromanyetik alan ve güç dağılımının elde edilmesi amacıyla FDTD yönteminin kullanılmasında büyük gelişmeler kaydedildi. Bu yöntemi kapalı bir mikrodalga ısıtma sistemine ilk uygulayan De Pourcq (1985)’dır. De Pourcq, sonu kısa devre edilmiş ve içine kayıplı dielektrik malzeme yerleştirilmiş transmisyon borusuna FDTD yöntemini uygulayarak malzeme üzerindeki güç dağılımlarını elde etmiştir. Daha sonra kısmi olarak kayıplı bir dielektrik malzeme ile yüklenmiş çok modlu mikrodalga rezanatöre uygulanmıştır. FDTD yöntemiyle elde edilen sonuçların hassasiyeti yeterince iyiydi ve koşturma zamanı çok kısaydı. Son zamanlarda Iskander (1994), Pierre-Yves (1994), Lui-Kun ve William (1995) FDTD yöntemini, farklı dielektrik ısıtma problemlerini de, ısıtma profilini elde etmek için kullandılar ve bu yöntemin bu problemlerdeki etkinliğini ortaya koydular. Sundberg ve arkadaşları (1996), FDTD yöntemini endüstriyel mikrodalga fırınına uygulamışlar ve malzeme üzerindeki güç dağılımını çıkartarak çeşitli mikrodalga fırınları için tasarım prensipleri oluşturmuşlardır. Zhao ve Turner (1996) üç boyutlu bir fırın içinde mikrodalga ışımasını FDTD kullanarak modellemişler. Liu ve ark. (1996) rijit dalga kılavuzu kullanarak polimerlerin birleştirilmesinde meydana gelen olayları anlamak için Maxwell Denklemlerini ısı denklemiyle birlikte FDTD kullanarak çözmüşlerdir. Seramik yüklü bir cavity için (1996) Clemens ve Saltiel sayısal olarak incelemiştir. Dielektrik malzemeyle kısmen doldurulmuş bir mikrodalga aplikatörde, mikrodalga ısıtmanın iki boyutlu sayısal analizi 1998’de S. Tada, R. Echigo ve ark. tarafından yapılmıştır. Oktay ve Akman (2001) bir dalga kılavuzu içinde seramik malzemelerin kurutulmasının modellenmesi için Maxwell denklemleri, ısı transfer ve kütle difüzyon denklemlerini birlikte FDTD yöntemiyle çözmüşler ve kurutma karakteristiğini çıkarmışlardır. Silviu Frandos (2003) dalga kılavuzunda mikrodalga ısıtma uygulamasını FDTD yöntemiyle yapmıştır. Yuichi Funawatashi ve Tateyuki Suzuki (2003) fırın içindeki dielektrik malzemenin mikrodalga ısıtmayla sayısal analizini FDTD yöntemiyle incelemişlerdir. Maxwell denklemlerinin farklı yapılardaki transmisyon borularının alan bileşenlerinin hesaplanmasında kullanılması suretiyle bu alan bileşenlerinin hesaplanması borunun boyutu arttıkça hesaplamalar analitik olarak zorlaşmakta hatta imkânsız hale gelmektedir. Bu alan bileşenlerini benzetmede kullanılan FDTD (sonlu zamanlar yöntemi), FEM(sonlu elemanlar yöntemi), iletim hat Matrisleri bize büyük 6 kolaylık sağlar. Bu yöntemler arasında işlem yoğunluğu en az olan FDTD yöntemidir. FEM yöntemi kapalı bir çözüm sunduğundan daha kesin çözüm verirken daha fazla program koşturma zamanı, daha güçlü donanım ve daha karmaşık yazılım gerektirir. Bunun yanında FDTD yöntemi, iteratif çözüm sunduğundan her zaman yakınsamama tehlikesi mevcut iken daha az program koşturma süresi, daha basit programlama olanağı sunar. 7 3. MATERYAL VE YÖNTEM Mikrodalga ısıtma teknolojisi günümüzde, birçok endüstriyel süreçte artan bir hızla etkin olarak kullanılmaktadır. Bir mikrodalga ısıtma sisteminin kalitesi, içerisine yerleştirilmiş yük üzerinde mümkün olduğunca düzgün sıcaklık dağılımı sağlayabilmesi ile değerlendirilir. Bu özelliğin sağlanması ile endüstriyel süreçlerde daha fazla verim ve daha kaliteli ısıtma söz konusu olacaktır. Belirli bir malzeme (yük) için, ısıtma sisteminin tasarımı ve optimizasyonu bu özelliğin sağlanmasında önemlidir. Tasarımın istenilen özellikler doğrultusunda yapılabilmesi de mikrodalga ısıtma esnasında oluşan fiziksel olayların (elektromanyetik alan, ısı transferi...) analizi ile mümkündür. Karmaşık ve zincirleme olaylar dizisini bünyesinde barındıran mikrodalga ısıtma sistemlerinde, ısıtma sürecinin kontrolü ve optimizasyonu, ancak meydana gelen fiziksel olayların matematiksel ifadelerinin birbirleri ile olan ilişkileri de göz önüne alınarak çözülmesiyle mümkün olacaktır. Bu sebeple bu bölümde, bir sonraki bölümde değinilecek mikrodalga ısıtmanın sayısal modellenmesine temel sağlayacak matematiksel bağıntılar ve bunların anlamları üzerinde durulacaktır. 3.1 Mikrodalga Isıtmanın Elektromanyetizma Temelleri Mikrodalga gücü, bir hacim içerisinde yayıldığı zaman, bu hacmin yüzeyinden itibaren belirli bir derinliğe kadar olan kısımda, malzeme ile elektromanyetik alan arasında karşılıklı bir etkileşim olur. Malzemenin elektrik özellikleri, bu etkileşime bağlı olarak değişir. Burada, bu etkileşimin en belirgin özelliği, elektromanyetik enerjinin ısıya dönüşmesidir. Malzeme ile dalga arasındaki etkileşim, elektromanyetik olaylarla oluştuğundan bu oluşumu inceleyebilmek için, Maxwell Denklemlerini göz önüne almak gerekir. Bu denklemler elektromanyetik dalga büyüklükleri arasındaki bağıntıları gösteren denklemlerdir. 8 Şekil 3.1: Elektrik ve Manyetik alanların herhangi bir t anında yayınımı Maxwell denklemleri şu iki temel kurala dayanmaktadır: 1. Zamana göre değişen elektrik alan, bir manyetik alan yaratmaktadır. 2. Zamana göre değişen manyetik alan, bir elektrik alan yaratmaktadır. Maxwell dalga denklemleri, elektromanyetik dalganın iki bileşeni olan E elektriksel alanı ve H manyetik alanı arasındaki bağıntıları eşitlik (3.1)’den (3.4)’ye kadar göründüğü gibi: r r ∂B ∇ × E = − (Maxwell Faraday Denklemi ) (3.1) ∂t r r r ∇ × H ∂D= + J s (Maxwell Ampere Denklemi ) (3.2)∂t r ∇.D = ρ (Maxwell Gauss Denklemi ) (3.3) r ∇.B = 0 (Maxwell Gauss Denklemi ) (3.4) r r r r Burada E , H , D , B sırasıyla elektrik ve manyetik alan şiddetleri, elektrik ve manyetik r akı yoğunluklarıdır. J s ve ρ sırasıyla kaynak tarafından sağlanan elektrik akımı ve yük yoğunluklarıdır. (Pozar,D.M) Ortam ile elektromanyetik alanlar arasında şu bağıntılar da vardır: r r r v r r B=µ.H J =σ.E D=ε.E (3.5) 9 Böylece eşitliklerimiz: r ∂H 1 r = − (∇× E) (3.6) ∂t µ r ∂E σ r 1 r = − E + (∇× H ) (3.7) ∂t ε ε Denklemleri x,y,z koordinatlarına bağlı olarak Elektrik ve Manyetik alanlarını kayıplı dielektrik içeren ortamda açık bir şekilde yazabiliriz. Manyetik alan bileşenleri: ∂H x ∂Eµ y ∂E= − z (3.8.a) ∂t ∂z ∂y ∂H µ y ∂E z ∂E= − x (3.8.b) ∂t ∂x ∂z ∂E µ ∂H z ∂E= x − y (3.8.c) ∂t ∂y ∂x Elektrik alan bileşenleri ise: ∂Ex E ∂H ∂H ε +σ x = z − y (3.9.a) ∂t ∂y ∂z ∂E ε y ∂H E ∂H+σ y = x − z (3.9.b) ∂t ∂z ∂x ∂Ez ∂H ∂Hε +σE y x ∂t z = − (3.9.c) ∂x ∂y σ =ωε0ε " r ε =ε ε ε =ε ' 0 r r r − jε '' r (3.10) Bir malzemenin elektrik enerjisini yutma özelliği, iki parametre ile tespit edilmiştir. Bunlardan birisi, σ iletkenlik katsayısı, diğeri de ε = ε ' − jε ''r r r ile verilen dielektrik katsayısıdır. ε ' ortamın dielektrik sabiti, ε ''r r ise ortamın dielektrik kayıp faktörüdür. Bu iki parametrenin değerlerine göre ortamdaki bir elektromanyetik ışımanın ısıya dönüşüp dönüşmediğini karakterize etmek mümkündür. Makroskobik boyutta incelendiğinde, malzeme içerisinde mevcut olan serbest yüklerin (iyonlar veya 10 çok az miktarda bulunabilecek serbest elektronlar) elektrik alandan etkilenmeleriyle yapacakları osilasyonları neticesinde ve polar moleküllerin (dipolar momentleri sıfırdan farklı olan moleküllerin) elektrik alanla yapacakları dönme hareketinden dolayı oluşan dielektrik relaksasyon ile olmak üzere iki farklı mekanizma malzemenin ısınmasına katkıda bulunacaktır. Bu mekanizmalardan ilki, ortamın iletkenlik sabiti ile ikinci de ortamın dielektrik kayıp faktörü ile temsil edilir. Ortamda kayıplardan dolayı ısıya dönüşen güç: 1 r rP = − ∫∫Re(E×H *)∂s (3.11)2 S ifadesinin reel kısımları tarafından elde edilen değerdir. Buradan, eşitlik(3.14)de birim hacim başına düşen güç yoğunluğu(Oktay ve Akman,2000): 1 r 2P = σ E (3.12) 2 şeklinde verilir. 3.2 Isı Transfer Denklemi Klasik ısı tekniğinde, kızılötesi lambalarından, rezistanslardan çıkan ısı bir malzemeye geçerken, önce yüzeyden başlayarak malzemenin içine doğru difüzyonla yayılır. Isı enerjisi veya ısı miktarı en genel hali ile maddenin moleküllerinin titreşimlerine bağlı düzensiz küçük dağılımların bir toplamıdır. Bir madde veya cisim üzerinden enerjinin geçişi bu molekül titreşimlerin en yakın moleküle transferi ile gerçekleşir. Cismin veya maddenin sıcaklığını artırdığımızda bu titreşimler büyük olur. Bazı durumlarda bu titreşimler çok büyük olduğu zaman cisimde bir faz değişikliği veya kimyasal bir olay meydana gelir. Bir cismin sıcaklığı ∆t artırıldığı zaman bu cisme verilen enerji miktarı: θ = m.c p∆t (3.13) 11 ile ifade edilir. Burada, Cp özgül ısı, m ise kütledir. Şekil 3.2’deki gibi bir hacim üzerinde, bir P gücü gelsin. Hacmin içinde bir T1 noktası ve dışında bir T2 noktası alınsın. Burada dikkat edilmesi gereken, hacmi V olan bir cismin elektromanyetik alan altındayken, malzemenin içinde bir noktadaki T1 sıcaklığı zamanla belirli değişim gösterir. Bu hacmin yüzeyi S ise, hacim S yüzeyi boyunca ikinci ortamla etkileşim halindedir. Bu ikinci ortamın karakteristikleri sıcaklık, basınç, nem gibi büyüklükleri zamanla değişecektir. (2) V (1) P .T1 P . T2 Şekil 3.2: Elektromanyetik alanın malzemeye etkisi Bu malzemenin herhangi bir andaki sıcaklığı T1, ısı denklemiyle ortaya çıkan, malzemenin içine dalan ısı enerjisinin ve sınır şartlarının belirlediği bir büyüklüktür. Bu büyüklük Isı Denklemi ile ifade edilir. Isı denklemi, ısı iletiminin ortaya çıkardığı bir formülasyondur. Isı denklemi en basit haliyle: r r J = −k∇T (3.14 ) r şeklinde verilir. Burada, J ısı iletimini ifade eder ve birimi joule.s/m ‘dir. k, ısı r iletkenlik katsayısını ifade eder. Isının aktığı yönde sıcaklık artar. Bu da J ile karakterize edilir. r ∫∫J .n r. s dθ∂ + ∫∫∫PK .∂v = ∫∫∫ .∂v (3.15) s v v dt 12 Pk gücü temsil eder, yani poynting vektörüdür. Denklem (3.13)’deki ifadeye hacim integrali uygulandığında ifade: ∂θ ∫∫∫ . ∂T ∂v = ∂t ∫∫∫ ρ .C p . .∂v (3.16)v v ∂t şeklini alır. ρ (yoğunluk) ve Cp (özgül ısı) basınca bağlı büyüklüklerdir. Eşitlik (3.15)’deki ifadede bulunan yüzey integralini hacim integraline çevirirsek, ifade: r r ∫∫J.n r.∂s = ∫∫∫(∇.J ).∂v = −∫∫∫∇.(k.∇T).∂v (3.17) S V V şeklini alır. Denklem (3.16)’da elde edilen hacim integrali, (3.14)’deki ısı denklemine uygulandığında, mikrodalga içerisinde bulunan malzemeye ait sıcaklık dağılımı: .(k. T (m,t)) ρ.C .∂T (m,t)∇ ∇ − p = −P(m,t) (3.18)∂t yukarıdaki ısı transferi denklemi tarafından verilir (Fu ve Metaxas, 1995). Bu ifade, mikrodalga ısıtma işlemlerinde yaygın olarak kullanılır. Burada dikkat edilmesi gereken nokta, P(m,t) kaynak tarafından verilen güçtür. P(m,t) pozitifse, dışarıdan malzemeye ısı verildiği anlaşılır. Negatifse, ısının kaybolduğu anlaşılır. TE10 modu için güç dağılımının, y yönünde düzgün bir biçimde dağıldığı ve iki boyutlu ısı transferi modelinin, x ve z yönünde mikrodalga ısıtma talimatlarını belirlemek için yeterli olacağı kabul edilir. Eğer iletim sıcaklığı sadece x ve z yönünde akarsa denklem (3.18)’ deki ifade aşağıdaki gibi elde edilir: ρ.C .∂T ∂ 2T ∂2k k Tp = + −P(m,t) (3.19)∂t ∂x2 ∂z 2 Mikrodalga gücü P(m,t), eşitlik (3.12)’deki ifadede w yerine 2πf konulmasıyla denklem (3.20)’deki biçime dönüşür. 13 P = π . f .ε .ε "0 r . E 2 (3.20) Bu güç, malzeme ısıtıldıkça değişir. Gücü hesaplayabilmek için, alan dağılımını bilmek gerekir. Malzeme içindeki gücün dağılımı, fırının geometrik yapısına bağlıdır. Malzemenin her tarafının, aynı miktarda ısınması için, malzemenin konulduğu yerdeki alanın düzgün dağılımlı olması gerekir. Bir malzemede, sıcaklık değişiminin bulunması için, ısı denklemiyle birlikte sınır koşullarının da sağlanması gerekir. Sınır koşulları, 1.ortam ile 2.ortam arasındaki ısı alışverişi ile bulunur. Yüzeydeki ısı alışverişi, 3 değişik işlemle gerçekleşir. Bunlar: 1. İletim ile ısı taşınması (Konduksiyon) 2. Işıma ile ısı iletimi (Radyasyon) 3. Taşınma ile ısı iletimi (Konveksiyon) 3.2.1 İletim İle Isı Taşınması (Konduksiyon) Kaynaktan çıkan ısı enerjisi, maddenin yüzeyindeki molekülleri titreştirir. Bu moleküller, titreşimi komşu moleküllere, onlar da diğerlerine şokla iletirler. Konduksiyon ile ısı iletimi çok yavaştır ve belli noktada statik dengeye ulaşır. Isı taşınımı gerçekleşirken, ısı belli bir oranda iletilir. Bu oran sabiti ısı iletkenliği adını alır. İletimle ısı taşınması, matematiksel olarak aşağıda verilen Fourier Yasası ile ifade edilir. r r J = −k∇T (3.21) r r J ısı iletimi, k ısı iletkenliği ve T ise sıcaklık dağılımıdır. Mikrodalga açısından düşünüldüğünde, dielektrik malzemenin ısı iletkenliğinin ölçümü zor olduğundan, böyle bir işlem pek kullanılmamaktadır. 14 3.2.2 Işıma İle Isı İletimi (Radyasyon) Herhangi bir sıcaklıktaki bir cisim termal radyasyonla enerji salar. Bu enerji elektromanyetik dalgalar vasıtasıyla taşınır. Isınmakta olan cismi çevreleyen yüzey bu salınan enerjiyi emer. Dolayısıyla cisim ve onu çevreleyen yüzey arasında bir ısı taşınması olur. Bu tip ısı taşınması daha çok yüksek sıcaklıklarda etkilidir. Düşük sıcaklıklarda ihmal edilebilir. Radyasyonla ısı iletimi tamamen kızıl ötesi bir ışımanın sonucudur. Mikrodalga sistemlerinde çok az kullanılan bir işlemdir. Eğer ışıma yoluyla bir alışveriş varsa, bu durum, Stephan- Boltzmann yasasına göre gerçekleşir. Çok yüksek sıcaklıklarda kullanılır. P = E.γ .T 4 (3.22) Burada, E malzemenin ısı yayma faktörü, γ stephan sabiti, P ise ısı gücüdür. 3.2.3 Taşınma ile Isı İletimi (Konveksiyon) Konveksiyon ile ısı taşınması iki mekanizma ile gerçekleşir. Konveksiyonla ısı taşınımı, rastgele moleküler hareketlerle taşınabileceği gibi makroskobik olarak toplu hareketle de taşınabilir. Bunlardan ilki cismin ısıtıldığı ortamdaki hava moleküllerinin cismin yüzeyine çarparak ısıyı taşıması; ikincisi ise bir moleküller topluluğu olarak hareket eden bir sıvı veya gaz akışının cisimden ısıyı taşımasıdır. Konveksiyonla ısı iletimi matematiksel olarak Newton Yasası ile ifade edilir: r nr.J = HC (T −T0 ) (3.23) Şeklinde verilir. Burada, nr yüzey birim vektörü, HC konveksiyon katsayısıdır. Mikrodalga ısıtma sistemlerinde en çok kullanılan işlem konveksiyonla ısı iletimidir. Isı transfer denkleminin çözümünde sınır koşulu olarak kullanılır. 15 3.3 Mikrodalga Isıtmanın Sayısal Modellenmesi Mikrodalga ısıtma problemlerini, karmaşık yapıları ve birbirleri ile ilişkili olarak gerçekleşen olaylar dizisini içermeleri sebebiyle analitik olarak çözümlemek çok güçtür. Dolayısıyla bu tür sistemlerin sayısal olarak çözümlenmesi, elde edilecek sonuçların değerlendirilerek sistemin tasarlanması açısından büyük önem taşır. Özellikle son yıllarda bilgisayar teknolojisindeki gelişmeler sayesinde sayısal çözümlerin hassasiyetlerinin arttırılmasıyla daha doğru sonuçlar elde edilmektedir. Bu kısımda kayıplı bir dielektrik malzeme ile yüklenmiş bir mikrodalga uygulayıcısı içindeki elektromanyetik alan dağılımının ve yük üzerindeki sıcaklık dağılımının bulunması amacıyla Maxwell denklemlerinin ve ısı transfer denklemlerinin FDTD- Sonlu Farklar Zaman Domeni Yöntemi ile çözümlerine değinilecektir. 3.3.1 FDTD Yöntemi (Sonlu Farklar Zaman Domeni Yöntemi) Sonlu farklar yöntemi sayısal analizde geniş kullanılma alanına sahiptir. Hatta sayısal işlemlerin pek çoğu için temel teşkil ederler. Bu yöntemde kapalı bir bölge kare ya da dikdörtgen şekilli elemanlara bölünür. Problemin koşullarının bilinmesi ve bölgenin kapalı olması gerekir. Düğüm noktalarına ait sonlu fark denklemleri yazılır. Bilinmeyen sayısı çoğaldıkça problemin çözümü zorlaşmaktadır. Sonlu farklar yöntemi kullanılarak çözülen alan hesaplarında yapılan ilk iş, problemin çözüleceği bölgeye uygun biçimde bir ızgara (grid) oluşturulmasıdır. Dolayısıyla çözümün doğruluğu, kullanılan ızgaranın göz boyutlarına, ızgaranın düzenine ve ızgaranın dış düğüm noktalarının sınırlarına uyumuna bağlıdır. Alan dağılımının bulunmasında, hesaplamanın yapıldığı bölgeye ilişkin sınır koşullarının verilmesi gereklidir. Sınır koşulları verilmezse bilinmeyen sayısı, denklem sayısından fazla olur, çözüme ulaşılamaz. Açık bölgelerin incelenmesinde, sistemin çözüm bölgesini etkilemeyecek şekilde hayali ve sınır koşulu verilen sınırlar ile kapatılması gerekir. 3.3.1.1 FDTD Yönteminde Türev İfadelerinin Taylor Serileri İle Analizi Diferansiyel denklemi verilen bölge için belirlenmiş sınır koşullarına bağlı, analitik olarak çözüldüğünde, sonuçta oluşan çözüm bölgesinde her noktada 16 diferansiyel denklemi tanımlanabilir. Problem analitik olarak çözümlenemediğinde veya sayısal hesaplama çok zor olduğunda analitik çözüm anlaşılmaz olur ve genellikle çözüm için sayısal tekniğe başvurulur. Sonlu fark yaklaşımı kullanıldığında, problem alanı ayrıklaştırılabilir. Böylece, sadece bölge civarındaki her nokta yerine, düğüm noktasının sonlu sayısında bilinmeyen bağımlı değişkenin değerleri dikkate alınır. Eğer N düğümleri seçilirse, N cebirsel denklemleri, diferansiyel eşitliklerinin ve problem için sınır koşullarının ayrıklaştırılması ile geliştirilir. Bu problem bölgesinde, olağan veya kısmi diferansiyel denklemlerinin çözümünün problemleri, uygun algoritmik çözümlerle bir takım cebirsel denklemlerin gelişmiş görevlerine dönüştürülür. Bu görünüşteki temel yaklaşım, cebirsel denklemlerin sonuçta oluşan düzenin yapısı parabolik, eliptik veya hiperbolik olan fiziksel sorun yönünden, kısmi diferansiyel denklemlerin karakterine bağlı olduğu olgu tarafından karmaşıklaşacaktır. Bundan başka, birçok ayrıklaştırma yöntemleri vardır, bu nedenle problemin yapısı için en uygun olanı hangisiyse onu seçmeliyiz. Taylor serilerinin formülasyonu için, türevin sonlu fark gösterimi; F(x,y) fonksiyon türevinin x=x0 ve y=y0 ‘da tanımlanması yoluyla aşağıdaki gibi ifade edilir. ∂F lim F(x0 + ∆x, y= 0 )− F(x0 , y0 ) (3.24) ∂x ∆X→0 ∆x Açık bir biçimde eğer fonksiyon F(x,y) sürekliyse, yukarıdaki (3.24) eşitliğinin sağ tarafı ∂F /∂x ’e makul yeterli bir küçüklükle yaklaşır ki, bu sonlu ∆x’ dir.Gelişmiş sonlu fark yaklaşımından türevler yerine bir temel form olan Taylor serileri kullanılabilir. Fonksiyon F(x)’in Taylor serisinin x0 noktasındaki açılımı göz önüne alınarak aşağıdaki ileri eşitlik (3.25) ve geri eşitlik (3.26) yönergeleri verilmiştir. 2 2 3 3 f (x0 +∆x) = f (x ) ∂f 0 + ( x) ∂ f (∆x) ∂ f (∆x) ∆ + + +... (3.25) ∂x ∂x2 30 0 2! ∂x 0 3! f (x x) f (x ∂f 2 2 3 3 0 − ∆ = 0 ) − (∆x) ∂ f (∆x) ∂ f (∆x) + 2 − 3 + ... (3.26)∂x 0 ∂x 0 2! ∂x 0 3! 17 Sonlu farklar yaklaşımı için bu iki temel formül df/dx’ in x0 noktasındaki birinci dereceden türevidir. Yeniden düzenlenmiş eşitlik (3.25) ve (3.26) ileri ve geri ilk türevi için sonlu farklar yaklaşımı, sırayla şu şekli alır: ∂f f (x0 +∆x)− f (∆x)= +0(∆x) (ileri) (3.27) ∂x 0 ∆x ∂f f (x0 )− f (x0 −∆x)= +0(∆x) (geri) (3.28) ∂x 0 ∆x 0(∆x) notasyonu, kırpma hatası (truncation error) olarak sonlu farklar yaklaşımı ile birleştirilerek tanımlanır. Bu notasyon türev ile sonlu fark gösterimi arasındaki farkı gösterir, geçişteki hatanın oranıdır. Verilen durum için şöyle yazılabilir: ∆x (∆x)20(∆x) = f "(x0 )+ f '' '(x )+ ... (3.29)2 6 0 eşitlik (3.25)’den (3.26)’ı çıkararak, merkezi farklar yaklaşımı şu şekilde belirlenebilir: ∂f f (x = 0 +∆x)− f (x0 −∆x) −0(∆x)2 (merkezi) (3.30) ∂x 0 2∆x (∆x) 2 50(∆x) 2 = f ' ' ' (x ) (∆x)0 + f ' ' ' ' ' (x0 )... (3.31)6 120 Kırpma hatası (truncation error) dikkatlice incelenirse merkezi farklar yaklaşması ileri ve geri farklardan daha kesin yaklaşmalarla hesaplanır. Merkezi farkların ileri ve geri fark yöntemlerine göre daha yaklaşıklıkla türev ifadelerine yaklaşım yapmaktadır. Çünkü 0(∆x)2 olduğu için hassasiyet artmaktadır. 3.3.1.2 Birinci Dereceden Türevin Sonlu Farklar Yaklaşımı x0’da örgü noktasında i seçilir. Sonra notasyon i+1 ve i-1’e başvurularak, sırasıyla örgü noktaları x0 +∆x ve x0 –∆x alınırsa, benzer biçimde notasyon i +2 ve 18 i –2’e başvurarak sırasıyla örgü noktaları x0 +2∆x ve x0 –2∆x, şeklinde yazılır. Bu notasyonu kullanarak, birinci dereceden türev için iki-üç ve dört-nokta formülasyonu sunulabilir. İki –Nokta Formülasyonu: (ileri) (3.32.a) (geri) (3.32.b) (Merkezi) (3.32.c) Burada h= ∆x olursa yukarıdaki üç formül tek eşitlik halinde yazılabilir: (3.33) Burada; İleri için; Merkezi için; Geri için; Üç-Nokta Formülasyonu: (3.34.a) (3.34.b) Dört nokta formülasyonu (3.35.a) (3.35.b) (3.35.c) 19 Yaklaşımın doğruluğunu geliştirmek için, sınırda bulunan i noktasındaki düğümün birinci dereceden türevini göstermede, sınırın bir yanındaki ikiden fazla örgü noktalarını kullanarak yukarıda verilen üç ya da dört-nokta formülasyonu kullanmak yararlıdır.(Özışık,1994) 3.3.1.3 İkinci Dereceden Türevin Sonlu Farklar Yaklaşımı İkinci dereceden türevin sonlu farklar yaklaşımını geliştirmek için eşitlik (3.25) ve (3.26) ile verilen Taylor Serileri açılımlarından faydalanılır. İkinci dereceden türevin merkez sonlu farklar yaklaşımını elde etmek için, Eşitlik (3.25) ve (3.26) toplanır, sonuçtaki ifade ( d2f/dx2)’den çözülür ve sonuçta kısaltılmış notasyon şu şekilde yazılabilir: (merkezi) (3.36) İkinci dereceden türevin ileri ve geri sonlu-farklar yaklaşımını geliştirmek için, fonksiyon f(x0 + 2∆x) ve f(x0 - 2∆x) Taylor serilerine açılabilir. Fonksiyon f’(0), f(x0 - 2∆x) ve eşitlik (3.25) arasından çıkarılabilir ve sonuç ifadesi (d2f/dx2)’den çözülebilir. İleri sonlu-farklar yaklaşması ikinci dereceden türevler için şöyle belirlenir: (İleri) (3.37) Benzer biçimde, Fonksiyon f’(0), f(x0 - 2∆x) açılımı ve eşitlik(3.26) de verilen açılım arasından çıkarılabilir ve sonuç ifadesi ( d2f/dx2)’den çözülebilir. İkinci dereceden türevin geri sonlu-farklar yaklaşımı şöyle belirlenir: (geri) (3.38) 20 Yukarıda verilen ikinci dereceden türev için sonlu-farklar yaklaşımı üç örgü noktasını kullanır. Üç noktadan daha fazlası için aşağıda bazı yaklaşımlar verilmiştir(Özışık,1994) (3.39) 3.4 Elektromanyetik Alanların FDTD Yöntemi ile Analizi Elektromanyetik Teorem için Maxwell Denklemlerine Sonlu Fark (3.40) Zaman Domen Yöntemi uygulanır. Elektromanyetikte en çok kullanılan sayısal yöntemlerin başında gelir. Yöntemun aslı, Maxwell denklemlerinde alan büyüklüklerinin zamana ve konuma bağlı türevleri yerine, bu büyüklüklerin Taylor serisine açılmasıyla elde edilen fark denklemlerinin kullanılmasıyla, sürekli olan elektromanyetik alanların sonlu hacimlerle ve eş zaman aralıklarıyla örneklenmesine dayanır. Neticede türev ifadelerine sahip sürekli denklem, ayrık fark denklemlerine dönüşür. Bu fark denklemlerinin çözümü ile elde edilen değerler, sürekli alanın ilgili zaman ve konumdaki değerine karşılık gelir. 3.4.1 FDTD Yönteminin Uygulanmasında Temel İlkeler Maxwell Denklemlerinin belli bir ortamda FDTD kullanarak çözebilmek için, ortamın birim hücrelere bölünmesi gerekir. Bu birim hücrelerin x- yönünde genişlikleri ∆x, y- yönünde genişlikleri ∆y ve z- yönünde genişlikleri ∆z’ dir. x, y, z eksenlerinde kaç adıma karşılık geldiğini belirtir. Birim hücre, koordinat sisteminde orijine indirgendiğinde, orijin noktasına denk gelen köşeyi numaralayan üç indis ile belirtilir. Bu indisler, Yee tarafından önerilen birim hücre Şekil 3.3’de verilmektedir. 21 Şekil 3.3: Yee Birim Hücresi Yee, alan bileşenlerini farklı noktalara yerleştirerek programlama açısından bazı zorlukları ortadan kaldırmıştır. Yee hücresinde, elektrik alanın x- bileşeni, x- yönünde uzanan kenarın, y- bileşeni, y- yönünde uzana kenarın, z- bileşeni, z- yönünde uzanan kenarın orta noktasına yerleştirilmiştir. Magnetik alanın x- bileşeni, y-z düzlemindeki yüzeyin, y-bileşeni, x-z düzlemindeki yüzeyin ve z-bileşeni, x-y düzlemindeki yüzeyin orta noktasına yerleştirilmiştir. Yee notasyonuna göre (i,j,k) noktası, FDTD ağında (i∆x, j∆y, k∆z) koordinatlarıyla belirtilen konuma yerleştirilmiştir. (∆x, ∆y, ∆z) üç yöndeki konum adımını gösterir. Zaman ve konumda değişen herhangi bir fonksiyon Fn(i, j, k)= F(i∆x, j∆y, k∆z, n∆t) (3.41) 22 şeklinde ifade edilir. Elektromanyetik alan büyüklükleri ∆t adımlarıyla örneklenirler. Elektrik alanları ∆t’ nin tam katlarında, manyetik alanlar ∆t’ nin yarım katlarında hesaplanırlar. Maxwell denklemlerinin hesaplanmasında kullanılan Yee birim hücresi, sonlu farklar yönteminde (FDTD) temel alınır. Bu yöntemin genel felsefesi, boşluk örgüsü ve eşit uzaklıktaki noktalar arasında belirgin her iki noktada, boşluğun sonlu durumunda sürekli alanın elektrik ve manyetik alanlarını örneklemeyle gösterilebilmesidir. Dalga yayılması, her bir örgü noktasına Maxwell denklemlerinin sonlu farklar yöntemi tarafından defalarca uygulanır. FDTD yönteminde, önemli olan Maxwell denkleminin alan ıraksaklık durumunu yaklaşık bir şekilde gösterebilmesidir. Mikrodalga sorunlarını, sayısal olarak çözmek genelde zordur. Çünkü birçok durumda sadece yalın geometrileri kolayca hesaplayabiliriz. Kapalı mikrodalga sistemi gibi, karmaşık üç boyutlu yapıda elektromanyetik alan probleminin çözümü zordur. FDTD ağının her düğümünde sinüsoidal genliğin ve elektromanyetik alanın aşamaları hesaplanır. Yöntemin çok önemli avantajlarından ikisi, kolayca uygulanılması ve hiçbir matris tersi içermemesidir. Yöntem bir, iki ve üç boyutlu problemler için aynı zamanda kolayca ayarlanabilir. Bu yöntem, iki ve üç boyutlu Maxwell denklemlerini çözmede kullanılan mikroşerit, boşluk ve dalga kılavuzu problemlerine uygulanabilir. Bu yöntemin riskli bir tarafı hesaplamalarda iteratif yani basamak basamak çözüme yaklaşmasıdır. Bu yüzden kesin sonuçtan uzaklaşılabiliyor. Yani (n+1)’ inci hesap için n’ inci hesap sonucunun bulunması gerekiyor. Sonuçta bir adım önce hesaplanmış olan veriler bir adım sonraki hesaplamada kullanılmaktadır. 3.4.2 FDTD Yönteminin Maxwell Denklemlerine uygulanması Bu uygulamada Yee notasyonlarını takip ederek (i, j, k) noktalarında FDTD ağının 3 kartezyen koordinatı olan (i∆x, j∆y, k∆z) ile sınırlarız. Burada i, j, k doğal tam sayılardır. Konum arttırımı bu üç yönde ∆ olarak arttıracağız. Böylece ifademiz eşitlik (3.41)’deki Fn(i, j, k)=F(i∆x, j∆y, k∆z, n∆t) gibi olacaktır. Elektrik ve manyetik alanların hesaplanış şekli olarak şunları söyleyebiliriz; elektrik alanlar zaman adımlarının tam katlarında t=(T+n) manyetik alanlar yarım katlarında t=(T+n+0.5) hesaplanır. Gerekli açılımların yaparak elektrik ve manyetik alanları eşitlik denklemlere göre denklemleri sırasıyla yazarsak, manyetik alanlar eşitlik (3.11)’deki ifadeleri, 23 elektrik alan için; eşitlik (3.12)’deki ifadeler kullanılarak alan değerlerini hesaplayabiliriz. (3.11.a)’daki denklemi ele alalım. Bu denklemde zaman ve konuma göre türevlere merkezi farkları uygulayalım. Bu işlemi yaparken Yee hücresini göz önünde tutacağız. Yee hücresine göre Hx manyetik alan bileşeni (i, j+0.5, k+0.5) noktasında tanımlıdır. ∂Hx H n+1/ 2 x (i, j +0.5,k +0.5)−H n−1/ 2 x (i, j +0.5,k +0.5)= (3.42) ∂t ∆t Türevin ayrıklaştırılmasını eşitlik (3.11.a)’da tanımlayarak, zamanda merkezi farklar, elektrik alanda konumda ileri farklar izlenir. Buna göre: H n+1 / 2x (i, j + 0.5, k + 0.5) − H n−1 / 2 x (i, j + 0.5, k + 0.5)µ 0 =∆t E n (i, j + 0.5, k + 1) − E n (i, j + 0.5, k ) y y (3.43) ∆z E n (i, j + 1, k + 0.5) − E n (i, j, k + 0.5) − z z ∆y H n+1/ 2x (i, j + 0.5, k + 0.5) = H n−1/ 2 x (i, j + 0.5, k + 0.5) ∆t + [E n (i, j + 0.5, k +1) − E n (i, j + 0.5, k ) (3.44) µ 0∆z y y ∆t − [E nZ (i, j +1, k + 0.5) − E n Z (i, j, k + 0.5)]µ 0∆y R ∆t ∆t ∆tx = , Ry = , Rz = (3.45)µ0∆x µ0∆y µ0∆z Olmak üzere (3.11.a) denkleminin ayrıklaştırılmasında izlenen yol (3.11.b) ve (3.11.c) denklemlerine uygulanırsa: 24 H n+1 / 2 (i, j + 0.5, k + 0.5) = H n−1 / 2x x (i, j + 0.5, k + 0.5) +R nz [E (i, j + 0.5, k + 1) − E n (i, j + 0.5, k )] (3.46 ) y y −R y [ E n Z (i, j + 1, k + 0.5) − E n Z (i, j , k + 0.5)] H n+1 / 2 (i + 0.5, j, k + 0.5) = H n−1 / 2Y Y (i + 0.5, j, k + 0.5) +R X [E n (i +1, j, k + 0.5) − E n (i, j, k + 0.5)] (3.47) Z Z −R [E nZ X (i + 0.5, j, k +1) − E n X (i + 0.5, j, k )] H n+1/ 2Z (i + 0.5, j + 0.5, k ) = H n−1/ 2 Z (i + 0.5, j + 0.5, k ) +R [E nY X (i + 0.5, j +1, k ) − E n X (i + 0.5, j, k )] (3.48) −Rx [E n Y (i +1, j + 0.5, k ) − E n Y (i, j + 0.5, k )] Böylece üç manyetik alan bileşenine ait diferansiyel denklemler ayrıklaştırılmış oldu. Yee hücresinde Ex alan bileşeni (i+0.5, j, k) noktasında tanımlıdır. Elektrik alan bileşenlerine ait diferansiyel denklemlerinin ayrıklaştırılması için (3.12.a) denkleminden başlayalım. Bu denklemde zamana ve konuma göre türevlere merkezi farklar uygulayalım: ∂E En+1/2x x (i+0.5, j,k)−E n−1/2 x (i+0.5, j,k)= (3.49) ∂t ∆t E n+1 / 2x (i + 0.5, j, k ) − E n−1 / 2 x (i + 0.5, j, k )ε + σ .E nX (i + 0.5, j, k ) =∆t H nZ (i + 0.5, j + 0.5, k ) − H n Z (i + 0.5, j − 0.5, k ) (3.50) ∆y H n (i + 0.5, j, k + 0.5) − H n (i + 0.5, j, k − 0.5) − Y Y ∆z Buradan sonra yapacağımız işlemler, manyetik alanlar için yaptığımızdan farklı olacak. FDTD ile Maxwell denklemlerinin çözümünde elektrik alanlar ∆t’ nin katlarında, manyetik alanlar ise ∆t(n+0.5) zamanlarında hesaplanır. Bu özellikle bilgisayar 25 ortamında yapılacak işlemler için kolaylık taşır. Bu sebeple elektrik alanın ayrıklaştırdığımız (3.50) denkleminde zaman adımını (n) 0.5 kaydırmamız gereklidir. Bu durumda manyetik alanlar (n+0.5) zaman adımına kayacaktır. E n+1x (i + 0.5, j, k) − E n x (i + 0.5, j, k)ε +σ .E n+0.5 (i + 0.5, j, k) = ∆t X H n+1/ 2Z (i + 0.5, j + 0.5, k) − H n+1/ 2 Z (i + 0.5, j − 0.5, k) (3.51) ∆y H n+1/ 2Y (i + 0.5, j, k + 0.5) − H n+1/ 2 Y (i + 0.5, j, k − 0.5)− ∆z (3.51) denklemine bakıldığında E n+0.5x (..) terimi problem yaratmaktadır. Çünkü elektrik alanlar ∆t’nin tam katlarında hesaplanır. Alanların sürekliliğini göz önüne alarak bu terim için, çevresindeki iki noktadaki değerlerin ortalaması alınır. n 1+ E n+1 nE 2 (i, j, k ) X (i, j, k ) + E X (i, j, k )X = (3.52 )2 Bu durumda elektrik alanın (3.51) denkleminin ayrıklaşmış hali şöyle yazılabilir: ε [En+1(i+0.5, j,k)−En(i+0.5, j,k)] σ+ [En+1 nx x x (i+0.5, j,k)−Ex (i+0.5, j,k)]=∆t 2 1 [Hn+1/ 2(i+0.5, j+0.5,k)−Hn+1/ 2Z Z (i+0.5, j−0.5,k)] (3.53)∆y 1 − [Hn+1/ 2(i+0.5, j,k+0.5)−Hn+1/ 2Y Y (i+0.5, j,k−0.5)]∆z Elektrik alanları n ve n+1 olarak parantezine alırsak: ⎛ ε σ ⎜ + ⎞ ⎟En+1x (i+0.5, j,k) ε σ −⎛ − ⎞⎜ ⎟Enx (i+0.5, j,k) = ⎝∆t 2 ⎠ ⎝∆t 2 ⎠ 1 [Hn+1/2(i+0.5, j+0.5,k)−Hn+1/2Z Z (i+0.5, j−0.5,k)] (3.54)∆y 1 − [Hn+1/ 2Y (i+0.5, j,k+0.5)−Hn+1/ 2Y (i+0.5, j,k−0.5)]∆z 26 n zaman adımındaki bileşeni eşitliğin diğer tarafına alırsak: 2ε +σ∆t En+1(i 0.5, j,k) 2ε −σ∆tx + = E n (i+0.5, j,k) 2∆t 2∆t x 1 [Hn+1/ 2Z (i +0.5, j +0.5,k)−Hn+1/ 2(i +0.5, j −0.5,k)] (3.55)∆y Z 1 − [Hn+1/ 2(i+0.5, j,k +0.5)−Hn+1/ 2Y Y (i+0.5, j,k −0.5)]∆z En+1(i 0.5, j,k) 2ε −σ∆t+ = Enx (i +0.5, j,k)2ε +2σ∆t x 2∆t [Hn+1/ 2Z (i+0.5, j +0.5,k)−Hn+1/ 2Z (i +0.5, j −0.5,k)] (3.56)(2ε +2σ∆t)∆y 2∆t − [Hn+1/ 2(i+0.5, j,k +0.5)−Hn+1/ 2(i+0.5, j,k −0.5)] (2ε +2σ∆t)∆z Y Y Katsayılara değişkenleri atarsak; C 2ε (i, j,k) −σ (i, j,k)∆ta = (3.57)2ε (i, j,k) +σ (i, j,k)∆t C 2∆tbx = (3.58.a)(2ε(i, j,k)+2σ(i, j,k)∆t)∆x C 2∆tby = (3.58.b)(2ε (i, j,k) + 2σ (i, j,k)∆t)∆y C 2∆tbz = (3.58.c)(2ε (i, j, k ) + 2σ (i, j, k )∆t)∆z (3.56)’daki eşitliğe bu katsayılar atanırsa: En+1x (i+0.5, j,k) =Ca(i+0.5, j,k)E n x (i+0.5, j,k) + Cby(i+0.5, j,k)[Hn+1/ 2Z (i+0.5, j+0.5,k)−Hn+1/ 2Z (i+0.5, j−0.5,k)] (3.59) −Cbz(i+0.5, j,k)[Hn+1/ 2Y (i+0.5, j,k+0.5)−Hn+1/ 2Y (i+0.5, j,k −0.5)] 27 Diğer elektrik alan bileşenlerine ait diferansiyel denklemlerde aynı yolla ayrıklaştırılırsa: En+1y (i, j +0.5,k) =C (i, j +0.5,k)E n a y (i, j +0.5,k) +Cbz(i, j +0.5,k)[Hn+1/ 2x (i, j+0.5,k +0.5)−Hn+1/ 2x (i, j+0.5,k −0.5)] (3.60) −C (i, j +0.5,k)[Hn+1/ 2(i+0.5, j +0.5,k)−Hn+1/ 2bx Z Z (i−0.5, j +0.5,k)] En+1Z (i, j,k +0.5) =Ca (i, j,k +0.5)E n Z (i, j,k +0.5) +C n+1/ 2bx(i, j,k +0.5)[Hy (i +0.5, j,k +0.5)−H n+1/ 2y (i −0.5, j,k +0.5)] (3.61) −C (i, j,k +0.5)[H n+1/ 2 (i, j +0.5,k +0.5)−H n+1/ 2by x x (i, j −0.5,k +0.5)] Elektrik alanlar da bulunmuş olur 3.4.3 FDTD Yönteminde Sınır Koşulları FDTD örgüsünün kenarlarındaki alanları hesaplarken problemle karşılaşmaktayız. Örgü kenarındaki yansımalar, örgü içindeki alanları bozacaktır. Biz böyle bir yansımanın gözlenilen alana ulaşmadan önce zaman-adımlı prosedürü sonlandırabiliriz ya da örgüyü çok çok geniş yapabiliriz. Fakat bunlar hesaplama açısından uygulanabilir yollar değildir. Bu yüzden örgü kenarında alan bileşenlerini güncellemede karşılaşılan soruna ekstra özen verilmelidir. Örgü kenarında karşılaşılan bu güncelleme sorunun en pratik çözümü sınır koşullarını kullanmaktır. Bir boyutlu durumda, istenilen koşul basit ve kesindir. Çünkü düzlemsel dalga örgü kenarına normal olarak gelmektedir. Bu yüzden, basit propagasyon gecikmesi kullanılabilir. İki ve üç boyutlu durumda problem çok zordur. Çünkü örgü kenarında dalga normal olarak gelmemektedir. Son on yılda birçok sınır koşulu geliştirilmiştir. Bunlardan bazıları MUR, PML ve DBC’dir. MUR Gerrit MUR adlı kişi tarafından geliştirilen, ismiyle anılan bir yöntemdir. MUR için, birinci ve ikinci dereceden olmak üzere iki tip yöntem vardır. 28 PML(Perfectly Matching Layer) ve DBC(Dispersive Boundary Condition) adlı üç yöntem bize hesaplamalarımızda büyük kolaylıklar sağlar. Bu üç yöntemin ortak amaçları ve işlem yaptıkları noktayı araştıracak olursak, bilgisayarda benzetimde sonsuz uzun alamadığımız boru için borunun sonuna veya başına konulan yutucu yüzeyin, işlemlerimizde borunun içine dağılan alanın tespitinde yardımcı olan çok önemli yöntemlerdir. Yutucu yüzeye gelen dalga geri yansımayacağı için, bu yansımama olayı aldığımız verilerle bir genelleme yapmamızı sağlıyor. Yutucu yüzeye gelen bir dalganın yansımaksızın yutuluşu Şekil 3.4’de görülmektedir. Şekil 3.4: Yutucu maddeye gelen dalganın durumu Bu üç yöntem aynı işlemi yaptıkları halde aralarında da bazı farklar vardır. MUR yöntemi bize en uygun olanıdır; çünkü MUR sadece ilk hücrede işlem yapar. Bu da demektir ki, bizim işlem yapacağımız dalga kılavuzu boyutlarında ilk hücrenin yutucu kesit olarak alınmasıdır. PML ve DBC yöntemleri ise 7 veya 8 hücreyi kapsayan bir uzunluktaki hücrede yutucu yüzey olarak işlem yapmaktadır ki bu da bizim kesiti daraltmakta veya işlem yapacağımız dalga kılavuzunun boyutlarını arttırmamıza neden olmaktadır. PML ve DBC yöntemlerinin bilgisayar işlem süreleri de fazlaca sürmektedir. Yani işlem yoğunluğu fazla yöntemlerdir. MUR bu konuda bize epey avantajlar sunmaktadır. Çünkü işlem süresini oldukça fark edilir düzeyde azaltmaktadır. Yalnız MUR’daki verilerimiz diğer iki yönteme göre daha az hassasiyet göstermektedir. Hata payımız artmaktadır. Özellikle kesim frekanslarında hassasiyet minimuma inmektedir.(Kunz,1993) 29 3.4.3.1 Birinci dereceden MUR yutucu sınır koşulu Yönteminin Matematiksel ifadesi FDTD yöntemi yoluyla birleştirilmiş Maxwell denklemlerinin çözümü, herhangi bir alan bileşeni için dalga denkleminin çözümüne eşittir. MUR’un türetimi daha çok genelleştirme yöntemlerinde olduğu gibi dalga denklemi ile başlar. Özel alan bileşenleri E ya da H ile işlem yapar. Tek alan bileşen W için MUR’ un çözümünü izleyerek dalga denklemi olarak yazılabilir. ( ∂ x2 + ∂ y2 + ∂ z2 – c -20 ∂ t2 )W=0 (3.62) ∂ x2=∂ /∂ x2 ’ dir ve c0 yayılım hızıdır. FDTD problem boşluğu (uzayı) veya ağı x≥0 da konumlandırılmıştır. Burada x=0’ da sınır vardır ki, bu dışsal sınır koşulları zorlanmadıkça dağılmış dalga ulaşacaktır ve yansıtır durumda olacaktır. Bu dışsal sınır koşulları yukarıdaki dalga denkleminden bulunur. Dağılmış dalga, hız bileşenleri vx , vy , vz ‘ye sahip olur, öyle ki v 2+v 2x y +v 2z =c 20 olur. Aynı zamanda Sx, Sy, Sz ters hız bileşenlerini de tanımlamak mümkündür.(Kunz,1993) S ∂x ∂ /∂x 1 ∂y ∂zX = = = SY = S∂t ∂ /∂t ∂x /∂t ∂t Z = (3.63) ∂t S 2 + S 2 S 2 1X Y + Z = 2 (3.64)c0 Dağılmış dalga, düzlemsel dalgalar tarafından yerel olarak yakınlaşmaktadır. W=Re[ψ(t+Sxx+Syy+Sz.z)] (3.65) ve Sx’i (1/c 20 -S 2y -S 2)1/2z olarak yerine koyarsa dalga denklemi şu hale gelecektir: ⎧ 1 ⎫ W Re ψ ⎡= ⎨ ⎢t + (c2 2 2 ⎤20 − SY − SZ ) x+ SY y + SZ z⎥⎬ (3.66) ⎩ ⎣ ⎦⎭ 30 birinci dereceden sınır koşulu x=0 sınırında denklem: ⎧ ∂ 1 −c−1⎨ 0 [1−(c 20SY ) − (c0SZ )2 ]2 ∂ ⎫⎬WX=0 = 0 (3.67) ⎩∂x ∂t⎭ eşitlik (3.66) ve(3.67) eşitliklerinde görülen W, (-x) yönünde giden dalgadır. Sy ve Sz çözümü için üstteki birinci derece sınır koşuluna bakarak karar verilmez. Birinci ve ikinci dereceden yaklaşımlarını takip ederek incelenirse bu belirsizliği ortadan kaldırabilir ve x=0 sınırında bulabilmek için W’ ya izin verir. Birinci dereceden yaklaşımı: [ 11− (c0S )2Y − (c0S 2 2Z ) ] =1+ 0[(c0SY )2 − (c S )20 Z ] .1 (3.68) Böylece birinci dereceden sınır koşuluna, birinci dereceden aşağıdaki eşitlik tarafından yaklaşılmaktadır: ⎛ ⎜ ∂ 1 ∂ ⎞⎜ − ⎟∂x c ⎟ WX =0 = 0 (3.69) ⎝ 0 ∂t ⎠ ikinci dereceden yaklaşımla birinci dereceden sınır koşulu: ⎡ ∂ 1 ⎧ 1 ⎤ ⎢ − ⎨1− [(c S )2 ∂0 Y −(c S )2 ]⎫ W∂x c 2 0 Z ⎬∂t ⎥ X=0 = 0 (3.70)⎣ 0 ⎩ ⎭ ⎦ veya zaman domenindeki türev eşitliğini alırsak: ∂ 1 ∂2 1 ⎧ − + c ⎡ ∂ 2 ∂ 2⎤ ∂ ⎫2 2 0⎨⎢ (SY ) − (SZ ) ⎥ ⎬WX=0 = 0 (3.71)∂x c0 ∂t 2 ⎩⎣ ∂t ∂t ⎦ ∂t⎭ FDTD ağının dışsal çevre uzunluğunda yayılma durumu, sonsuzluğun hesaplama bölgesine uzanmasıyla benzetilebilir. Yansıyan dalga veya dağılmış dalga yutucu yüzeye uygulanarak emilebilir. Açık ortamlarda veya yarı açık ortamlarda FDTD’ yi kullanarak Maxwell denklemlerinin çözümü, programlama açısından çok 31 büyük hafıza gerektirir. Bu hacimde açık kısımlara yutucu yüzey denilen sınır koşulları uygulanır. Yutucu yüzeyler, üzerine gelen dalgayı geri yansıtmazlar yani dalga sonsuza doğru ilerlemektedir. Yutucu sınır koşulu, sınır düzleminde elektromanyetik alanın TE10 modunun faz hızı ile z koordinat azalması yönünde dalga yayılımının içerdiği olguyu ifade eder. Bu tanımlamalardan sonra yutucu yüzey z=0’da FDTD’ye MUR’ u uygularsak: Şekil 3.5 MUR yönteminde açık uçlara yutucu yüzey yerleştirilmesi Kaynak tarafında z=0 noktasında yutucu sınır koşulları uygulanırken, ⎛ ⎞ ⎜ ∂ 1 ∂⎜ −∂z c ∂t ⎟ ⎟E x = 0 (3.72 .a)z=0 ⎝ 0 ⎠ ⎛ ⎜ ∂ 1 ∂ ⎞⎜ − ⎟⎟E y = 0 (3.72 .b) ⎝ ∂z c z=00 ∂t ⎠ c0 TE10 modunun faz hızıdır. Eşitlik (3.72.a)’ i açarsak, En+1x (i +0.5, j,0) c ∆t −∆x = En (i +0.5, j,1)+ 0x [En+1x (i +0.5, j,1)−Enx (i +0.5, j,0)] (3.73)c0∆t +∆x Eşitliğin sağ tarafındaki (n+1)’inci E alan ifadesinin n zaman cinsinden eşitlik (3.72.a) yerine koyarak çözümü kolaylaştırmış oluruz. Aynı işlemi eşitlik (3.72.b)’ da yaparsak: 32 En+1(i, j 0.5,0) En (i, j 0.5,1) c0∆t −∆yy + = y + + [En+1y (i, j +0.5,1)−Eny (i +0.5, j,0)] (3.74)c0∆t +∆y bulunmuş olur. 3.4.3.2 Dielektrik yüzey üzerinde sınır koşullarının FDTD Yöntemi ile Analizi Bir mikrodalga ısıtma sisteminde, Maxwell denklemlerinin sayısal çözümünde sınır şartlarının da çözüme ilave edilmesi gerekir. Isıtma sisteminin iletken duvarları üzerindeki elektrik alan, r rn×E = 0 (3.75) s şartını sağlamalıdır. Bu, elektrik alanın iletken duvarlar üzerindeki teğetsel bileşenin sıfır olması anlamına gelir. Çözüm esnasında bu duvarlar üzerinde yer alan elektrik alanın teğetsel alan bileşenleri doğrudan sıfıra eşitlenerek sınır şartı sağlanır. Mikrodalga ısıtma sistemi dielektrik malzeme içerdiğinden dolayı dielektrik malzemenin yüzeyindeki sınır şartında, iki farklı dielektrik ortam arasındaki yüzeyde elektrik alan, r rn .(ε r r 1 E1 ) = n .(ε E ) (3.76 )s 2 2 şartını sağlamalıdır. Bu şart, iki farklı dielektrik ortam arasındaki yüzeyde elektrik elektrik alanın normal bileşeninin süreksiz, elektrik alanın teğetsel alanının dielektrik yüzeye teğet olan bileşenleri süreklidir. FDTD formülasyonunda elektrik alanların teğetsellik şartı kullanılacaktır. Bunun için, Şekil 3.6’da verilen iki dielektrik ortam arasındaki sınır şartını elde edelim. Bu düzleme, elektrik alanın x ve y bileşenleri teğettir ve süreklidir.(Pozar,D.M) 33 Şekil 3.6: İki dielektrik ortam arasındaki sınır şartı Birinci ortamdaki elektrik alanın x-bileşeni ve ikinci ortamdaki elektrik alanın x- bileşeni önceki kısımlarda verilen eşitlik (3.12.a)’daki elektrik alanın x-bileşenine ait diferansiyel denklemi ayrı ayrı sağlar. ∂E ε x1 σ E ∂H+ = z ∂H 1 1 x1 − y (3.77) ∂t ∂y ∂z ∂Ex2 ∂Hz ∂Hε y2 +σ∂t 2 Ex2 = −∂y ∂z Her iki denklemi toplar ve ortalamasını alırsak ve elektrik alanın teğetsel bileşenin sınır üzerinde sürekli olduğunu göz önüne alırsak: ε1 +ε 2 ∂Ez σ1 +σ 2 ∂HE y ∂H+ = − x (3.78.a) 2 ∂t 2 z ∂x ∂y ε1 +ε2 ∂Ex σ1 +σ2 E ∂H ∂H + z y 2 ∂t 2 x = − (3.78.b) ∂y ∂z FDTD yönteminde, yalıtkanın geometrisi, katsayılar Ca ve Cb’ i tanımlayarak kolaylıkla modellenebilir. Ancak, iki dielektrik sınırı arasında, teğet halindeki elektrik alan bileşenleri için ε ve σ denklemlerde kullanılmış olması gerektiğinden gerekli sınır tanımlamaları yapılmalıdır. 34 Katsayıları da bu sınır koşuluna göre yeniden yazarsak, 1 σ 1 +σ− 2 R Ca (i + 0.5, j, k ) ε = 1 + ε 2 (3.79) 1 σ+ 1 +σ 2 R ε1 + ε 2 R Cb(i+0.5, j,k) = a (3.80) ε1 +ε2 σ1 +σ+ 2 R 2 2 R (c= 0∆t )2 R ∆ta = (3.81)∆x 2ε0 Ey ve Ez, yalıtkan ara birime teğet halindeki bileşenler olduğunda benzer denklemler bulunur. Normal elektrik alan bileşenleri ve manyetik alan bileşenleri için denklemler, orijinal biçimleriyle kullanılmazlar ve düzeltilmezler. Aynı şekilde FDTD yöntemine uygulayarak eşitlik (3.74)’e benzeyecek biçimde eşitliğimizi yazabiliriz. Kenar ve köşelerde sınır koşullarına göre hesaplama yaparken alan bileşenlerin sabit ε ve σ değerlerinde Şekil 3.7’den de görüldüğü gibi bu alanın farklı sabitlere sahip hücrelerle komşu olması neticesinde, bu sabitlerin toplamının dörde bölümüyle gerekli ε ve σ değerleri hesaplanabilir. Karşılaştırmalı olarak yapılan bu işlem komşu dört hücrelerin ε ve σ değerleri farklıysa bu yola başvurulur. Şekil 3.7: Farklı dielektrik değerlere sahip hücrelerde E alanı 35 ε ε1 +ε 2 +ε 3 +ε 4T = (3.82)4 σ σ1 +σ 2 +σ 3 +σ 4T = (3.83)4 şeklinde alınarak ifade edilebilir. 3.4.4 Kararlılık İçin Hücre Boyutunun Belirlenmesi Arzu edilen kesin sonuçlar için temel kesitleme hücre boyutunun en küçük dalga boyundan daha küçük olmuş olması gerektiğidir. Açık soru “Ne kadar daha az? ”dır ve bu soruya eklenmesi gereken bir diğer soru “ Doğru sonucun nasıl olmasını istiyorsunuz?” dur. Her dalga boyu başına 10 hücre olmalıdır. Bu demektir ki, her hücrenin kenarı 1/10 λ olmalı veya ilgili en yüksek frekanstan (en küçük dalga boyu) daha az olmalıdır. Kabul edilebilir bir sonuç için en az her dalga boyu için dört hücre gerekir. Eğer hücre boyutu bu Nyquist örnek sınırından daha küçükse, λ=2∆x makul sonuçlar ele etmek için önemli olan nyquist sınırı üstünde sinyal bileşenleridir. Hesaplanan sayısal boşluğun bazı kısmı nüfuz edici malzeme ile doldurulursa maksimum hücre boyutunu belirlemek için malzeme için de dalga boyunu kullanmamız gerekir. Eğer düzgün dağılımlı (tek biçimli) hücre boyutu bu gücü başından sonuna kadar kullanırsa problem boşluğundaki hücre göreceli olarak küçük olacaktır. Bu hücre sayılarının sayısını çokça arttırmak gerekecektir. Herhangi bir dikkate değer zaman adımında FDTD örgüsünün alan dağılımının ayrık uzaya ait örneği olduğunu hesaba kattığımızdan dolayı hücre boyutunun neden bir dalga boyundan çok daha küçük olmuş olması gerektiğini anlamış olmalıyız. Nyquist örnekleme teoreminden, uzaya ait bilginin uygun bir şekilde örneklenmesi için, dalga boyu başına en az iki örneklem gereklidir. Çünkü bizim örneklememiz kesin değildir ve en küçük dalga boyumuz tam olarak belirli değildir. Dalga boyu başına ikiden fazla örneklem gereklidir. Diğer bir ilintili düşünce, örgü dağılma hatası (grid dispersion error) dır. FDTD deki doğal yaklaşımlardan, farklı frekanslardaki dalgalar ağ boyunca biraz farklı hızlarda yayılacaktır. Yayılma hızındaki bu fark, ağa ilişkin yayılma yönüne 36 bağlıdır. Yanlışsız ve kararlı sonuçlar için örgü dağılma hatası (grid dispersion error) kabul edilebilir seviyeye indirgenmeli ki bu da hücre boyutunu azaltarak başarılabilir. Başka bir hücre boyutu düşüncesi, sorun geometrisinin önemli özelliklerinin doğru olarak, modellenmiş olması gerektiğidir. Bundan daha küçük bazı özel geometri özelliklerinde, ilginin tepkisini belirlemede etkiler olmadıkça normal olarak bu 1/10 λ veya daha küçük hücreler yaparak otomatik olarak karşılanabilir. Örnek olarak ince tel antenleridir; çünkü tel kalınlığını 1/10 λ’ dan 1/20 λ’ ya değiştirirsek anten empedansını etkileyecektir. Başka bir örnek olarak, dikdörtgensel hücreler ile pürüzsüz yüzey modellemenin merdiven (staircase) etkileridir. Belirgin hatalara yol açabilen pürüzsüz hedeflerden düşük düzeyli saçılmanın hesaplanmasıdır. İyi sonuçlar için bu ve buna benzer durumlarda, çok küçük hücreler gereklidir veya aynı boyuttaki Yee hücrelerinden daha iyi olan gerçek geometriye yaklaşık olan özel ağlar ya da (sub-cell) alt hücre modellenmesi gibi alternatif ölçümler gereklidir. Hücre boyutu belirlenince, nesne ile her bir boyuttaki dışsal sınır arasındaki kabul edilir boşluk uzay miktarı ve nesneyi modellerken gerekli, hücre sayısı bulunur. Buradan FDTD boşluğu için gerekli toplam hücre sayısı belirlenir(Kunz,1993). 3.4.5 Kararlılık İçin Zaman Adımı Boyutu Hücre boyutu belirlenince, ∆t maksimum zaman adımı boyutu Courant kararlılık kriterine uymak zorundadır. Courant kararlılık kriterini anlamak için, bir FDTD ağı boyunca propagasyon yapan bir düzlemsel dalgayı düşünelim. Bu dalga üzerinde herhangi bir nokta bir zaman adımı içinde, bir hücreden fazlasına geçiş yapmamalıdır; çünkü FDTD ‘nin bir zaman adımı boyunca dalga, sadece bir hücreden yakınındaki hücreye yayılım yapabilir. Bu zaman adımı kısıtlamasını belirlemek için, biz bir düzlemsel dalga doğrultusu seçeriz. Böylece düzlemsel dalga alan nokta konumlarında çok daha hızlıca propagasyon yapar. Bu doğrultu, FDTD ağının örgü düzlemlerine dik olacaktır. “d” boyutlu bir ağ için (d=1,2 veya 3) ve bütün hücre kenarları ∆u ‘ ya eşit iken, problem içindeki herhangi bir ortamda “ v” maksimum propagasyon hızını buluruz. Serbest uzayda, genellikle ışık hızında kararlılık için: v∆ t ∆u≤ (3 .84 ) d 37 3-Boyutlu dikdörtgensel ağ içi daha çok genellersek: −1 ⎛ ⎞ ∆t ⎜c 1 1 1≤ ⎜ 2 + ⎟2 + 2 ⎟ (3.85) ⎝ (∆x) (∆y) (∆z) ⎠ Yapılan deneyimler gösterdi ki; Eşitlik (3.84) ve (3.85)’ de verilen ∆t değerinin gerçek hesaplamalar kesin sonuçlar üretecektir ve pek çok durumlarda daha kesin sonuçlar ∆t nin daha küçük değerlerini kullanarak elde edilemeyecektir. Eşitlik tutunca, ayrıklaşmış dalga gerçek dalga propagasyonuna son derece yaklaşır ve ağ saçılma hatası ( grid dispersion error) minimuma iner. Ancak burada kural dışı istisna durumlar ortaya çıkabilir. Malzemenin iletkenliği, sıfırdan çok büyük olduğu zaman, zaman adımı eşitliği (3.85) ile ilgili indirgenmiş olması gereken bir durum ortaya çıkacaktır. İletken malzemeler için (σ>0) kararlı hesaplamalar Courant sınırından daha küçük zaman adımları gerekecektir. Bu genellikle problem oluşturmaz; çünkü pek çok hesaplamalarda, zaman adımı boyutu serbest uzaydaki ışık hızı tarafından ayarlanır. İletken malzemede hız, serbest uzayda daha küçük olacağından, hem serbest uzay hem de iletken malzemeler içeren FDTD hesaplamasında, zaman adımı öyle olacaktır ki, Courant sınırı her yeri tatmin edecektir. Ancak yüksek iletken malzemeler içindeki kısa dalga boyları, çevreleyen serbest uzay bölgelerinden daha küçük FDTD hücrelerine gerek duyabilir. Zaman adımının Courant sınırının altına indirilmesini gerektiren diğer bir durum doğrusal olmayan malzemelerde karşımıza çıkar(Kunz,1993). 3.5 Isı Transfer Denkleminin Fdtd Yöntemi İle Analizi Sıcaklık profilini, malzemenin karakteristiklerine bağlı olarak, modelize edebilmek için çeşitli malzeme ve yöntemler kullanılır. Bu modelizasyonda, işlemi yöneten denklem sistemleri kendi aralarında ilişkili ve birçok değişken parametreye bağlıdır. Böyle bir işlemi, incelemek çok karmaşık olduğundan, çoğu zaman geliştirilmiş matematiksel yöntemlere ihtiyaç vardır. Mikrodalga ısıtma işlemlerinde kullanılan ısı denkleminin ayrıklaştırılmasında Açık (Explicit), Kapalı (Implicit) ve Crank Nicolson yöntemleri kullanılabilir. Yapılan araştırmalarda; uzun sürecek mikrodalga ısıtma işlemlerinde Implicit ve Crank Nicolson yöntemi Explicit’e göre daha avantajlı olduğu 38 bilinmektedir. Implicit ve Crank Nicolson yöntemi karşılaştırıldığında, Crank Nicolson yönteminin uzun ısıtma sürelerine çıkıldıkça, Implicite göre daha doğru sonuç vereceği görülmüştür. Bunun sebebi ise, Crank nicolson yönteminin, Implicite göre uzun program koşturma süresi gerektirmesine rağmen, hata payı az olan sonuçların bulunması için daha avantajlı olduğu yapılan çalışmalarla görülmüştür. 3.5.1 Crank-Nicolson yöntemi ile ısı Denkleminin Ayrıklaştırılması Crank Nicolson tarafından 1947’de gerçekleştirilen yöntemde yanlışsız daha kesin sonuçlar elde edilmektedir. Bu yöntemde zaman türevi için(Özışık,1994): ∂T Tn+1 n = (i,k) −T(i,k) (3.86) ∂t ∆t konum türevi için ise: ∂2T 1 ⎡T n+1 n+1 n+1 n n n = (i−1,k ) − 2T(i,k ) +T(i+1,k ) T+ (i−1,k ) − 2T(i,k ) +T(i+1,k ) ⎤ ⎢ ⎥ (3.87) ∂x2 2 2 2⎢⎣ (∆x) (∆x) ⎥⎦ eşitlik (3.19)’deki ısı denklemine yukarıda elde edilen türev denklemlerinin ayrıklaştırılmış ifadesi uygulanırsa ısı denklemi (3.88)’deki gibi çıkar. T n+1 n(i ,k ) −T(i ,k ) k ⎡T n+1 − 2T n+1 n+1 n n n ρ .C = (i−1,k ) (i ,k ) +T(i+1,k ) T(i−1,k ) − 2T(i ,k ) +T+ (i+1,k ) ⎤ P ⎢ ⎥∆t 2 2 2⎢⎣ (∆x) (∆x) ⎥⎦ k ⎡T n+1 − 2T n+1 n+1 n(i−1,k ) (i ,k ) +T(i+1,k ) T(i−1,k ) − 2T n n + (i ,k ) +T(i+1,k ) ⎤ ⎢ ⎥ + P A2 (∆z)2 (∆z)2 (i ,k ) (3.88) ⎢⎣ ⎥⎦ 39 Şekil 3.8: Crank-Nicolson yöntemi için Sonlu Fark Molekülleri Şekil 3.8’de A düğüm noktasına etki eden sıcaklıklar görülmektedir. Burada 3 bilinmeyen T n+1i değerlerine karşılık, 3 tane bilinen T ni değeri vardır. Bilinmeyen T n+1i sıcaklık değerleri bilinen T ni değerlerinden bulunur. Kare örgü durumunda ∆x=∆z=δ, ∆t zaman adımı PA(i,k) poynting vektörü ve r Kararlılık kriteri: P A ∆t(i,k ) = P(i,k ) (3.89)ρCP ρC δ 2 ∆t ≤ P (3.90) k r k∆t= 2 (3.91)2ρCPδ Bu denklemleri kullanarak eşitlik(3.88)’in alacağı yeni durum: T n+1 n(i,k ) −T(i,k ) = r[T n+1(i−1,k ) +T n+1 n+1 n+1 n+1(i,k−1) −4T(i,k ) +T(i+1,k ) +T(i,k+1) ] (3.92) + r[T n n n n n A(i−1,k ) +T(i,k−1) −4T(i,k ) +T(i+1,k ) +T(i,k+1) ]+ P(i,k ) Bu denklemde (n+1). Değerler eşitliğin sol tarafına çekilirse: (1+ 4r )T n+1 − r[T n+1 + T n+1 + T n+1 n+1( i ,k ) ( i−1,k ) ( i ,k−1) ( i+1,k ) + T( i ,k+1) ]= (3.93) (1− 4r )T n( i ,k ) + r[T n(i−1,k ) + T n n n A(i ,k−1) + T( i+1,k ) + T(i ,k+1) ]+ P( i ,k ) 40 (3.93)’deki denklemin çözümünde Gauss-Seidel çözüm yönteminden yararlanılacaktır. Bu denklemde çözüm uzayımızın dışına çıkan noktaların yerine, sınır koşularından elde ettiğimiz değerler kullanılacaktır. 3.5.2 Isı Transfer Denkleminde Sınır Koşulu Isı denkleminde sınır koşullarının, sonlu farklar biçimindeki gösterimi aşağıdaki gibidir. Şekil 3.9’da malzemenin yüzeyindeki sınırlar ve ısının dağılım yönü gösterilmiştir. Şekil 3.9: Malzeme yüzeyindeki sınırlar ve ısının dağılım yönü nr r .J = H(T −T0) (3.94) r (J =−k∇.T) Yukarıdaki konveksiyon denklem referans alınarak sınırlardaki denklemler: x=0 için: k∂T =H(T−Text) (3.95.a)∂x x=d için: k∂T− =H(T−Text) (3.95.b)∂x 41 z=0 için; k∂T =H(T−Text) (3.95.c)∂z z=h için: −k∂T =H(T−Text) (3.95.d)∂z Burada, x=0, x=d, z=0 ve z=h değerleri sınır yüzeylerini ifade eder. Bu ifadelere, konum türevleri için merkezi farklar uygulanırsa: x=0 için: T n(1,k) −T n k (−1,k) =H(T n(0,k) −Text) (3.96.a)∆x T n Hδ n n(−1,k) =− (Tk (0,k) −Text)+T(1,k) x=d için: Tn −Tn −k (M+1,k) (M−1,k) =H(Tn −T ∆x (M,k) ext ) (3.96.b) Tn Hδ n n(M+1,k) =− (Tk (M,k) −Text)+T(M−1,k) z=0 için: T n −T n k (i,1) (i,−1) = H(T n(i,0) −Text) (3.96.c)∆z T n Hδ(i,−1) = − (T n (i,0) −T n k ext )+T(i,1) z=d için: T n n − k (i ,N +1) −T(i ,N−1) = H (T n(i ,N ) −T∆z ext ) (3.96.d ) T n Hδ= − (T n −T ) +T n(i ,N+1) k (i ,N ) ext (i ,N −1) 42 elde edilir. Burada β=Ηδ/k biçiminde ifade edilebilir. Sınır koşulları da yazıldıktan sonra, ısı denklemi malzemenin sınırlarında yazılsın. i=0 ve k=0 için: (1+ 4r)T n+1 n+1 n+1 n+1 n+1(0,0) − r[T(−1,0) + T(0,−1) +T(1,0) + T(0,1) ]= (3.97.a) (1− 4r)T n n n n n A(0,0) + r[T(−1,0) +T(0,−1) +T(1,0) + T(0,1) ]+ P(0,0) sınır değerler T(-1,0) ve T(0,-1) yerine konulursa: [1+ 2r(2+ β )]T n+1 − 2r(T n+1 +T n+1(0,0) (1,0) (0,1) ) = [1− 2r(2+ β )]T n n n(0,0) + 2r(T(1,0) +T(0,1) ) (3.97.b) +4βrT Aext + P(0,0) i=0 ve k=N için, sınır koşullarının da eklenmiş haliyle ısı denklemi: [1+2r(2+β)]T n+1 n+1 n+1(0,N) −2r(T(0,N−1) +T(1,N) ) = [1−2r(2+β)]T n n n(0,N) +2r(T(0,N−1) +T(1,N) ) (3.98) +4βrT Aext +P(0,N) i=M ve k=0 için aynı yol uygulandığında: [1+2r(2+β)]T n+1 −2r(T n+1 +T n+1(M ,0) (M−1,0) (M ,1) ) = [1−2r(2+β)]T n n n(M ,0) +2r(T(M−1,0) +T(M ,1) ) (3.99) +4βrText +P A (M ,0) i=M ve k=N için de aynı yol uygulandığında: [1+ 2r(2+ β )]T n+1 n+1 n+1(M ,N ) − 2r(T(M ,N−1) +T(M −1,N ) ) = [1− 2r(2+ β )]T n(M ,N ) + 2r(T n(M ,N−1) +T n(M −1,N ) ) (3.100) +4βrT Aext + P(M ,N ) ifadeleri bulunur. 43 (3.97)’den (3.100)’e kadar bulunan eşitlikler, malzemenin köşelerinde geçerli olan matematiksel ifadelerdir. Malzemenin yüzeylerinde oluşan ısı dağılımı için ise ısı denklemi sonlu farklara uygun biçimde benzetilsin. i=0 ve k≠0, k≠N için, sınır koşulları eklendiğinde: [1+ r(4+β)]T n+1(0,k ) − r(2T n+1 +T n+1 n+1(1,k) (0,k−1) +T(0,k+1) ) = [1− r(4+β)]T n + r(2T n(0,k) (1,k) +T n n+1(0,k−1) +T(0,k+1) ) (3.101) +2βrT Aext +P(0,k) denklemi elde edilir. i=M ve k≠0, k≠N için [1+ r(4 + β )]T n+1 − r(2T n+1 n+1 n+1(i ,0) (i ,1) +T(i−1,0) +T(i+1,0) ) = [1− r(4 + β )]T n(i ,0) + r(2T n n(i,1) +T(i−1,0) +T n+1(i+1,0) ) (3.102) +2βrText + P A (i ,0) bulunur. i≠0, i≠M ve k=0 için: [1+ r(4+ β)]T n+1(M ,k) − r(2T n+1 n+1(M−1,k) +T(M ,k−1) +T n+1(M ,k+1) ) = [1− r(4+ β)]T n n n n+1(M ,k) + r(2T(M−1,k) +T(M ,k−1) +T(M ,k+1) ) (3.103) +2βrText +P A (M ,k) i≠0, i≠M ve k=N için: [1+ r(4+β)]T n+1 − r(2T n+1 n+1 n+1(i,N ) (i,N−1) +T(i−1,N ) +T(i+1,N ) ) = [1−r(4+β)]T n + r(2T n +T n +T n+1(i,N ) (i,N−1) (i−1,N ) (i+1,N ) ) (3.104) +2βrText +P A (i,N ) en genel halde FDTD ile ayrıklaştırılmış ısı denklemi: [1+r(4+β)]Tn+1 −r(Tn+1 n+1 n+1 n+1(i,k) (i−1,k) +T(i,k−1) +T(i+1,k) +T(i,k+1) )= [1−r(4+β)]Tn +r(Tn +Tn n n(i,k) (i−1,k) (i,k−1) +Tii+1,k) +T(i,k+1) ) (3.105) +PA(i,k) bulunur. 44 3.6 Gauss-Seidel Yöntemi Doğrusal cebirsel denklemlerinin çözümü sayısal hesaplamalarda en çok rastlanan problemlerden biridir. Böyle bir problem: 1. İki noktalı sınır değer probleminin sonlu farklar yöntemi ile çözümünde, 2. Kısmi diferansiyel denklemlerinin sonlu farklar yöntemi ile çözümünde, 3. Doğrusal programlama ile yapılan optimizasyon hesaplarında, 4. Sonlu elemanlar yöntemi ile yapılan çözümlerde, olduğu gibi çeşitli uygulamalarda çıkabilir. Bir doğrusal cebirsel denklemin çözümünde kullanılabilecek bir yöntemin seçiminde bu yöntemin çözüm sırasında gerektirdiği işlem sayısı, özellikle çarpma ve bölme işlem sayısı, önemli rol oynar. Bu arada yöntemin kolay programlanabilir olması, mümkün olduğu kadar az yuvarlatma yanlışlarına yol açması ve sonuçların güvenilir olması da yöntemin seçiminde göz önünde tutulacak diğer noktalardır. Homojen denklemlerin çözümü için kullanılabilecek yöntemler dolaysız ve dolaylı olarak ikiye ayrılır. Dolaylı yöntemler hem algoritmanın kolayca programlanabilir olması hem de yuvarlama hatalarının az ve iterasyonlar yapıldıkça birikme olmamasından dolayı kullanılırlar. Hızlı olması bakımından Gauss-Seidel iterasyonu dolaylı yöntemler arasında en çok kullanılan yöntemdir(Özışık,1994). Nokta iteratif yöntemine dayanan bu yöntem, büyük ve seyrek (içinde çok 0 olan) sistemlerin çözümünde kullanılır. Çok basit bir yöntemdir. Yöntemde izlenilecek yol aşağıdaki gibidir. 1. Ana köşegenin bilinmeyen değerleri için her bir eşitlik çözülür. 2. Bütün bilinmeyenler için bir önceki hesaplamalar yapılır. 3. Birinci adımda ana köşegendeki bilinmeyenleri başarıyla çözebilmek için yaklaşımlar hesaplanır. Bu şekilde her bir hesaplamalar pozitifse en çok bu aralarda iterasyonun birinci yaklaşımı kullanılır. 4. İterasyonun birinci yaklaşımından değerler belirlenir ve ne zaman pozitifse hesaplanan değerler iterasyonun ikinci yaklaşımını tamamlamada kullanılır. 5. Açıkça belirtilen yakınsama kriterlerine kadar devam eden bu yöntem bütün bilinmeyenler için geçerlidir. 45 Gauss-Seidel yöntemi ile ilgili denklemler aşağıdaki biçimde verilsin. a 11 T 1 + a 12 T 2 + .... + a 1 n T n = d 1 ( 3 .106 .a ) a 21 T 1 + a 22 T 2 + .... + a 2 n T n = d 2 ( 3 .106 .b ) .......... ..... a n 1T 1 + a n 2 T 2 + .... + a nn T n = d n ( 3 .106 .c ) Burada i=1’ den n’ye kadar aii≠0’ dır. 3x3 booyutlu bir matris incelendiğinde, ana köşegendeki bilinmeyen değerler aşağıdaki biçimde yazılır. T 11 = ( d 1 − aa 12 T 2 − a 13 T 3 ) (3 .107 .a ) 11 T 12 = ( d 2 − a 21 T1 − a 23 T 3 ) (3 .107 .b )a 22 T 13 = ( d 3 − a 31 T1 − a 32 Ta 2 ) (3 .107 .c ) 33 İlk atadığımız değerler T (0) , T (0)1 2 ve T (0)3 seçilsin. Bu T değerlerini eşitlik(3.107)’ de yerine koyarsak: T (1) 11 = (d 1 − a ( 0 ) ( 0 ) a 12 T 2 − a13 T3 ) (3 .108 .a ) 11 T (1) 12 = (d − a T (1) 2 21 1 − a T ( 0 ) ) (3 .108 .b ) a 23 322 T (1) 1 (1) (1)3 = (d 3 − a 31T1 − a 32 T 2 ) (3 .108 .c )a 33 Birinci değerleri bulabiliriz. Yukarıdaki eşitliklerde görüldüğü gibi eşitlik(3.108.a)’ da bulunan T (1)1 değeri, eşitlik (3.108.b)’ de yerine konur ve ilk atadığımız T (0)3 değeri, korunur. İşlemlerin her biri aynı sırayı takip eder ve bu şekilde bütün köşegen değerlerini bulmamız kolaylaşır. Eşitlik (3.108)’ i genellersek: T ( n+1) 1= (d ( n ) ( n )1 1 − a12 T2 − a13T3 ) (3.109 .a )a11 T ( n+1) 1= (d − a T ( n+1)2 2 21 1 − a 23T ( n ) 3 ) (3.109 .b )a 22 T ( n+1) 13 = (d − a T ( n+1) − a T ( n+1)3 31 1 32 2 ) (3.109 .c )a 33 46 denklemleri elde edilir. Daha düzgün ve genel bir ifadeyle yukarıdaki eşitlikler aşağıdaki biçimde belirtilir: i−1 M T (n+1) 1 ⎧ ⎫ i = ⎨di −∑a (n+1)ijT j −∑aijT (n)j ⎬ (3.110)aii ⎩ j=1 j=i+1 ⎭ Bu işlemler (convergence) yakınsaklık kriteri sağlanana kadar devam eder. Yakınsaklık kriteri: T n +1i − T n i n +1 ≤ ε (3 .111 )Ti 47 3.7 Malzeme Parametreleri Bir malzemenin mikrodalga ısıtma sürecinde sıcaklık dağılımını elde etmek için bu malzemenin elektriksel ve termal özelliklerinin bilinmesi gerekir. Mikrodalga ısıtma ve kurutma işlemlerinde, fırın içindeki malzemeye ait dielektriksel özellikler ve nem içeriği sıcaklığa bağlı olarak değişir. Malzemenin sıcaklıkla değişen parametreleri; Cp (özgül ısı), k (kappa, ısı iletkenlik katsayısı) ve ε " (dielektrik kayıp faktörü) değerleridir. Diğer parametreler Hc (konveksiyon katsayısı) ve ρ (yoğunluk). Bu tezde incelenen malzemeler tuğla (Al2O3 -SiO2), motor yağı ve sudur. Tuğla için dielektrik sabiti 3-4 arasındadır; dielektrik kaybının sıcaklıkla değişimi, Şekil 3.10’da görülmektedir. 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 100 200 300 400 500 600 700 800 900 1000 sicaklik (C) Şekil 3.10: Tuğla için dielektrik kayıp sabitinin sıcaklıkla değişimi (Goodson, C,1997) Çizelge 3.1: Tuğlanın termo-fiziksel malzeme özellikleri (Metaxas,1996) Özgül ısı J/kgK 840 Isı iletkenlik 0,69 W/mK Yoğunluk 1600 Kg/m3 Nem oranı %10 er" 48 Kullanılmış motor yağının Dielektrik sabiti 2.1-3 ;dielektrik kayıp faktörü 0.2- 0.3 değerleri arasındadır. Suyun dielektrik sabiti 77 dielektrik kayıp faktörü 12 alınmıştır. Şekil 3.11’de su için sıcaklığa bağlı kayıp faktörü ( dielectric loss foctor) verilmiştir. (Arıcı ve Özkaynak,1978) Çizelge 3.2: Motor yağının termo-fiziksel malzeme özellikleri Sıcaklık 20 ºC 40 ºC 60 ºC 80 ºC 100 ºC 120 ºC Özgül ısı 1851 1934 2018 2102 2186 2269 J/kgK Isı iletkenlik 0,144 0,143 0,142 0,141 0,140 0,138 W/mK Yoğunluk 871 858 845 832 820 807 Kg/m3 Çizelge 3.3: Suyun termo-fiziksel malzeme özellikleri Sıcaklık 20 ºC 40 ºC 60 ºC 80 ºC 100 ºC 120 ºC Özgül ısı 4183 4178 4191 4216 4216 4232 J/kgK Isı iletkenlik 0,598 0,627 0,651 0,669 0,682 0,685 W/mK Yoğunluk 998 992 983 972 Kg/m3 Şekil 3.11: Suyun Sıcaklığa Bağlı Dielektrik Kayıp Faktörü ε" 49 3.8 DENEYSEL ve FDTD SONUÇLARI Mikrodalga ısıtmada, zaman arttıkça malzemenin sıcaklığı da doğru orantılı olarak artacaktır. Sıcaklık dağılımlarının bilgisayar ortamında FDTD ile çözülen denklemlerin sonuçlarında elde edilen sayısal sonuçlar kesin olmasa da hata payı az olan, yaklaşık sonuçlar elde etmemizi sağlayacaktır. Sıcaklık malzemenin her yerinde aynı olmayabilir. Çünkü malzemeye bir elektrik alan uygulandığı zaman malzeme bu alanın etkisiyle ısınmaya başlayacaktır. Malzemenin sıcaklığı, elektrik alanın uygulandığı noktada daha fazla olur. Malzemenin kenarlarına doğru gidildikçe ısıda düşme olur. FDTD yönteminin uygulamasında, dalga kılavuzu x ekseninde 20 yee hücresine, y ekseninde 10, z ekseninde 40 yee hücreye bölünmüştür. Kaynak z ekseninde 10.uncu hücreye yerleştirilmiştir. Konum adımı 4.3(mm) alınmıştır. İncelenen malzeme 20(mm)x20(mm) ebatlarında uygulanan frekans 2.45GHz konum adımı 4.3(mm) alınmış hücre sayısı 5 olmuştur. 86(mm) x43(mm) boyutlarındaki sonu kısa devre edilmiş bir dalga kılavuzu aşağıda Şekil 3.12’de görülmektedir. Şekil 3.12: Deney Düzeneği 50 Deney için kullanılan aletler 2.45GHz frekansına sahip güç kaynağı (power supply), magnetron, magnetronu yansıyan güçten koruyan sirkülator, sirkülatorun yönlendirdiği su yükü (water load), gelen ve yansıyan güçleri ölçmek için bir yönlü kuplör, uygunlaştırıcı 3-Stub Tuner, WR340 dalga kılavuzu ve kısa devre pistonu. Kılavuzun kesit boyutları 86(mm) x 43(mm) ve uzunluğu 172(mm)’dir. z=0 düzlemine yutucu yüzey, z=172(mm) düzlemine kısa devre yerleştirilmiştir. Uyarma z=43(mm) düzleminden zaman faktörü f=2.45GHz frekanslı bir sinüsoidal bir işaret olan TE10 modunun elektrik alanı ile yapılmıştır. Kaynak terimi: Ey (x,t) E π = 0 sin xG(t) (3.112)a Formunda verilen TE10 modunun elektrik alanının y bileşenidir. Burada E0 alanın genliği ve a dalga kılavuzunun geniş olan boyutudur. G(t) fonksiyonu ise alanın zamana bağımlı olduğunu gösterir. Zaman faktörü bir gauss darbesi seçilebileceği gibi tek ton sinüsoidal bir işaret de seçilebilir. Mikrodalga jeneratörlerinde çıkış genellikle güç boyutunda verildiğinden bu güce bağlı olarak alanın genliği: E 2 4P µ 10 = 0 (3.113)ab ε λ 20 1−⎛ ⎞⎜ 0 ⎟ ⎝ 2a ⎠ şeklinde belirtilir. Burada P jenaratörün çıkış gücü, a ve b dalga kılavuzunun boyutları ve λ0 serbest uzayın dalga boyudur.(Oktay ve Akman,2003) Deneysel uygulamada dielektrik malzeme olarak ilk önce su analiz edilmiştir. Deneysel düzenek kurulduktan sonra kılavuz içine magnetron vasıtasıyla güç verilmiştir. Dielektrik üzerindeki sıcaklık bir ölçülmüş, zamanla sıcaklığın arttığı görülmüştür. Değişen güç miktarlarına göre sıcaklık-zaman eğrileri elde edilmiş ve FDTD yöntemiyle bu sıcaklık zaman eğrisi program vasıtasıyla bulunmuştur. Kılavuz içine verilen 100Watt için oluşturulan sıcaklık zaman eğrileri Şekil 3.13’ de verilmiştir. 51 100 90 80 70 60 50 40 ** deneysel 30 -- FDTD 20 10 0 20 40 60 80 100 120 zaman(sn) Şekil 3.13: P=100Watt için suyun deneysel ve FDTD sonuçları FDTD yöntemi uygulanarak alınan simulasyon sonuçlarından, başlangıçta hızlı bir artış görülmektedir, sonraki sıcaklıklarda ise bu artış daha yavaş olmaktadır. Deneyden alınan sonuçlarda ise başlangıçtan itibaren doğrusal bir artış görülmektedir. Deneyde sudan farklı dielektrik özelliklere sahip bir malzeme olarak, kullanılmış motor yağı analiz edilmiştir. Aşağıda simulasyon sonuçları ile deneyden elde edilen sonuçlar birlikte verilmiştir. Şekil 3.14’de farklı güçler için deneysel verilerle elde edilen motor yağının mikrodalga güç altında zamana bağlı sıcaklık grafiği görülmektedir. 140 P=180W 120 100 80 60 P=72W 40 20 0 200 400 600 800 1000 1200 1400 zaman(sn) Şekil 3.14: P=180w ve P=72w için motor yağının deney sonuçları SICAKLIK (C) SICAKLIK (C) 52 Deneyde motor yağı analiz edilirken, önce düşük güçte mikrodalga verilmiş ve sonuç gözlemlenmiştir. Düşük güçlerde motor yağının çok geç ısındığı görülmüştür. Daha yüksek güç verildiğinde sıcaklığının daha hızlı bir şekilde arttığı gözlemlenmiştir. Daha düşük güçlerde deney düzeneğindeki cihazlardan kaynaklanan kayıplar çok fazladır. 140 120 100 80 60 ** Deneysel -- FDTD 40 20 0 50 100 150 200 250 300 350 400 zaman(sn) Şekil 3.15: Kullanılmış motor yağının deneysel ve FDTD sonuçları 60 55 50 45 40 35 ** Deneysel -- FDTD 30 25 20 0 200 400 600 800 1000 1200 1400 zaman(sn) Şekil 3.16: P=72W için motor yağının deneysel ve FDTD sonucu SICAKLIK (C) SICAKLIK (C) 53 Şekil 3.15 ve Şekil 3.16’da hem deneysel hem de simulasyon sonucu birlikte görülmektedir. Düşük güçlerde magnetronun kararsız olması, simulasyon sonucu ile deney sonuçları arasında fark oluşturmaktadır. Ayrıca deney düzeneğinde oluşan kayıpları da göz önüne almak gerekir. 120 110 100 90 80 70 60 * Deneysel - FDTD 50 40 30 20 0 10 20 30 40 50 60 70 80 90 t(sn) Şekil 3.17: P=360W için tuğla (Al2O3SiO2) deneysel ve FDTD sonucu Laboratuarda analiz edilen diğer malzeme ise, su ve yağdan daha düşük bir kayıp faktörüne sahip olan Işıklar Holding AŞ’den temin edilen tuğla olmuştur. Tuğlanın bu sahip olduğu dielektrik özellikten dolayı gücü emme (absorbe etme) ve ısıya dönüştürme yeteneği daha az olduğu için yüksek güç verilmiş ve gözlemler sonucu içerdiği nem dolayısıyla sıcaklık birden artmış kurumaya geçtikten sonra ise sıcaklık bir dengeye ulaşmış ve sabitlenmiştir sıcaklık artması azalmıştır. Dielektrik kayıp faktörü sıcaklıkla artıyorsa ısıtma hızla artacaktır. Dielektrik kayıp faktörü düşük olan malzemelere verilen güç miktarı fazla olmalıdır. FDTD ile simulasyonu sonucunda da benzer sonuç elde edilmiştir. SICAKLIK(C) 54 4. ARAŞTIRMA SONUÇLARI VE TARTIŞMA Bir elektromanyetik alan altında bulunan bir dielektrik malzemenin davranışı incelendiğinde, malzemenin yapısına bağlı olarak elektromanyetik enerji ya bu malzeme üzerinde toplanır ya da dielektrik kayıplarından dolayı ısı enerjisine dönüşür. Malzeme ile elektromanyetik alan arasındaki etkileşimde oluşan termik olayları anlamak için, malzemenin hem elektromanyetik alan altındaki davranışı, hem de ısınmasından dolayı ortaya çıkan olayları incelemek gerekir. Mikrodalga frekanslarındaki elektromanyetik enerjinin ortamla etkileşmesinde, elektrik ve termik mekanizmalar ortaya çıkar. Bu mekanizmaları, matematiksel model altında belirli bir yaklaşımla incelemek için, hem elektromanyetik olayları hem de ısı olaylarını ifade eden matematiksel modelleri oluşturmak ve bunlar arasındaki ilişkileri ortaya koymak gerekir. Yüksek frekanslı radyo ya da en genel halde elektromanyetik dalgaların, dielektrik kayıplarla ısıya dönüşme olayını kullanan endüstriyel olaylarda, işlemin davranışını ve kontrolünü ancak bu modellerle benzetebiliriz. Dielektrik kayıplardan yararlanılarak malzemelerin elektromanyetik enerji ile ısıtılması süreçlerinde karmaşık ve zincirleme olaylar dizisi oluşur. Bu olaylar ortamda meydana gelen etkileşmeler neticesinde ortaya çıkar. Bu etkileşmeler mikrodalga ısıtma sürecinde iki şekilde ortaya çıkar. Birincisi; elektromanyetik alanlar, sıcaklık ve ortamın fiziko-kimya reaksiyonlarıdır. Malzemenin ısınması, ortamda kütle transferi ve bazı kimyasal reaksiyonların oluşmasına sebep olacaktır. Bu da malzemenin elektriksel ve fiziksel özelliklerinin değişmesine sebep olacağından ortamdaki elektromanyetik alan dağılımının da değişmesine neden olacaktır. Elektromanyetik alan dağılımındaki değişim malzeme üzerindeki güç dağılımını da etkileyeceğinden ısınma süresinin de değişmesine sebep olacaktır. Bu karşılıklı etkileşme süreç boyunca devam edecektir. İkincisi; malzeme, fırın ve ortama gücü sağlayan kaynaktan oluşan elektromanyetik alan ile malzeme arasındaki ilişkidir. Malzemenin ısıtılması neticesinde dielektrik özelliklerinde meydana gelecek olan değişme ortamdaki elektromanyetik alan dağılımını değiştirmesine ve dolayısıyla fırının yapısına bağlı olarak ortamın mikrodalga kaynağına verdiği tepkinin değişmesine sebep olacaktır. Bu değişiklik de yine ısıtma üzerinde etkili olacaktır. Böylece zincirleme bir çevrim başlayacaktır. 55 Kayıplı bir dielektrik malzeme ile yüklenmiş bir mikrodalga uygulayıcısı içindeki elektromanyetik alan dağılımının ve yük üzerindeki sıcaklık dağılımının bulunması amacıyla Maxwell denklemlerinin ve ısı transfer denklemlerinin FDTD- Sonlu Farklar Zaman Domeni Yöntemi ile çözümlerine değinilmiştir. Bu amaçla, elektromanyetik alanların ve ısı transfer denkleminin FDTD ile modellenmesi yapılmıştır. FDTD yöntemi, Maxwell Denklemlerindeki kısmi türev operatörlerinin merkezi farklara dayalı sonlu farklar karşılıkları ile değiştirilip, doğrudan zaman ve konum domenlerinde sayısallaştırılmasına dayanır. Bu yöntem zaman adımlı yaklaşıma dayanmaktadır. Elektromanyetik Alanların modellenmesinde; zamana bağlı Maxwell denklemlerinde konum ve zaman türevleri yerine Taylor serisi açınımlarından elde edilmiş fark denklemleri kullanılarak, sürekli konum ve zamanda tanımlı Maxwell denklemlerinin ayrıklaştırılması yapılmıştır. Bu ayrıklaştırmada, konum türevleri yerine merkezi, zaman türevleri yerine de ileri fark denklemleri kullanılmıştır. Ayrıca, gerekli sınır koşullarının matematiksel ifadeleri elde edilerek, bu ifadelerde FDTD yöntemi uygulanmıştır. Mikrodalga ısıtma işlemi için, ısı transfer denklemi FDTD ile çözümlenmiştir. FDTD ile ısıtma sistemlerini analiz etmek için pek çok yaklaşım vardır. Bunlardan biri Explicit (açık) Yöntemidir. Bu yöntem, çok basit hesaplamalar içindir. Zaman adımının büyük olması durumunda kararlılık kriteri tarafından sınırlandırılır. Hesaplamalar ise, büyük zaman dilimlerinde gerçekleşir. Bu yüzden adımların sayısının gerekenden büyük olması hesaplamaların sayısına engel teşkil eder. İmplicit (kapalı) ve Crank- Nicolson diğer yöntemlerdendir. Bu yöntemler Explicit’e göre daha uzun ısıtma süreçlerinde daha avantajlıdır. Isı Transfer Denklemleri’ni çözümlemekte Crank- Nicolson yöntemi, hata payı daha az ve doğru sonuçlar vereceğinden bu yöntem esas alınmıştır. Crank-Nicolson Isı Transfer Denklemi’ni sayısal olarak çözerken, zamanda ve konumda ikinci dereceden yaklaşım sağlamaktadır ve çözüm koşulsuz kararlıdır. Fark denklemlerini elde etmek için, ikinci dereceden konum türevleri ve zaman türevleri ileri sonlu fark ifadeleri kullanılır. Fark denklemlerini çözerken matris tersi gerektiğinden, bunu sayısal olarak çözümlemek çok güçtür. Bu sebepten bu denklemleri iteratif yöntemlerden olan Gauss-Seidel yöntemi ile çözümlenmiştir. 56 Dalga kılavuzu içindeki malzemeye elektromanyetik alan uygulanmıştır. Bu alanın malzeme üzerindeki etkisi incelenmiş ve elde edilen sonuçlar değerlendirilmiştir. Sıcaklık dağılımının simulasyonunu yapmak için dalga kılavuzunda tanımlanan FDTD yöntemi kullanılmıştır. Çalışmada kullanılan malzeme dalga kılavuzu içinde, elektriksel alanın maksimum olduğu konumda yerleştirilmiştir. Böyle bir sistemde, malzemenin mikrodalga ısıtma hesaplamalarını yapabilmek için ısı transfer denkleminin çözülmesi gerekmektedir. Isı transfer denklemini çözerken, elektromanyetik alan dağılımının hesaplanması zorunludur. Magnetrondan verilen gücün ve yayılan modun belli olmasından, boş bir dalga kılavuzunun merkezindeki elektrik alan yoğunluğunu hesaplamak mümkündür. TE10 dalga kılavuzunda elektrik alan dağılımı y-koordinatında değişiklik göstermez ve kılavuz y-koordinatında düzgün dağılımlı olarak dağıldığından, ısı transfer denklemi iki boyutta çözülmüştür. Modellemede Elektrik alan ve Manyetik alanları hesap edilir ve buradan güç yoğunluğu hesaplanır. Elektromanyetik alanlar ve güç yoğunluğu, kararlı hale ulaşana kadar tekrar tekrar hesaplanarak güncellenir. Malzeme üzerindeki her bir hücredeki sıcaklık artışı 1sn periyotla hesaplanmıştır. Bu sıcaklıklar kullanılarak, dielektrik sabitinin ve dielektrik kayıp faktörünün yeni değerleri hesaplanır ve bu değerler elektrik ve manyetik alanlarının ve güç yoğunluğunun tekrar hesaplanmasında kullanılır. Hesaplamalar her ∆t adımında yapılır ve çevrim istenilen ısıtma zamanına ulaşılana kadar tekrarlanır. Isıtma süreci, malzemenin fiziksel ve elektriksel özellikleriyle doğrudan ilişkilidir. Bu nedenle dielektrik sabiti, dielektrik kayıp faktörü, ısı iletkenliği, özgül ısı gibi karakteriksel parametreler, sıcaklıkla önemli ölçüde değişmektedir. Dielektrik sabiti malzemenin mikrodalga enerjisinin iç kısımdan geçmesi sırasında yutma (absorbe etme) yeteneğinin bir ölçüsüdür. Kayıp faktörü ise giren mikrodalga enerjisinin malzeme içerisinde ısı olarak tüketilmesiyle kayıp olma miktarını vermektedir. Bu nedenle yüksek kayıp faktörlü bir malzeme mikrodalga enerjiyle kolaylıkla ısıtılabilmektedir. Farklı malzemeler için mikrodalga emilimi, dolayısıyla mikrodalganın etkileri farklıdır. Yüksek dielektrik kayıplara sahip malzemeler su ve alkol gibi, mikrodalga radyasyonla kolay bir şekilde ısınırlar. Seramiklerin çoğu mikrodalga enerjiyi geçirimli veya zayıf yutan (absorbe eden) malzemelerdir. 57 Isıtma sırasında dielektrik katsayısı ve dielektrik kayıp faktörü sıcaklıkla değişir ve bunların değişim bilgisi işlem kontrolü için önemlidir. Kısacası, seramik malzemesinin mikrodalga ile ısıtılması frekans ve sıcaklığa bağlı olan kayıp faktörü ve dielektrik sabiti ile ilişkilidir. Enerji malzemeler tarafından emilir ve sıcaklık olarak atılır. Elektrik alan enerjisi malzemeye geçtiğinde enerjide azalma meydana gelir. Negatif ve pozitif kutuplara sahip polar moleküllü yapılar alternatif yüksek elektrik alanı içinde şiddetle titreşirler. Böylece moleküller hareketin ve çekimin direncinin üstesinden gelirler. Sürtünmenin yol açtığı ısı nedeniyle de malzeme sıcaklığı yükselir. Dielektrik sabiti malzemeye bir güç verildiğinde, malzemenin depolayacağı elektrik alan enerji miktarının tespit edilmesini sağlar. Elektrik alan kaynağa, kılavuzun boyutlarına ve malzemenin dielektrik özeliklerine bağlıdır. Malzemelerin dielektrik özellikleri de sıcaklığa bağımlı olduğu için, kılavuzdaki elektrik alan ve malzemenin sıcaklığı birbirine bağımlıdır sonucu çıkmaktadır. Bu karmaşık bileşik ilişki, mikrodalga kılavuz içinde ısıtılan malzemelerin sıcaklık hesaplamalarında dikkate alınmalıdır. Dalga kılavuzuna yerleştirilen dielektrik malzeme için elektrik alanı ve sıcaklık dağılımı, verilen mikrodalga güç boyunca FDTD yöntemiyle sayısal olarak analiz edilip incelenmiştir. Bu incelemelerden dalga kılavuzundaki elektrik alan içinde bulunan dielektrik malzemenin dielektriksel parametrelerine bağlı olduğu görülmüştür. Dielektrik malzemenin ısıtma süreci dielektrik malzemenin dielektriksel ve termal parametrelerine bağlıdır. Malzemenin dielektriksel (εr, tanδ) özelliklerinin değişimini hesaba katmak mikrodalga ısıtma açısından önemlidir. Dielektrik sabiti ve dielektrik kayıp faktörü sıcaklığa ve neme bağlıdır. Deneysel çalışmada malzeme olarak, tuğla, kullanılmış motor yağı ve çeşme suyu test edilmiştir. Deneydeki temel düzenekte bir dalga kılavuzu ile üzerinde sıcaklığı ölçmek bir açıklık, 2.45GHz frekansına sahip mikrodalga güç kaynağı, ayarlayıcı aygıt (tuner), mikrodalga gelen ve yansıyan güçleri ölçen güç metre bulunmaktadır. Ayarlanabilir kısa devre bir pistona sahiptir. Bu piston yardımıyla kılavuzun boyu ayarlanır. Kılavuzun boyunun ayarlanması, kılavuzun rezonansa gelmesini sağlamak ve malzeme üzerine maksimum elektrik alanın düşmesini sağlamaktadır. Malzeme elektrik alanın maksimum olduğu yerde konumlandırılmıştır. Kılavuzun boyutları, mikrodalga güç kaynağının frekansı verildiğinde kılavuzdaki elektrik alan bulunur. 58 Dielektrik malzeme tarafından yutulan güç direkt olarak ölçülemez. Yönlü kuplör, kılavuzun duvarları ve dielektrik tarafından emilen gücü dolaylı olarak ölçer. Ölçülen güç ikisi tarafından yutulan güçtür. Malzeme tarafından güç emildikçe sıcaklık da artmaktadır. Ölçüm malzemenin ortasından alınmaktadır. Burada dielektrik kayıp faktörü daha yüksektir ve kenarlara doğru gittikçe yani soğuk olan yerlere gittikçe azalır. Deneysel çalışmada önce su analiz edilmiş ve ısıtma sürecinde sıcaklık artışı gözlenmiştir. Suyun dielektrik sabiti ve dielektrik kayıp faktörü diğer dielektriklerden yüksektir. Bu nedenle su diğer malzemelerden daha kuvvetli ve mükemmel olarak mikrodalgaları emer ve hızla buharlaşır. Nem iyi yayılmasa bile mikrodalga belli bir uzaklıktan malzemeye nüfuz eder. Çok nem çok mikrodalga absorbsiyonudur. Yağ ve su çok farklı dielektrik özelliklere sahiptir. Bir yağ için tipik değer 2,1-3 iken su için bu değer 77’dir. Bir başka deyişle, suyun depoladığı enerji yağdan birim hacim başına daha fazladır. Sıcaklık artışı da buna bağlı olarak yağdan daha hızlıdır. Analiz edilen diğer malzeme olan tuğlanın mikrodalga enerjiyi emme kapasitesi diğerlerinden daha azdır. Dielektrik kayıp faktörü düşük olan malzemelere verilen güç miktarı fazla olmalıdır. Tuğlanın mikrodalga enerji ile ısınması sonucunda sıcaklığında birden artış gözlemlenmiştir. Sıcaklıktaki bu kontrol edilemez yükselmenin sebebi dielektrik kayıp faktörünün pozitif eğimidir. Hızlı sıcaklık gradyanı sonucunda malzemenin kimyasal özelliklerinin değişmesine neden olacaktır. Hızlıca dilektrik kaybının düşmesiyle yani mikrodalga enerjiyi emme yeteneğinin düşmesi sonucu sıcaklık artmayacak dengeye gelecektir. Burada aynı şartlar altunda en geç ısınan malzeme tuğla olmalıdır. Çünkü dielektrik kayıp faktörü diğer malzemelerden daha küçüktür. Enerjiyi absorbe etme yeteneği daha azdır. Değişik güç seviyelerinde farklı malzemelere uygulanan mikrodalga güç sonucunda sıcaklık zaman eğrileri elde edilmiştir. Simulasyon sonuçları ile deneysel sonuçlar arasında farklar görülmüştür. Bu farkların pek çok nedenleri olabilir. Bu nedenler arasında deneysel düzenekten kaynaklanan kayıplar gösterilebilir. Mikrodalga güç kaynağı verilen düşük güç seviyesinde kararlı değildir. Diğer bir neden de malzemelerin ısınmayla birlikte dielektrik özelliklerinin ve termal özelliklerinin 59 sıcaklıkla değişen parametreler olması yüzünden malzemenin ısınmasıyla ortaya çıkan kararsızlığı gösterilebilir. Bu tez çalışmasında farklı dielektrik özeliklere sahip 3 malzeme incelenmiştir. Dielektrik özelliklerinin elektrik alana ve sıcaklığa etkisi hem sayısal olarak hem de deneysel ortamda analiz edilerek sonuçları ortaya konulmuştur. Mikrodalga enerji ile işleme, daha hızlı süreç, daha düzgün dağılımlı ısıtma, daha kesin kontrol imkanı ve endüstride uygulamada daha tutarlı ürün kalitesi sağlamaktadır. Diğer taraftan Mikrodalga Isıtmanın yaygın olarak kullanılmasına engel olan en büyük etmen, ısıtılan malzemenin dielektrik özelliklerinin sıcaklığa bağlı olarak değişmesidir. Malzemenin dielektrik kayıp faktörü, malzemenin mikrodalga enerjiyi yutma (absorbe etme) ölçüsüdür. Eğer dielektrik kayıp faktörü sıcaklıkla artıyorsa ısıtma hızla artacaktır. Birçok endüstriyel malzemeler seramik ve polimerler gibi, tam bir ısıtma için kusursuz işleme sıcaklıkları gerektirir. Fazla ısıtma seramiklerde çatlamalara, polimerlerde erimeye neden olacaktır. Yüksek sıcaklıklarda bu malzemelerin kararsızlığı endüstriyel mikrodalga ısıtmada tüm erim için büyük bir sınırlamadır. 60 5. KAYNAKLAR ADU, B.and OTTEN, L.1996. Modelling Microwave Heating Characteristics of Granular Hygroscopic Solids. Journal of Microwave Power and Electromagnetic Energy, vol 31, no.1, p.35-42 AL-RIZZO, H.M., TRANQUILLA, J.M.AND MA,F. 2000. Incorporation of the Waveguide Feed and Cavity Wall Losses in to a Cartesian/Cylindrical Hybrid Finite- Difference Time Domain Analysis of a Geometrically- Composite Microwave Applicator. Journal of Microwave Power and Electromagnetic Energy, vol 35, no.2, p.110-118 ARICI, ÖNER ve ÖZKAYNAK, TANER. 1978. Katı Sıvı ve Gazların Termodinamik Özellikleri, Termas Yayınları. CHENG, DAVİD K., 1989. Field and Wave Electromagnetics. Addison-Wesley Publishing Company. CRAVEN, M., CROSS, T., and BİNNER, J.,1996 Enhanced Computer Modeling for High Temperature Microwave Processing of Ceramic Materials, in: Microwave Processing of 134 Materials V, eds. M. Iskander, J. J.O. Kiggans, and J. Bolomey, Vol. 430, pp. 351-356,MRS, Pittsburgh. COLLIN, R.E. Foundations for Microwave Engineering. The IEEE Press series on Electromagnetic wave Theory. CURTİS, J., 1999.Experimental Verification for Microwave Processing of Materials in a Single Mode Rectangular Resonant Cavity,M.S. Thesis, Virginia Polytechnic Institute and State University. 61 DE POURCQ, M. 1985. Field and power-density calculations in closed microwave systems by three-dimensional finite differences.IEE Proceedings 132, Pt. H 6:p.360-368 DIBBEN, D.C. and METAXAS, A.C. 1994. Finite Element Time Domain Analysis of multimode applicators using edge elements. Journal of Microwave Power and Electromagnetic Energy, vol 29, p.242-251 DUCHEZ, W. , 1996. Role of Electric Field Profiles in Continuous Microwave Processing of Thermal Runaway Materials, M.S. Thesis, Virginia Polytechnic Institute and State University. FRANDOS, S.2003. Mathematical model for indirect excitation of a wave-guide with electromagnetic waves in TE10 mode,9th International Conference on Microwave and RF Heating,Loughborough, U.K. FRANDOS, S.2003.FDTD numerical 3D application for microwave heating with a waveguide and a multimode cavity, 9th International Conference on Microwave and RF Heating,Loughborough, U.K. HAALA and WIESBECK, W.2000. Simulation of Microwave, Conventional and Hybrid Ovens Using a New Thermal Modelling Technique. Journal of Microwave Power and Electromagnetic Energy, vol 35, no.1, p.34-42 FU, W.B. and METAXAS, A.C.1994. Numerical prediction of three-dimensional power density distributions in a multi-mode cavity. Journal of Microwave Power and Electromagnetic Energy, vol 29, p.67-75 GOODSON, C.1997, Simulation of Microwave Heating of Mullite Rods, M.S. Thesis, Virginia Polytechnic Institute and State University. 62 HARMS, P.H., CHEN, Y., MITTRA, R., SHIMONY, Y. 1996 Numerical Modelling of Microwave Heating System. Journal of Microwave Power and Electromagnetic Energy, vol 31, no.2, p.114-121 ISKANDER, M., 1991 Computer Modeling and Numerical Techniques for Quantifying Microwave Interactions With Materials," in: Microwave Processing of Materials II, eds. W. Snyder, W. Hutton, M. Iskander, and D. Johnson, Vol. 189, pp. 149-171, MRS, Pittsburgh. JURGENS, T.G., TAFLOVE, A., UMASHANKAR, K.and MOORE, T.G.1992. Finite –difference time-domain of modelling of curved surfaces. IEEE Trans. Antennas Propagat. 40, p.357-366 KUNZ, K.S. and LUEBBERS, J. 1993. The finite difference time domain method for electromagnetics. CRC Press, Boca Raton, FL. LIU, F. , TURNER, I. and BIALKOWSKI, M.E. 1994. A Finite Difference Time Domain Simulation of Power Density Distribution in a dielectric Loaded Microwave Cavity. Journal of Microwave Power and Electromagnetic Energy, vol 29, no.3, p.138- 148 LIU, F., TURNER, I., SIORES, E. and GROOMBRIDGE, P. 1996. A Numerical and Experimental Investigation of the Microwave Heating of Polymer Materials Inside a Ridge Waveguide. Journal of Microwave Power and Electromagnetic Energy, vol 31, no.2, p.71-82 LYE, S.W.,SIORES,A.,TAUBE,A.and MORRISON,R.2000. Numerical Modeling of the Temperature Distribution Within a Bonded Paper Web When undergoing Microwave heating in a Waveguide. Journal of Microwave Power and Electromagnetic Energy, vol 35, no.2, p.125-131 63 MCCONNELL, BRIAN G. 1999. M.S Thessis. A Coupled Heat Transfer And Electromagnetıc Model For Sımulatıng Mıcrowave Heatıng Of Thın Dıelectrıc Materıals In A Resonant Cavıty, Virginia Polytechnic Institute and State University. METAXAS, A.C. 1996. Foundations of Electroheat. John Wiley&Sons Press, England MARTIN, D.MATEESCU,D.and JIANU,A.2000. The Influence of Microwave Heating on the Characteristics of Polyelecrolytes. Journal of Microwave Power and Electromagnetic Energy, vol 35, no.4, p.216-224 OKTAY, A and AKMAN,A.1999 Thermal simulation of microwave heating for material processing.7th Ampere HF&MW Heating Conference, Sept,1999,Valencia,Spain. OKTAY, A and AKMAN,A.2000.FDTD Modelling of dielectric heating Process: Application on ceramic materials. 7th Int. Conference on Optimization of Electrical and Electronic Equipment,pp.95-99,Brasov. OKTAY, A and AKMAN,A.2000.Numerical and experimental analysis of dielectric heating process. Journal of Electrical Engineering OKTAY, A and AKMAN,A.2003. An Analysis of the FDTD Method for Modeling the Microwave Heating of Dielectric Materials within 3D Cavity System. Journal of Electrical Engineering, Vol. 3. OKTAY, A and AKMAN,A.2001. Microwave Heat and Mass Transfer characteristics of dielectric materials. International Seminar HIS-2001, Padova,Italy. ÖZIŞIK, M.N. 1994. Finite difference methods in heat transfer. CRC Press, Boca, Raton, FL. 64 ÖZIŞIK, M.N. 1985.Heat Transfer a Basic Approach. Mc Graw-Hill Book Company,New York POZAR, D.M. 2005 Microwave Engineering. 3rd ed. – Hoboken, John Wiley & Sons . TRANQUILLA, J.M., MA,F. ,0AL-RIZZO, H.M.and CLARK, K.G. 1999. A Cartesian-Cylındrical Hybrid FDTD Analysis of Composite Microwave Applicator Structures. Journal of Microwave Power and Electromagnetic Energy, vol 34, no.2, p.97-105 ZHAO, H. and TURNER, I.W. 1996. An Analysis of the Finite Difference Time Domain Method for Modelling the Microwave Heating of Dielectric Materials Within a Three-Dimensional Cavity System. Journal of Microwave Power and Electromagnetic Energy, vol 31, no.4, p.199-214