OFİS ORTAMINDA TÜRBÜLANSLI AKIŞ MODELLERİNİN KARŞILAŞTIRILMASI VE ÇEŞİTLİ ISITMA SİSTEMLERİNİN ISIL KONFOR İNDİSLERİ VE GAGGE MODELİ KULLANILARAK SAYISAL SİMÜLASYONU Bahadır Erman YÜ CE T.C. BURSA ULUDAĞ ÜNİVERSİTESİ FEN BİLİMLERİ ENSTİTÜSÜ OFİS ORTAMINDA TÜRBÜLANSLI AKIŞ MODELLERİNİN KARŞILAŞTIRILMASI VE ÇEŞİTLİ ISITMA SİSTEMLERİNİN ISIL KONFOR İNDİSLERİ VE GAGGE MODELİ KULLANILARAK SAYISAL SİMÜLASYONU Bahadır Erman YÜCE Doç. Dr. Erhan PULAT (Danışman) DOKTORA TEZİ MAKİNE MÜHENDİSLİĞİ ANABİLİM DALI BURSA – 2019 U.Ü. Fen Bilimleri Enstitüsü, tez yazım kurallarına uygun olarak hazırladığım bu tez çalışmasında;  tez içindeki bütün bilgi ve belgeleri akademik kurallar çerçevesinde elde ettiğimi,  görsel, işitsel ve yazılı tüm bilgi ve sonuçları bilimsel ahlak kurallarına uygun olarak sunduğumu,  başkalarının eserlerinden yararlanılması durumunda ilgili eserlere bilimsel normlara uygun olarak atıfta bulunduğumu,  atıfta bulunduğum eserlerin tümünü kaynak olarak gösterdiğimi,  kullanılan verilerde herhangi bir tahrifat yapmadığımı,  ve bu tezin herhangi bir bölümünü bu üniversite veya başka bir üniversitede başka bir tez çalışması olarak sunmadığımı beyan ederim. 27/06/2019 Bahadır Erman YÜCE ÖZET Doktora Tezi OFİS ORTAMINDA TÜRBÜLANSLI AKIŞ MODELLERİNİN KARŞILAŞTIRILMASI VE ÇEŞİTLİ ISITMA SİSTEMLERİNİN ISIL KONFOR İNDİSLERİ VE GAGGE MODELİ KULLANILARAK SAYISAL SİMÜLASYONU Bahadır Erman YÜCE Bursa Uludağ Üniversitesi Fen Bilimleri Enstitüsü Makine Mühendisliği Anabilim Dalı Danışman: Doç. Dr. Erhan PULAT Bu çalışmada iki ve üç boyutlu bir oda geometrisinde alttan ısıtma, radyatör ile ısıtma ve klima ile ısıtma durumunda ısıl konfor Hesaplamalı Akışkanlar Dinamiği metodu ile incelenmiştir. Oda geometrisi için Bursa Uludağ Üniversitesi Makine Mühendisliği Binası ofis odalarından biri kullanılmıştır. Oda içerisinde oturma pozisyonunda bir manken modellenmiştir. Sayısal çalışmada kullanılan yöntemlerin doğruluğu test etmek için deneysel verilerin mevcut olduğu teyit amaçlı doğrulama çalışmaları yapılmıştır. k-ε modelleri (Standard, Re-Normalisation Group, Realizable), k-ω modelleri (Standard, shear stress transport, baseline) ve literatürde daha önce yapılmış BGS modelleri birbiri ile kıyaslanmıştır. Doğal, zorlanmış ve birleşik taşınım akışları ayrı ayrı test edilmiştir. Yapılan çalışmalar sonucunda Standard k-ε türbülans modeli deneysel verilerle uyumluluğu nedeniyle tez çalışmasında kullanılmak için seçilmiştir. Oda içinde akış zamana bağlı, sıkıştırılabilir (ideal gaz), birleşik ve türbülanslı olup korunum denklemleri ANSYS Fluent 18.0 yazılımı kullanılarak çözülmüştür. Manken yüzey sıcaklığını elde etmek için manken yüzeyi on altı farklı parçaya bölünmüştür. Manken cilt sıcaklığı, Gagge ısıl duyarlılık modeli baz alınarak geliştirilen C programlama dilindeki kod yardımıyla her bir zaman adımı için dinamik olarak hesaplanmıştır. Manken yüzey sıcaklığı çevre hava sıcaklığı ve hızına duyarlıdır. Her ısıtma sistemi için sıcaklık, hız, PMV, PPD ve PD dağılımları elde edilmiştir. Klima ile ısıtma durumunda diğer durumlardan farklı olarak dokuz farklı üfleme sınır şartında çözümler elde edilmiştir. Klimanın üfleme hızı ve sıcaklığındaki değişimlerin odanın ısınmasına ve konfor indislerinin değişimine olan etkisi incelenmiştir. Anahtar Kelimeler: Hesaplamalı akışkanlar dinamiği, türbülans, ısıl konfor, Gagge modeli 2019, x + 177sayfa. i ABSTRACT PhD Thesis COMPARISON OF TURBULENT FLOW MODELS IN OFFICE MEDIUM AND NUMERICAL SIMULATION OF VARIOUS HEATING SYSTEMS BY USING THERMAL COMFORT INDICES AND GAGGE MODEL Bahadır Erman YÜCE Bursa Uludağ University Graduate School of Natural and Applied Sciences Department of Mechanical Engineering Supervisor: Assoc. Prof. Erhan PULAT In this study thermal comfort was investigated by the Computational Fluid Dynamics method in two and three-dimensional room geometry, and in case of underfloor heating, radiator heating and heating with air conditioner. One of the office rooms of Bursa Uludag University Mechanical Engineering Building was used for the geometry. A thermal manikin is modeled in the sitting position in the room. In order to test the accuracy of the methods used in the numerical study, some benchmarks studies which contain experimental data, were performed. k-ε models (Standard, Re-Normalisation Group, Realizable), k-ω models (Standard, shear stress transport, baseline) and BGS models which previously performed to benchmark studies in the literature, were compared with each other. Natural, forced and combined convection were tested separately. As a result of the studies, Standard k-ε turbulence model was chosen to be used in the thesis study because of its compatibility with experimental data. The flow in the room is time- dependent, compressible (ideal gas), mixed turbulent and the conservation equations are solved using ANSYS Fluent 18.0 software. To obtain the manikinsurface temperature, the manikinsurface is divided into sixteen different parts. The model skin temperature was calculated dynamically for each time step with the help of the code in C programming language developed based on Gagge thermal sensitivity model. Manikinsurface temperature is sensitive to ambient air temperature and speed. Temperature, velocity, PMV, PPD and PD distributions were obtained for each heating system. In the case of heating with air conditioner, different solutions have been obtained in nine different inlet boundary conditions. The effect of air conditioning blowing speed and temperature changes on air temperature in the room and comfort indices were investigated. Key words: Computational fluid dynamics, turbulence, thermal comfort, Gagge model 2019, x + 177 pages. ii TEŞEKKÜR Bana Hesaplamalı Akışkanlar Dinamiği ile tanışma fırsatı veren, bilgi ve deneyimi ile çalışmalarımı yönlendiren, değerli danışmanım Doç. Dr. Erhan PULAT’a, teşekkürlerimi sunarım. Doktora tez çalışmam süresince desteklerini esirgemeyen Doç. Dr. Nurullah ARSLANOĞLU ve Öğr. Gör. Dr. Mustafa MUTLU’ya teşekkürlerimi sunarım. Doktora tezimi hazırlarken geçirdiğim yoğun çalışma temposuna rağmen benden desteklerini esirgemeyen, hep yanımda olan ve motive eden değerli eşim Nurten YÜCE’ye ve aileme içten teşekkürlerimi sunarım. Bahadır Erman YÜCE 27/06/2019 iii İÇİNDEKİLER Sayfa ÖZET.................................................................................................................................. i ABSTRACT ....................................................................................................................... i TEŞEKKÜR ..................................................................................................................... iii İÇİNDEKİLER ................................................................................................................ iv SİMGELER ve KISALTMALAR DİZİNİ ....................................................................... v ŞEKİLLER DİZİNİ ........................................................................................................ viii ÇİZELGELER DİZİNİ ..................................................................................................... x 1.GİRİŞ…….… ................................................................................................................ 1 2. KAYNAK ARAŞTIRMASI ......................................................................................... 6 3. MATERYAL VE YÖNTEM ...................................................................................... 10 3.1. HAD Tarihi .............................................................................................................. 10 3.2 Korunum Denklemleri .............................................................................................. 11 3.3. Türbülans ve Türbülans Modelleme ........................................................................ 13 3.3.1 k-ε Türbülans Modelleri ......................................................................................... 16 3.3.2 k-ω Türbülans Modelleri ........................................................................................ 24 3.4. Yakın Duvar Yaklaşımı............................................................................................ 33 3.5. Isıl Konfor ................................................................................................................ 38 4.BULGULAR ................................................................................................................ 54 4.1. Validasyon Çalışmaları ............................................................................................ 54 4.1.1. Zorlanmış Taşınım ................................................................................................ 56 4.1.2. Doğal Taşınım ....................................................................................................... 65 4.1.3. Birleşik Taşınım .................................................................................................... 70 4.2. Bir Ofis Odasında Alttan Isıtma Sisteminin Isıl konfor Aaçısından İncelenmesi ... 78 4.3. Bir Ofis Odasında Radaytör ile Isıtma Durumunun Isıl konfor Aaçısından İncelenmesi ..................................................................................................................... 87 4.4. Bir Ofis Odasında Klima ile Isıtma Durumunun Isıl konfor Aaçısından İncelenmesi ......................................................................................................................................... 94 5. SONUÇ VE TARTIŞMA ......................................................................................... 107 KAYNAKLAR ............................................................................................................. 110 EKLER .......................................................................................................................... 117 ÖZGEÇMİŞ .................................................................................................................. 176 iv SİMGELER ve KISALTMALAR DİZİNİ Simgeler Açıklama A Alan (m 2) A DuBois yüzey alanı (m 2) D a İvme (m/s2) b Katsayı (mm) C Gagge modelinde taşınım ile olan ısı kaybı (W/m2) c Sabit basınçta özgül ısı (J/kgK) p c Vücudun özgül ısısı (kJ/kgK) p,b c Kanın özgül ısısı (kJ/kgK) p,bl C Solunum ile olan taşınımla ısı kaybı (W/m 2) res CSIG Kordan gelen soğuk sinyal (boyutsuz) cr CSIG Deriden gelen soğuk sinyal (boyutsuz) sk c Sabit hacimde özgül ısı (J/kgK) v d Çap (m) df Serbestlik derecesi D Hidrolik çap (m) H E Buharlaşma ile olan ısı kaybı (W/m 2) E Difüzyon ile olan buharlaşma ısı kaybı (W/m 2) dif E Maksimum buharlaşma ile olan buharlaşma ısı kaybı (W/m 2) max E Solunumla olan buharlaşma ısı kaybı (W/m 2) res E Terleme ile olan buharlaşma ısı kaybı (W/m 2) rsw E Deriden toplam buharlaşma ısı kaybı (W/m 2) sk F Fisher katsayısı (boyutsuz) f Giysi faktörü cl g Yerçekimi ivmesi (m/s2) Gr Grashof sayısı h Taşınımla ısı geçiş katsayısı (W/m 2K) c h Suyun buharlaşma gizli ısısı (kJ/kg) fg l Giysi yalıtımı (clo) cl i Vücut parçası numarası (boyutsuz) j Giysi katmanı numarası (boyutsuz) k Giysi tabakaları arasındaki havanın ısıl iletkenliği (mm.W/m2K) k Türbülans kinetik enerjisi (m2/s2) K İç vücut ile deri arasındaki etkin iletim sayısı (W/m 2K) l Boy (m) l Uzunluk ölçeği (m) LR Lewis oranı (°C/kPa) v m Kütle (kg) m Kütlesel debi (kg/s) M Toplam metabolik ısı üretimi (W/m 2) M Aktiviteye bağlı metabolik ısı üretimi (W/m 2) act m Kan akışı (kg/m 2s) bl m Birim zamanda birim cilt yüzey alanında üretilen ter miktarı (kg/m 2s) rsw MS Beklenen varyans (boyutsuz) M Titreme ile oluşan metabolik ısı üretimi (W/m 2) shiv nl Toplam giysi tabakaları sayısı (boyutsuz) Nu Nusselt sayısı (boyutsuz) p Basınç (Pa) P Ortam havasının su buharı basıncı (kPa) a Pr Prandtl sayısı (boyutsuz) P Ortalama deri sıcaklığındaki su buharı basıncı (kPa) sk ,m P Derideki su buharı basıncı (kPa) sk ,s q Kordan deriye olan toplam ısı transferi (W/m 2) cr ,sk r Yarıçap (m) Ra Rayleigh sayısı (boyutsuz) R Dış hava tabakası ısıl direnci (m 2°C/W) a R Giysi tabakalarının arasındaki havanın ısıl direnci (m 2°C/W) al Re Reynolds sayısı (boyutsuz) R Dış hava tabakası buharlaşma direnci (m 2kPa/W) e,a R Giysi tabakalarının arasındaki havanın buharlaşma direnci e,al (m2kPa/W) R Kumaşın buharlaşma direnci (m 2kPa/W) e, f R Toplam buharlaşma direnci (m 2kPa/W) e,t R Kumaşın ısıl direnci (m 2°C/W) f R Toplam ısıl direnç (m 2°C/W) t S Deride depolanan ısıl enerji (W/m 2) sk T Sıcaklık (K, °C) t Zaman (s) T Ortam sıcaklığı (°C) a T Vücut ortalama sıcaklığı (°C) b,m T Vücut ortalama nötr sıcaklığı (°C) b,n t Giysi yüzey sıcaklığı (°C) cl T Kor tabakası sıcaklığı (°C) cr T Ortalama kor tabakası sıcaklığı (°C) cr ,m T Nötr kor tabakası sıcaklığı (°C) cr ,n T Operatif sıcaklık (°C) o vi T Işınım sıcaklığı (°C) r T Deri tabakasının sıcaklığı (°C) sk T Vücut ortalama deri tabakası sıcaklığı (°C) sk ,m T Nötr deri tabakası sıcaklığı (°C) sk ,n u X yönünde hız değeri (m/s) v Y yönünde hız değeri (m/s) V Hava hızı (m/s) w Deri ıslaklığı (boyutsuz) w Z yönünde hız değeri (m/s) w Difüzyon kaynaklı deri ıslaklığı (boyutsuz) dif w Terin buharlaşması için gerekli olan deri ıslaklığı (boyutsuz) rsw W Yapılan iş (W/m2) WSIG Vücuttan gelen ılık sinyal (boyutsuz) b WSIG Kordan gelen ılık sinyal (boyutsuz) cr WSIG Deriden gelen ılık sinyal (boyutsuz) sk x Giysi katmanları arasındaki hava tabakası kalınlığı (mm) V Hacim (m3)  Isıl yayınım katsayısı (m2/s) a Deri bölgesinde bulunan vücut kütlesi (boyutsuz)  Isıl genleşme katsayısı (m2/s)  Değişkenin son durum ile başlangıç durumu farkı  Delta fonksiyonu  Türbülans yayınım oranı (m2/s3)  Dinamik viskozite (kg/ms)  Kinematik viskozite (m2/s)  Yoğunluk (kg/m3)  Kayma gerilmesi (Pa)  Özgül yayınım oranı (1/s)  Zaman (s) Kısaltmalar Açıklama HAD Hesaplamalı akışkanlar dinamiği PMV Tahmin edilen ortalama oy (Predicted Mean Vote) PD Yüzde memnuniyetsizlik (Percent Dissatisfaction) PPD Tahmin edilen yüzde memnuniyetsizlik (Predicted Percentage of Dissatisfied) vii ŞEKİLLER DİZİNİ Sayfa Şekil 3.1 HAD gelişimine etki eden faktörlerin kronolojik olarak gösterimi ................. 11 Şekil 3.2 Doğrudan sayısal simülasyon (DSS) modeliyle modellenen büyük ve küçük ölçekli girdaplar (Çengel ve Cimbala 2006) ................................................................... 13 Şekil 3.3 BGS simülasyonunda modellenen büyük ölçekli girdaplar (Çengel ve Cimbala 2006) ............................................................................................................................... 14 Şekil 3.4 Hızın zamana bağlı dalgalanması ve ortalama çalkantı terimlerinin şematik gösterimi (Çengel ve Cimbala 2006) .............................................................................. 15 Şekil 3.5. Y+ değerine bağlı olarak sınır tabaka katmanları (ANSYS Fluent Theory Guide, 2016) .................................................................................................................... 34 Şekil 3.6 Duvar fonksiyonu ve yakın duvar modeli çözümlerinin gösterimi (ANSYS Fluent Theory Guide 2016) ............................................................................................. 35 Şekil 3.7 PMV ölçeğine göre PPD değerlerinin dağılımı (ASHRAE 2009) .................. 41 Şekil 3.8 İnsan ve çevresinin eş merkezli silindir modeli (ASHRAE 2009) .................. 42 Şekil 3.9 Kullanılan HAD yazılımının ve hazırlanan kodu eş zamanlı çalışma ilişkisi. 53 Şekil 4.1 IEA Annex 20 odasının ölçüleri ve sınır şartları ............................................. 57 Şekil 4.2 Zorlanmış taşınım için a) x=H, b) x=2H, c) y=0,084m, d) y=2,916m kesitlerinde hız dağılımları ................................................................................................................. 57 Şekil 4.3 Akım çizgileri (Bu çalışma: ANSYS-Fluent 18.0; Rong ve Nielsen (2008): ANSYS-CFX 11.0; Dreau ve ark. (2012): STAR-CCM+; Pulat and Ersan (2015): ANSYS-CFX 13.0 Non-isothermal) ............................................................................... 59 Şekil 4.4 Boyutsuz hızın y=0.084 m’de x bileşeni for a) k-ε modelleri and b) k-ω modelleri ve Boyutsuz hızın y=2.916 m’de x bileşeni c) k-ε modelleri d) k-ω modelleri ......................................................................................................................................... 62 Şekil 4.5 Boyutsuz hızın x=H’da x bileşeni for a) k-ε modelleri, b) k-ω modelleri, c) Hız vektörleri ve x=H ve x=2H’da SST k-w modeliyle elde edilmiş U/U0 boyutsuz hız profili, x=2H’da boyutsuz hız profilleri a) k-ε modelleri and b) k-ω modelleri ......................... 63 Şekil 4.6 Nielsen (1990)’in, Davidson ve Nielsen (1996)’in SGS BGS tahminleri ve Emmerich ve McGrattan (1998)’ın zorlanmış taşınımda x=H’da boyutsuz hız profillerinin deneysel verilerle karşılaştırılması a) k-ε ve b) k-ω modelleri (Smagorinsky sabiti her BGS modeli için Cs = 0.14 alınmıştır.) ........................................................... 64 Şekil 4.7 Doğal taşınım için çalışılan geometrinin ölçüleri ve sınır şartları ................... 65 Şekil 4.8 Farklı konumlar için dört farklı ağ yapınsa sıcaklık dağılımları...................... 66 Şekil 4.9 Farklı konumlar için dört farklı ağ yapında hız dağılımları ............................. 66 Şekil 4.10 Doğal taşınımda farklı türbülans modellerinde elde edilen akım çizgileri ... 67 Şekil 4.11 Sıcaklık ve hız dağılımları a) y/H = 0.1 b) y/H = 0.5 c) y/H =0.9 ................. 69 Şekil 4.12. Doğal taşınımda y/H=0.1’de dikey hız profillerinin karşılaştırılması a) RANS modelleri b) BGS modelleri (Yang ve ark., 2015), (Vreman BGS: Vreman model of BGS, Lily BGS: Lily model of BGS, Wale BGS: Wall-Adapting Local Eddy-Viscosity SGS Model of BGS, Wmles BGS: Wall Modelled BGS) ....................................................... 70 Şekil 4.13 Zorlanmış taşınım için kullanılan geometrinin ölçüleri ve sınır şartları ........ 71 Şekil 4.14 Zorlanmış taşınımda farkı eleman sayıları için a) x/L konumunda sıcaklık dağılımı, b) y/L konumunda sıcaklık dağılımı, c) x/L konumunda hız dağılımı, d) y/H konumunda hız dağılımı sonuçları .................................................................................. 72 viii Şekil 4.15 a) Literatürde elde edilmiş hız vektörleri ve akım çizgileri b) k-ε modelleri ile elde edilmiş hız vektörleri ve akım çizgileri c) k-ω modelleri ile elde edilmiş hız vektörleri ve akım çizgileri. ............................................................................................................. 73 Şekil 4.16 x/L = 0.5’de a) k-ε ve b) k-ω modellerinde hız ve sıcaklık dağılımları. ...... 75 Şekil 4.17 x=L/2’de BSL k-ω’inde elde elde edilen hız profili ve hız vektörleri ........... 76 Şekil 4.18 y/H = 0.5’te birleşik taşınımda sıcaklık ve hız profillerinin karşılaştırılması a) k-ε and b) k-ω modelleri ................................................................................................. 77 Şekil 4.19 İki boyutlu alttan ısıtma çalışması için oluşturulan model ve sınır şartlarının atandığı bölgeler .............................................................................................................. 79 Şekil 4.20 Ofis odasının geometrik özellikleri ................................................................ 80 Şekil 4.21 Alttan ısıtmanın incelendiği I) iki ve II) üç boyutlu çalışmada, mankeni ortalayan kesitte (z=1,5 m) a) Sıcaklık, b) Hız, c) PMV, d) PPD ve e) PD dağılımları. 84 Şekil 4.22 Radyatör ile ısıtmanın incelendiği iki boyutlu ve üç boyutlu çalışmada, mankeni ortalayan kesitte (z=1,5 m), a) Sıcaklık, b) Hız, c) Karbondioksit hacim oranı, d) PMV, e) PPD ve f) PD dağılımları. ............................................................................ 89 Şekil 4.23 Alttan ısıtma ve radyatör ile ısıtma durumunda vücut parçalarının sıcaklıkları. ......................................................................................................................................... 91 Şekil 4.24 Klima ile ısıtmanın incelendiği I) iki ve II) üç boyutlu çalışmada 9 farklı durum için sıcaklık dağılımları. .................................................................................................. 96 Şekil 4.25 Alttan ısıtma, radyatör ve klima ile ısıtmanın incelendiği iki boyutlu ve üç boyutlu çalışmalarda elde edilen ortalama oda sıcaklıkları. ........................................... 98 Şekil 4.26 Klima, radyatör ve alttan ısıtmanın kullanıldığı iki boyutlu çalışmalarda 16 farklı vücut parçasının sıcaklıkları. ................................................................................. 98 Şekil 4.27 Klima ile ısıtmanın incelendiği I) iki ve II) üç boyutlu çalışmada 9 farklı durum için hız dağılımları. ....................................................................................................... 100 Şekil 4.28 Klima ile ısıtmada dokuz farklı durum için a) ortalama sıcaklık ve b) ortalama hız değerleri ................................................................................................................... 101 Şekil 4.29 Klima ile ısıtmanın incelendiği I) iki ve II) üç boyutlu çalışmada 9 farklı durum için PMV dağılımları. ................................................................................................... 102 Şekil 4.30 Klima ile ısıtmanın incelendiği I) iki ve II) üç boyutlu çalışmada 9 farklı durum için PPD dağılımları. ..................................................................................................... 103 Şekil 4.31 Klima ile ısıtmanın incelendiği I) iki ve II) üç boyutlu çalışmada 9 farklı durum için PD dağılımları. ....................................................................................................... 104 ix ÇİZELGELER DİZİNİ Sayfa Çizelge 3.1 16 farklı vücut parçası için alan ve ağırlık değerleri (Tanabe ve ark. 2002) ......................................................................................................................................... 44 Çizelge 3.2 16 farklı vücut parçası için kor sıcaklıkları (Tanabe ve ark. 2002) ............. 45 Çizelge 3.3 Vücut parçaları üzerindeki giysilerin özellikleri (Arslanoğlu 2015, Mccullough ve ark. 1989) ............................................................................................... 46 Çizelge 0.1 GCI değeri ve değere ulaşmak için yapılan hesaplamalar ........................... 83 Çizelge 0.2 Alttan ısıtma sisteminin incelendiği durumda vücut parçaları üzerindeki sıcaklık değerleri. ............................................................................................................ 85 Çizelge 0.3 AYI değeri ve değere ulaşmak için yapılan hesaplamalar ........................... 88 Çizelge 0.4 Alttan ısıtma sisteminin incelendiği durumda vücut parçaları üzerindeki sıcaklık değerleri. ............................................................................................................ 90 Çizelge 0.5 Klima ile ısıtma çalışmasında klimadan üflenen hava için sınır şartları. .... 94 Çizelge 0.6 AYI değeri ve değere ulaşmak için yapılan hesaplamalar ........................... 95 x 1. GİRİŞ Binalar çağlar boyunca değişken hava şartlarına (sıcak hava ve soğuk hava, rüzgar, yağmur, kar vb.) karşı koruma sağlamakta ve bina sakinlerine konforlu bir iç ortam sağlamaktadır. Binaların koruma kapasitesi ve becerisi ise geçmişten günümüze sürekli olarak artmış ve değişen ihtiyaçlara göre yenilenmiştir. Başlangıçta, iklimlendirme işlemi herhangi bir enerji tüketimi olmaksızın bina tasarımından elde edilen pasif soğutma ile sağlanıyordu (Samuel ve ark. 2013). Binayı yapan ve binada yaşayan kişi genelde aynı olduğu için konfora verilen önem de daha fazlaydı. Tarih öncesi çağlardan beri, kar ve buz soğutma için kullanılmıştır. Kış aylarında buz toplama ve yaz aylarında kullanım için depolama işi, 17. yüzyılın sonlarına doğru popüler bir meslek haline gelmiştir (Nagengast 1999). Literatürde iklimlendirmenin ardındaki ilk temel konseptin, sazların pencerelere asıldığı ve odaların damlayan suyla nemlendirildiği eski Mısır'da uygulandığı ifade edilmektedir. Pencereden içeri giren sıcak hava, damlayan suyun buharlaşması ile soğuyarak iç ortamı serinletmiştir. Antik Roma'da, su kemerleri aynı zamanda evleri de soğutması için evlerin etrafından dolaştırılmıştır (Bahadori 1978). 1840'larda, Florida’lı mucit Dr. John Gorrie, bölgesinin önemli sorunlarından biri olan yüksek sıcaklık için bir fikir geliştirdi. Gorrie, soğutmanın sıtmaya karşı bir tedavi ve hastaları rahatlatmanın bir yolu olduğunu düşünüyordu. Ancak hastanelere soğutmayı sağlayacak olan buz kuzeyin uzak dağlarından, nehirlerden ve göllerden sağlanıyordu. Bu işlem ise pahalı bir lojistik hizmeti gerektiriyordu. Gorrie, bu zorluk ve maliyetin üstesinden gelmek için, yapay soğutma fikri üzerinde çalışmaya başladı. Çalışmalarının sonucunda su, rüzgar, buhar veya hayvan gücü ile çalışan bir kompresör kullanarak buz oluşturan bir makine tasarladı ve 1851'de bu makinenin patentini aldı. Gorrie, mali destekçisinin ölümü üzerine bu teknolojiyi pazara sunmada başarısız olmuştur fakat icadı modern klima ve soğutma için temel oluşturmuştur (Lester, 2019). 1 1902 yılında, ilk modern elektrikli klima ünitesi Buffalo, New York'ta Willis Carrier tarafından icat edildi. Carrier, Brooklyn, New York'taki Sackett-Wilhelms Yayıncılık Şirketi için bir dergi kapaklarının kırışması ile alakalı bir mühendislik problemini çözmek için iklimlendirme denemelerine başladı. Carrier tarafından Buffalo'da tasarlanan ve inşa edilen ilk klima, 17 Temmuz 1902'de çalışmaya başladı. Bir dizi deneyle, Carrier, soğutma serpantini kullanarak “Apparatus for Treating Air” adını verdiği, nemi kontrol eden bir sistem tasarladı. Teknolojisini test ederek ve geliştirerek, tekstil fabrikalarında havanın nemini ve sıcaklığını düzenleyen bir otomatik kontrol sistemi tasarlayarak patentini aldı. Carrier, nem kontrolü ve iklimlendirmenin diğer birçok endüstriye fayda sağlayabileceğini fark etmişti ve sonunda Buffalo Forge'dan ayrılarak The Carrier Air Conditioning Company of America’yı kurdu. Bu cihaz ilk olarak kamu binaları ve ardından tiyatro ve sinemalarda kullanılmaya başlandı. Başlarda boyut ve maliyeti yüksek olan bu klima yıllar içinde geliştirilmeye devam etti. Soğutma teknolojilerindeki gelişmelere rağmen, bu sistemler evler için çok büyük ve pahalıydı. Frigidaire, 1929’da evler için kullanıma uygun ilk split tip klimayı duyurdu. Fakat sistem hala ağır, pahalı ve uzaktan kontrol edilen bir yoğuşma ünitesine ihtiyaç duyuyordu. 1947’de gelişen teknoloji ile daha kompakt hale getirilmiş ve birçok farklı şirketin üstünde çalıştığı klimalar önemli satış rakamlarına ulaşmıştır. 20. yüzyıldan itibaren modern bina inşaat teknolojisi tanıtılmıştır. Bu yeni teknolojilerin ve çelik gibi yapı malzemelerin kullanımı daha yüksek ve daha derin binaların yapımını da mümkün kılmıştır. Yapay ışıklandırma, yalıtım malzemeleri ve havalandırmayla birlikte artık dış ortam iklim şartlarından neredeyse tamamen bağımsız yeni bir iç mekan havası oluşturulmuştur. Artan teknolojik gelişmeler ve modern binaların sayısının artması yapı endüstrisinde uzmanlaşmaya yol açmıştır. Sonuç olarak binalar yalnızca bina sakinlerinin ihtiyacı için değil çok daha geniş kapsamlı işlevler için inşa edilmeye başlanmıştır. 2 Günümüzde pek çok insan zamanının çoğunu kapalı mekanlarda ve paylaşılan ortak alanlarda geçirmektedir. Havalandırma sistemlerinin devreye girmesi ve buna bağlı olarak iç mekan ikliminin artan kontrolü sayesinde, insanların iç mekan havası ve iklimine olan beklentileri artmıştır. Bununla birlikte, bu beklentilerin yerine getirilmesi, bina içindeki termal konfor koşullarının dar bir aralıkla (Povl Ole Fanger 1970) sınırlı olması ve bu aralığın her birey için farklılık göstermesiyle zorlaşmaktadır. Kullanılan yaygın iklimlendirme sistemlerinin varlığına rağmen, havalandırma hakkındaki şikâyetler hâlâ önemli bir problemdir. Akademik anlamda ısıl konfor, yirminci yüzyılın başlarında, Adolf Pharo Gagge'in işyerindeki stresli durumlardan kaynaklanan belirli problemleri çözmek için yaptığı çalışmalarla doğmuştur. 1970’lerde Povl Ole Fanger ve diğer araştırmacılar, çalışmayı gerçek bir disiplin haline getirmiştir (Fabbri 2015). Isıl konfor, ısıl ortamdan memnuniyeti ifade eden öznel bir ifadedir çünkü hem fizyolojik hem de psikolojik olarak, insandan insana büyük farklılıklar olduğu için, belirli bir ortamdaki tüm insanları tatmin etmek çoğu zaman mümkün değildir. Konfor için gerekli olan çevresel koşullar herkes için aynı değildir (Ashare Standard-55, 2004). Termal konfora ilave olarak, iç hava kalitesi de (Indoor Air Quality - IAQ) çözülmesi gereken bir mühendislik problemi haline gelmiştir. Bu sorun, enerji tasarrufunun diğer tüm ihtiyaçların önüne geçtiği 1970'lerin başındaki enerji krizinden kaynaklandı. Bina cepheleri bu yıllardan itibaren, hava sızıntısını önlemek için kademeli olarak daha iyi yalıtılmıştır. Doğal veya lokal mekanik havalandırma sistemleri merkezi mekanik havalandırma ile değiştirildi. Bir başka önlem olarak, hava akış hızı azaltıldı. Bu önlemler farklı bir mühendislik probleminin ortaya çıkmasına sebep oldu. Düşük havalandırma kapasiteli merkezi kontrollü bir sistem, kabul edilebilir koşulların sağlanmasında yeterli olmayabilir. Bu sistemlerde temiz hava doğrudan kişiye verilemez. Dahası, doğal havalandırma yoluyla (yani, bir pencere açarak) koşulların düzeltilmesi, kapalı cephe nedeniyle artık mümkün değildir. Bunlara ek olarak, binaların giderek daha da sızdırmaz hale getirilmesi ve temiz havanın azalması sonucu, modern yapı malzemeleri, mobilya ve temizlik ürünlerinden kaynaklanan gazların iç mekandaki yoğunluğu arttı. Bunun 3 sonucunda havalandırma kaynaklı sağlık problemlerinde giderek artan bir oran ortaya çıktı. Bina içerisindeki çalışanların üretkenliği ve verimliliği termal çevrenin koşullarıyla yakından ilişkilidir (Akimoto ve ark. 2010, Shin-ichi Tanabe ve ark. 2015). İyi bir iç hava kalitesine sahip çalışma ortamı, her modern ekonomide gayri safi milli hasılanın (GSMH) önemli bir bölümünün ofis binalarında çalışan insanlar tarafından kazanıldığını düşünülürrse önemi daha da artacaktır. Tasarımda temel zorluk, termal konfor ve iç hava kalitesinin önemi göz önüne alındığında, bir bina ve kişi için kabul edilebilir kapalı çevre koşullarının sağlanmasıdır. Böyle bir tasarım yaklaşımı, birey bazında, konfor odaklı tasarım gerektirirken mevcut durumda bunun zıttı, tekdüze bir tasarım anlayışı güdülmektedir. Oda yerine oturan kişiye odaklanmak, bina iç hava akışı yapısı ile ilgili detaylı bilgi gerektirir. Hava değişim hızı, temiz hava girişi ve temiz hava konsantrasyonunu etkileyen diğer gazların miktarı gibi bilgiler, ısıl konfor ve iç hava kalitesini değerlendirmek için yeterli değildir. Ayrıca iç ortam havası tek bir bölge olarak kabul edilmemelidir. Hava akışı ve buna bağlı olan tüm parametreler kontrol hacminin tamamında incelenmelidir. İç ortam hava akışı yapısının önemi hakkındaki bu bilinç göz önüne alındığında, tasarım aşamasında akış özelliklerini belirlemek ve tahmin etmek için güvenilir araçlara ihtiyaç duyulmaktadır. Bu araçlardan en önemlilerinden biri de Hesaplamalı Akışkanlar Dinamiği metodudur (HAD) (Nielsen 2015). İç hava akışını belirlemek ve tahmin etmek için mevcut olan seçenekler geçmişten günümüze sürekli gelişim göstermiştir. Bu gelişimin en önemli sebebi ise bilgisayar donanım teknolojilerinin hızlı gelişim sürecidir. Bir odadaki hava akışı genellikle izotermal olmayan, türbülanslı, üç boyutlu ve zamana bağlı olarak karakterize edilir. Bu karakterizasyon ise çözülmesi gereken daha fazla denklem ve bunun sonucunda yüksek kapasiteli bilgisayar donanımları ve zaman maliyeti ihtiyacı ortaya çıkarmaktadır. Günümüzde türbülans modelleme ve akış probleminin hesaplanmasında kullanılan çözücü tekniklerindeki gelişmeler simülasyon sürecini önemli ölçüde geliştirmiştir. 4 Bugünkü bilgisayarlarda gün geçtikçe artan işlemci hızı ve bellek kapasitesi, havalandırma ve diğer birçok mühendislik alanında HAD kullanımını önemli ölçüde kolaylaştırmış ve yaygınlaştırmıştır. Günümüzde sayısal çözümlerdeki tolere edilebilir hatalarla birlikte deneysel olarak ölçüm yapmanın zor ve/veya bazı durumlarda imkansız olduğu alanlarda HAD metodu Nielsen’in (1974) HAD’ı kullanarak iç ortam hava akışını incelenmesi deneysel çalışmalara büyük ve cazip bir alternatif sunmuştur. HAD, kapalı hacmin her tarafında yüksek doğrulukta ve detaylı matematiksel veriler vermektedir. HAD ile oda içerisindeki ısıl konfor indisleri, sıcaklık, basınç, hız, nem ve kirletici verileri ihmal edilebilir hata payları ile bulunabilmektedir. Bu çalışmada kapalı hacimlerde farklı ısıtma sistemlerinde gerçekleşen akışlar ısıl konfor dikkate alınarak incelenmiştir. Doğal, zorlanmış ve birleşik ısı transferinin meydana geldiği türbülanslı akışlar için kullanılan Fluent yazılımı literatürde var olan teyit amaçlı karşılaştırma (benchmark) çalışmaları ile doğrulanmıştır. Uygun türbülans modeli seçilerek incelenecek olan ısıtma sisteminde kullanılmıştır. 5 2. KAYNAK ARAŞTIRMASI Bu bölümde, tez kapsamında araştırılan konuyla ilişkili olarak literatürde yapılmış olan nümerik çalışmalar özetlenmiştir. Kapalı hacimlerde havalandırma, radyatörle ve alttan ısıtma sistemleri ısıl konfor odaklı olarak incelenmiştir. Bölüm sonunda literatürde eksik kalan konular ve bu çalışmanın literatüre katkısı vurgulanmıştır. 1950'lerden 1970'lere kadar, havacılık endüstrisi ve meteoroloji tahminlerinde HAD’ın kullanımı, HAD’ın erken gelişmesinde başlıca itici güç olmuştur ancak HAD’ın bina havalandırması gibi mühendislik problemlerinde kullanılmasıyla, son 30 yılda HAD gelişiminde ana itici güç mühendislik problemleri olmuştur. Nielsen (1973), kendi tanımı ile “muhtemelen” havalandırma çalışmalarına HAD’ı ilk uygulayan kişiydi. 1970'lerden bu yana, hem HAD’da (daha hızlı bilgisayarlar, daha hızlı sayısal yöntemler ve gelişmiş türbülans modellemesi) hem de havalandırmada (havalandırma yöntemleri ve havalandırma kontrolü) önemli gelişmeler olmuştur. Artık HAD akışkan çalışmalarının yapıldığı hemen hemen tüm disiplinlerde kullanılmaktadır. Tokyo'daki kentsel ısı adası analizi için 33 km x 33 km'lik bir alanda 5 milyar elemanlı HAD simülasyonları yapılmıştır (Ashie ve Kono 2011). Gan (1995) çalışmasında, oda hava dağıtım sistemlerinin tasarımını optimize etmek için iç ortamı etkileyen dört önemli unsuru (hava akışı, termal konfor, hava kalitesi ve enerji kullanımı) incelemiştir. İngiltere’nin kış ve yaz şartlarına göre yaptıkları on üç farklı simülasyonda, havalandırılan bir ofiste hava hareketinin, termal konforun, kirletici dağılımın ve enerji kullanımını HAD metodu kullanarak incelemiş ve birbirleri ile kıyaslamıştır. Çalışmasının sonucunda ısıtma işlemi için en verimli hava dağıtım sisteminin soğutma için aynı etkiyi göstermediğini ifade etmiştir. Isıtma ve soğutma sisteminin ayrık olarak kullanılmasının enerji tasarrufu açısından faydalı olacağını ancak bunun fiziksel olarak çoğu zaman mümkün olmadığını ifade etmiştir. Taze havanın odanın altından yukarı doğru hareket etmesinin hem oda sakinlerine daha fazla temiz havanın ulaşmasını sağladığını hem de enerji tasarrufu için daha uygun olduğunu belirtmiştir. 6 Chen ve Xu (1998) bulundukları dönemde HAD problemleri için önemli bir problem teşkil eden hesaplama maliyeti üzerine bir çalışma yapmıştır ve sıfır denklemi (zero- equation) türbülans modelini sundular. Standart k-ɛ türbülans modeli her ne kadar havalandırma çalışmalarında yaygın olarak kullanılsa da izafi olarak çok fazla bilgi işlem süresi ve yüksek işlem becerisine sahip bilgisayar ihtiyacı gerektirmektedir. Geliştirdikleri model, türbülans viskozitesinin, boyutsuz uzunluğun (length scale) ve lokal ortalama hızın bir fonksiyonu olduğunu varsaymaktadır. Bu model, kapalı hacimlerde doğal, zorlanmış ve birleşik taşınım akışlarında ve deplasmanlı havalandırmayı analiz etmek için kullanılmıştır. Sonuçları deneysel verilerle ve standart k-ɛ modeli ile elde edilen sonuçlarla karşılaştırmış ve makul derecede uyumlu olduklarını ifade etmiştir. Geliştirdikleri sıfır denklemli modelin çözülen denklem sayısındaki azalma sebebi ile daha az bilgisayar belleği kullandığını ve hesaplama hızının k-ɛ modeline göre en az on kat daha fazla olduğunu ifade etmişlerdir. Eleman sayısı ise diğer modellerin gerektirdiğinden daha az olduğu için o dönem bilgisayarlarında üç boyutlu HAD problemlerinin çözümünün çok daha kolay olacağını belirtmişlerdir. Chung ve ark. (1998) ise yaptıkları çalışmada üfleme ve emme menfezlerinin ısıl konfora ve havalandırma verimine olan etkilerini deneysel ve sayısal olarak incelemişlerdir. Türbülans modeli olarak Standart k-ɛ modelini kullanmışlardır. Sayısal sonuçlar için üfleme kısmında yüksek doğrulukta değerler elde ettiklerini ve çıkış bölümü için ise deneysel verilere kıyasla makul değerler elde ettiklerini ifade etmişlerdir. Üfleme ve emme menfezlerinin birbirine dik konumlandırıldığı sistemlerin paralel konumlanan sistemlere oranla daha verimli bir havalandırma sağladığını ve oda içerisine bölmeler eklemenin havalandırmaya çok büyük bir oranda etki ettiğini ifade etmişlerdir. Xing ve ark. ( 2001), Oda içerisindeki termal mankenin soluma çevresi etrafında ve oda içindeki diğer bölgelerde havanın tazeliğini, değişim indeksini ve havalandırma verimliliğini sayısal ve deneysel olarak incelemişlerdir. HAD metodu ile elde ettikleri verilerin deneysel verilerle oldukça tutarlı olduğunu ifade etmişlerdir. Birçok hastanede, bulaşıcı hatalık taşıyan hastalar için izolasyon odaları kullanılmaktadır. Bu özel odalarda, hava dağılımı ve bakteri yayılımı oda havalandırması ile yakından 7 ilgilidir. Shih ve ark. (2007), hareket eden bir kişinin ve otomatik bir sürgülü kapının açılıp kapanmasının hız, basınç ve kirletici parçacıklar dahil olmak üzere oda hava dağılımı üzerindeki etkisini araştırmak için HAD yöntemini kullanmıştır. Yürüyen kişinin ve sürgülü kapının hareketini simüle etmek için dinamik ağ yapısı kullanmışlardır. Sayısal sonuçlara göre, çalışmada hareketli cisimlerin odadaki hava dağılımı üzerindeki etkisi ele alınmıştır. Hız, basınç ve hava dağılımı kişi hareketinden etkilense de manken orijinal konumuna gelip durduğunda bu dağılımlar ilk haline geri dönmektedir. Ayrıca bu hareket kirletici partiküllerin üzerinde kayda değer bir etkisinin olmadığını belirtmişlerdir. Sürgülü kapının açılması ve kapanması ise basınç ve hız dağılımları üzerinde önemli etkilere sahiptir. Sürgülü kapının açılıp kapanması izolasyon odası ile dış oda arasında bir hava akışı oluşturmakta ve bu akım ani basınç değişikliklerine sebep olmaktadır. Ayrıca, belirtilen basınç farkı sayısal olarak ne olursa olsun, izolasyon odasındaki basınç değişimleri aynı eğilimlere sahiptir. Bununla birlikte, kirleticilerin, sadece açılıp kapanma esnasında oluşan basınç farkının büyüklüğünden etkilendiğini ifade etmişlerdir. Kapı kapandıktan sonra ise basınç farklılıkları kirletici hareketinde önemli bir etkiye sahip olmamaktadır. Araştırmacılar çalışmalarında ANSYS Fluent yazılımını kullanmışlardır ayrıca Standard k-ε türbülans modeli ve Standard duvar fonksiyonu yaklaşımını tercih etmişlerdir. Teodosiu ve ark. (2017) oda havalandırması çalışmalarında HAD kullanımının deneysel çalışmalarla olan tutarlılığını göstermek için bir çalışma yapmıştır. Yalıtılmış bir oda oluşturup elde ettikleri hava hızı değerlerini HAD metodu ile elde ettikleri sonuçlarla karşılaştırmışlardır. Çalışmalarında realizable k-ε türbülans modelini ve Güçlendirilmiş Duvar Yaklaşımı (Enhanced Wall Treatment) duvar yaklaşımını kullanmışlardır. Yapılandırılmamış ağ ile ağdan bağımsızlığa ulaşmışlardır. Bunlara ek olarak çalışmaları üç boyutlu ve sıkıştırılamaz akış kabülüne dayanmaktadır. Bu bağlamda bu tez çalışmasına benzer özellikler barındırmakta ve çalışmalarında HAD ile elde ettikleri verilerin deneysel sonuçlarla tutarlığı olduğunu belirtmişlerdir. Literatür incelendiğinde farklı türbülans modellerinin farklı ısı transfer mekanizmalarında verdiği sonuçların gerek deneysel verilerle gerekse literatürle ayrıntılı bir karşılaştırmasına rastlanmamıştır. Özellikle Büyük Girdap Modelleri (BGS) (Large Eddy 8 Simulation) ile yapılan detaylı karşılaştırmalar literatüre katkı sağlamıştır. Bununla birlikte, yazarın bilgisi dahilinde, literatürde Gagge modeli ile yapılan HAD çalışmalarında manken üzerindeki dinamik sınır şartı için oda içerisindeki ortalama hava sıcaklığının ve hızının kullanıldığı görülmüştür. Bu çalışmada ise on altı farklı vücut parçasının her birinin etrafındaki ortalama sıcaklık ve hız değerleri ayrı ayrı hesaplanmış ve model içinde kullanılmıştır. Böylece vücudun etrafında oluşan farklı sıcaklık ve hız dağılımlarına her bir vücut parçasının farklı tepkiler verebilmesi sağlanmıştır. 9 3. MATERYAL VE YÖNTEM Akışın ve ısı transferinin olduğu sistemlerin tasarımı ve analizinde iki farklı yöntem söz konusudur. Bu yöntemler deney yapma ve hesaplamadır (Çengel ve Cimbala 2006). Hesaplama yöntemi ise kendi içinde analitik ve sayısal olarak ayrıştırılabilir. Bu çalışmada sayısal yöntem yani HAD metodu kullanılmıştır. Momentum, ısı ve kütle transferini çözebilmek için momentum, enerji ve türbülans denklemleri ticari bir yazılım olan ANSYS Fluent 18.0 vasıtası ile çözülmüştür. 19.yy başlarında türetilmiş bu kısmi diferansiyel denklemler ve yalnızca sayısal olarak çözülebilmektedir. Kullanılan ticari yazılım koeunum denklemlernin ayrıklaştırılmasında sonlu hacim metodunu kullanmaktadır. Isıl konfor değerlendirmelerini yapmak için ise gerek akademik alanda gerek uluslararası standartlarda kendine yer edinmiş olan Ortalama tahmini oy (Predicted Mean Value - PMV), ortamdan memnun olmayanların yüzdesi (Predicted Percent Dissatisfaction - PPD) ve yüzde memnuniyetsizlik (Percent Dissatisfaction - PD) ısıl konfor indeksleri incelenmiştir. Vücut sıcaklığının ise vücut içi ve çevre sıcaklığı ve ortam hızı gibi pek çok farklı parametreye göre sürekli olarak zaman içinde değiştiği bir model olarak ele alınabilmesi için geçici zaman şartında Gagge Isıl Duyarlılık Modelini (Gagge Thermal Sensation Model) (Gagge ve ark. 1941) kullanılan paket programla beraber eş zamanlı çözüm oluşturması için C programlama dilinde kodlanmıştır. 3.1. HAD Tarihi HAD tarihsel olarak öncelikle doğrusallaştırılmış denklemleri çözmek için geliştirildi. 1930'larda, bir silindir etrafındaki akışkanın akışına uygun verileri sunabilen iki boyutlu (2B) yöntemler geliştirilmiştir (Thomson 1973). Modern HAD'a benzeyen en eski hesaplama türlerinden biri Lewis Fry Richardson tarafından yapılan hesaplardır, bu hesaplamalarda sonlu farklar yönetimini kullanmıştır. Önemli ölçüde başarısız olmalarına rağmen, bu hesaplamalar, Richardson'un " Weather prediction by numerical process" kitabıyla birlikte (Richardson ve Sydney 1965) modern 10 HAD ve sayısal meteorolojinin temelini oluşturmuştur. Bu yıllardan itibaren HAD gelişimine etki eden faktörler kronolojik olarak şekil 3.1’de gösterilmiştir. Bu yıllardan itibaren ise bilgisayarların sürekli olarak artan hesaplama kabiliyeti ile araştırmacıların HAD’a olan ilgisi artmıştır ve bu durum da HAD’ın gelişmesinde önemli bir rol oynamıştır. IV. RANS (1990’lar) +Viskoz Akış III. Euler Yaklaşımı (1980’ler) +Dönme Hareketi II. Nonnineer potansiyel akış +Nonlineerite I. Lineer potansiyel akış (1960’lar) Viskoz olmayan, türbülanssız akışlar Şekil 3.1 HAD gelişimine etki eden faktörlerin kronolojik olarak gösterimi 3.2 Korunum Denklemleri HAD yazılımlarının çözüm aşamasında kullandığı sürekli rejimde üç boyutlu sıkıştırılamaz akışkanlar için kartezyen koordinatlarda korunum denklemleri Denklem 3.1-3.4’de verilmiştir. HAD teorisi ile ilgili denklemler için Fluent Theory Guide (2016) Süreklilik Denklemi;  (3.1)  v   Sm t Denklem 3.1 süreklilik denkleminin genel şeklidir ve sıkıştırılabilir akışların yanı sıra sıkıştırılamaz akışlar için de geçerlidir. Sm ifadesi ise kaynak terimidir. 11 Momentum Denklemi (Batchelor 1967);   (3.2) v   .vv   p .  g  F t    Denklem 3.2’de p statik basınç,  ise aşağıda ifade edildiği gibi gerilme tensörüdür.  g yerçekimsel kütle kuvveti, F ise kaynak terimidir. Gerilme tensörü;   2  (3.3)    Tv v  .vI  3   Burada ki  , moleküler viskozitedir, I birim tensördür ve eşitliğin sağ tarafındaki ikinci terim, hacim genişlemesinin etkisidir. Enerji denklemi ise Denklem 3.4’te ifade edilmiştir.     (3.4)  (E) v E  p  k T h J   v eff j j  S t   eff  h    Denklem 3.4’te gösterilen E ifadesi, p v2 (3.5) E  h    2 Şeklinde ifade edilir. h ifadesi hissedilir entalpidir ve ideal gazlar için Denklem 3.6’da olduğu gibi tanımlanır: h Y h (3.6) j j Denklem 3.6’da gösterilen h j ifadesi ise Denklem 3.7’de gösterilmiştir. T (3.7) hj   c p , jdt Tref 12 3.3. Türbülans ve Türbülans Modelleme Türbülans, orta ve yüksek Reynolds sayılarında akışkanlarda gözlenen üç boyutlu, kararsız ve rastgele harekettir. Ayrıca türbülansın bu özellikleri sürekli olarak daimi değildir ve üç boyutludur. Bu sebeple türbülanslı akışların çözümü laminer akışlara göre çok daha zordur. Şekil 3.2’de bir akış alanı içindeki en büyüğünden (L) en küçüğüne (η) kadar girdaplar şematize edilmiştir. Şekil 3.2 Doğrudan sayısal simülasyon (DSS) modeliyle modellenen büyük ve küçük ölçekli girdaplar (Çengel ve Cimbala 2006) Türbülanslı bir akışta tüm farklı ölçeklerde akışın özelliklerini nümerik olarak çözebilmek için araştırmacılar Doğrudan Sayısal Simülasyon (Direct Numerical Simulation - DSS) metodunu geliştirmiştir. Bu metot ile akış içindeki farklı ölçeklerdeki girdaplar çözülebilmekte fakat bu çözümün büyük bir hesaplama maliyeti olmaktadır. DSS çözümleri için üç boyutlu ve çok sık ağ yapıları gerekmektedir. Bu tarz çözümler şu an kullanılan daha az hesaplama maliyeti gerektiren Reynolds Ortalamalı Navier Stokes (RONS) (Reynolds Averaged Navier Stokes Equations - RANS) modelleri için dahi yüksek performanslı işlem kabiliyetine sahip bilgisayarlar gerektirmektedir. Bunun yanı sıra Reynolds sayısı arttıkça girdaplar ve aralarındaki ölçek farkı önemli ölçüde artmakta bu durum DSS çözümlerini mevcut donanım kabiliyeti seviyesinde uygulanabilir olmaktan çıkarmaktadır. 13 Bu probleme bağlı olarak araştırmacılar Büyük Girdap Simülasyonu (BGS) (Large Eddy Simulation - BGS) adını verdikleri yeni bir yöntem geliştirdiler. BGS, DSS’ten farklı olarak olarak Şekil 3.3’te gösterilen büyük ölçekli girdapların çözümünü yaparken küçük ölçekli girdapları ise modeller. Bu sebeple daha az hesaplama maliyeti gerektirir fakat yine de gerektirdiği donanım ve zaman maliyeti çok yüksektir. Şekil 3.3 BGS simülasyonunda modellenen büyük ölçekli girdaplar (Çengel ve Cimbala 2006) Tüm bu hesaplama maliyetleri çoğu zaman ve çoğu problem için uygulanabilir olma niteliğini sağlayamamaktadır. Bu sebeple çoğu HAD uygulamasında RONS yaklaşımıyla türbülans modelleri kullanılmaktadır. Türbülans modelleri, büyük veya küçük ölçekte hiçbir girdabı doğrudan çözmez bunun yerine tüm girdapları modeller. Tek bir türbülans modelinin diğer modellere göre tüm problem türleri için üstün olarak nitelendirilmesi bugün için maalesef mümkün değildir. Türbülans modelinin seçimi, akışın fiziği, çözülmek istenen problemin türü, gerekli doğruluk seviyesi, mevcut hesaplama kaynakları ve simülasyon için gerekli zaman miktarı gibi hususlara bağlıdır. Her bir türbülans modelinin kendine göre farklı yetenekleri ve limitleri vardır. RONS modelleri akışı ve akışta oluşan dalgalanmaların özelliklerini inceler. Reynolds ortalamalı eşitliklerde karşımıza ekstra terimler çıkar. RONS modelinde, Navier Stokes denklemlerine ilave çözüm değişkenleri (zaman ortalamalı) olarak Şekil 3.4’te gösterilen ortalama ve dalgalanma hız bileşenleri eklenir. 14 Şekil 3.4 Hızın zamana bağlı dalgalanması ve ortalama çalkantı terimlerinin şematik gösterimi (Çengel ve Cimbala 2006) u  u u ' (3.8) i i i Burada u i ve u ' i ortalama hız ve çalkantıdır. Basınç ve diğer skaler birimler için ise Denklem 3.9’da ki ifade kullanılır.    ' (3.9) Bu ifadeler süreklilik ve momentum denklemlerine aktarıldıklarında sırasıyla Denklem 3.10 ve 3.11’de gösterildiği gibi elde edilir.    ui   0 t x (3.10) i   u   (3.11)   p  u j 2 u ui  u u      ii j     kij   u ' '   iu j  t x j xi x j   x j xi 3 xk  x j Yukarıdaki eşitlik RONS denklemleri olarak ifade edilir. Buradaki ekstra Reynolds  ' ' gerilmesi,  uiu j , terimlerinin modellenmesi gerekmektedir. Bu çalışmada türbülans modelleri olarak iki denklemli türbülans modellri şeklinde ifade edilen ε modelleri: 15  Standard k-ε,  RNG k-ε (Renormalizasyon Grubu Teorisi Tekniği- Renormalization Group Technique),  Realizable k-ε, ve ω modelleri:  Standard k-ω,  SST k-ω (Kayma Gerilmesi Transportu-Shear Stress Transport),  BSL k-ω (Baz Model-Baseline) kullanılmıştır. İki denklemli modeller hesaplama maliyeti açısından DSS, BGS gibi modellere karşılık daha kullanışlıdır. 3.3.1 k-ε Türbülans Modelleri İki denklemli modeller, endüstriyel HAD'da tarihsel olarak en yaygın kullanılan türbülans modelleridir. İki transport denklemini çözerler ve Girdap viskozitesi (Eddy Viscosity) yaklaşımını kullanarak Reynolds gerilmelerini modellerler. ANSYS Fluent'teki Standart k-ε modeli bu yaklaşımı kullanır ve Launder ve Spalding (1974) tarafından önerildiğinden beri hesaplamalı akışkanlar dinamiği uygulamalarında kullanılan modeller arasında en çok kullanılan model olmuştur. Çeşitli türbülanslı akışlar için kararlılığı, ekonomikliği ve makul orandaki doğruluğu, endüstri ve akademide bu yaygın kullanımının sebepleridir. k-ε modellerinin dezavantajları ise, diğer modellere göre ters basınç gradyanları ve sınır tabakası ayrılmasında sonuçların makul değerlere daha uzak olmasıdır. Bu modeller gözlemlere göre genellikle gecikmeli ve daha zayıf bir ayrılma öngörürler. Bu durumda, pürüzsüz yüzeylerden (aerodinamik gövdeler, difüzörler vb.) ayrılan akışlar için aşırı iyimser tasarım değerlendirmeleriyle sonuçlanabilmektedir. Bu nedenle k-ε modelleri dış aerodinamikte yaygın olarak kullanılmaz. Bu çalışma da ise benzer tasarım zorluklarına rastlanmamaktadır ve bu dezavantajın çalışmada elde edilen sonuçlara kayda değer bir etkisinin olmadığı görülmüştür. 16 k-ε Modeli Standart k-ε modeli, türbülans kinetik enerjisi (k) ve türbülans yayınım oranı (ε) transport denklemlerine dayanan bir modeldir (Launder ve Spalding 1972). Türbülans kinetik enerjisi için olan transport denklemi, tam denklemden türetilirken, türbülans yayınım oranı, fiziksel muhakeme kullanılarak elde edilmiştir ve matematiksel olarak kesin karşılığı ile benzerlik göstermemektedir. Türbülans kinetik enerjisi ve yayınım oranı için elde edilen transport denklemleri sırasıyla Denklem 3.12 ve 13’te gösterilmiştir:     u  k  (3.12) k  kui   t  t   G kG b Y M  Sk t xi x j   k  x j      u      2 (3.13)   ui   t  t   C1 (G kC 3Gb )C2  S t xi x j     x j  k k Bu denklemlerde, G k , hesaplanan ortalama hız gradyanlarına bağlı olarak türbülans kinetik enerji üretimini ifade eder. G b , kaldırma etkilerine bağlı türbülans kinetik enerjisi üretimidir. YM , sıkıştırılabilir türbülanslı akışta dalgalı dilatasyonun toplam yayınım hızına etkisini gösterir. C1 ve C2 ise sabittir.  k ve   , sırasıyla k ve ε için türbülanslı Prandtl sayılarıdır. Sk ve S kullanıcı tanımlı kaynak terimlerdir. Türbülans viskozitesi ise Denklem 3.14’da ifade edilmiştir: k 2 (3.14) t  C  C sabit bir değerdir. Model sabitlerinin değerleri aşağıda ifade edilmiştir: C1 =1,44, C2 =1,92, C =0,09,  k =1,0,   =1,3. 17 Bu değerler pek çok deneysel çalışmanın tetkik edilmesiyle elde edilmiştir ve pek çok farklı uygulama alanı için önemli doğrulukta çalışmaktadır. RNG k-ε Modeli: RNG k-ε modeli, renormalizasyon grubu teorisi denilen istatistiksel bir teknik kullanılarak elde edildi (Yakhot ve Orszag 1992). Standart k-ε modeline benzer, fakat üstünde bir takım iyileştirmeler yapılmıştır. RNG k-ε modelinde, türbülans yayınım hızı (ε) denkleminde, hızla uzayan akışların doğruluğunu artıran ek bir terim vardır. RNG teorisi, türbülanslı Prandtl sayısı için analitik bir formül sağlarken, standart k-ε modeli kullanıcı tarafından belirlenen sabit değerleri kullanır. Standart k-ε modeli yüksek Reynolds sayılarında kullanılan bir modeldir fakat RNG teorisi düşük Reynolds sayılarındaki etkileri gösteren efektif viskozite için analitik olarak türetilmiş bir diferansiyel formül sunmaktadır.     u k  (3.15) k  kui   t   keff  G kG b Y M  Sk t xi x j   k x j          2 (3.16)   ui   eff  C1 G kC3G b C2  R  S  t x  i x j  x j  k k Denklem 3.15 ve 16’daki ifadeler Standard k-ε modeli ile aynı olmakla beraber  k k için ters etkili Prandtl sayısını,  ise ε için ters etkili Prandtl sayısını ifade eder. R ise kaynak terime bağlı bir ifadedir. Denklem 3.17, girdap modifikasyonu, türbülanslı ve buna bağlı olarak dönüşlü akışın, efektif Reynolds sayısına göre nasıl değiştiğine dair doğru bir açıklama elde etmek için entegre edilmiştir; bu ifade, modelin düşük Reynolds sayılarında ve duvara yakın olan akışlarda daha makul çözümlerin elde edilmesini sağlar. 18   2k   (3.17) d   1.72 d    3  1Cv  (3.18) eff    Cv 100 Türbülans viskozitesi Denklem 3.15’te gösterildiği gibi hesaplanmaktadır: k 2 (3.19)   C t   Burada C =0.0845’dir ve RNG teorisinden türetilmiştir. Daha önce Standard k-ε modelinde bu değerin deneysel verilerden yola çıkarak 0.09 şeklinde ifade edildiği belirtilmişti. İki değerin farklı yaklaşımlarla türetilmelerine rağmen birbirine yakın olması şaşırtıcıdır. Türbülanstaki dönme etkisi, dönen akışlarda doğruluğunun artması için RNG modeline dahil edilmiştir. Bu işlem türbülans viskozitesi denkleminin modifiye edilmesi ile sağlanmıştır:  k  (3.20) t  t0 f  s ,,     Burada, t0 , türbülans viskozitesi kullanılarak girdap modifikasyonu olmadan hesaplanan türbülans viskozitesinin değeridir.  , ANSYS Fluent içinde değerlendirilen karakteristik bir girdap sayısıdır ve  s , akışta dönüşlerin baskın olduğu veya sadece hafif dönüşlerin olduğu durumlara bağlı olarak farklı değerler alan girdap girdap sabitidir. Hafif dönüşlerin olduğu durumlarda  s = 0,07 olurken güçlü dönüşlerde bu ifade daha yüksek değerler alabilir. 19 Ters etkili Parandtl sayısıları olan  k ve  RNG teorisinden analitik olarak türetilmiş olan Denklem 3.21 ile elde edilir: 0,6321 0,3679 (3.21) a 1,3929 a  2,3929   mol a0 1,3929 a0  2,3929 eff    Burada a0 =1,0’dır. Yüksek Reynolds sayılarında mol  1 ,a  a 1,393olur.   k     eff  RNG k-ε modelinin model sabitleri ise: C1 =1,42, C2 =1,68’dir. Realizable k-ε Modeli: Realizable k-ε modelinin Standard k-ε modelinden iki önemli farkı vardır, bunlar:  Realizable k- 𝜀 türbülans modeli, türbülans viskozitesini farklı bir formülle ifade eder.  𝜀 transport denklemi, türbülans dalgalanmasının ortalama karekök metodu ile türetilmiştir. “Realizable” terimi bu türbülans modelinin Reynolds gerilmelerindeki matematiksel kısıtlamaları sağladığını ve türbülanslı akışın fiziği ile uyumlu olduğunu ifade eder. Sırasıyla k ve 𝜀 için transport denklemleri Denklem 3.22 ve 23’te ifade edilmiştir.     u  k  k  ku tj       G kG b YM  Sk (3.22) t x j x j   k  x j      u     2  (3.23)   u j   t     C1S  C k2 C1 C3Gb  S t x j x j     x j  k  v 20    (3.24) C1  max 0.43,    5   Sk /  , (3.25) Bu denklemlerde, G k , ortalama hız gradyanlarına bağlı olarak türbülans kinetik enerjisi üretimini ifade eder. G b , kaldırma etkileri nedeniyle oluşan türbülans kinetik enerjisi üretimidir. YM dalgalanan dilatasyonun sıkıştırılamaz türbülansın toplam yayınım oranına etkisini gösterir. C2 ve C1 sabittir.  k ve   , sırasıyla k ve ε için türbülans Prandtl sayılarıdır. Sk ve S , kullanıcı tanımlı kaynak terimleridir. C 2 1,0 Türbülans viskozitesi Standard ve RNG k-ε’da olduğu gibi hesaplanır: k 2 t  C   Fakat C diğer iki modelde olduğu gibi sabit değildir: 1 (3.26) C  * U k A0  A  S * Bu ifadede U : ~ ~ (3.27) U *  S ijSij ij ij Ayrıca denklemin diğer bileşenleri: ~ (3.28) ij  ij  2 ijkk ij  ij  ijk (3.29) k A0  4.04 AS  6 cos (3.30) 21 1 (3.31)   cos1  6W  3 S S (3.32) ij jkSki W  S 3 1  u j u  (3.33) j Sij     2  xi xi  Şeklinde ifade edilir. Bundan hareketle C ’nin ortalama uzamanın, dönme hızının, sistemin dönme hareketinin açısal hızının ve türbülans bileşenlerinin (k ve ε) bir fonksiyonu olduğu görülebilir. Realizable k-ε için model sabitleri: C1 1,44 , C2 1,9 ,  k 1.0 ,   1.2 ’dir. k-ε Modellerinde Türbülans Üretiminin Modellenmesi: Türbülans kinetik enerji üretimini temsil eden Gk terimi, Standard, RNG ve Realizable k-ε modelleri için aynı şekilde modellenmiştir. u (3.34) ' ' jGk  i j xi Gk 'nın Boussinesq hipotezi ile tutarlı bir şekilde değerlendirilmesi için: G   2 (3.35) k tS S , denklem 3.36’de tanımlanmıştır: S  2S (3.36) ijSij 22 k-ε Modellerinde Kaldırma Etkileri: Yerçekimi alanı ve sıcaklık gradyanı aynı anda mevcut olduğunda, ANSYS Fluent'teki k-ε modelleri kaldırma etkileri sebebiyle oluşan türbülans üretimini (Gb ) aşağıdaki şekilde ifade eder:  T (3.37) Gb   g t i Prt xi Prt enerji için türbülanslı Prandtl sayısıdır ve gi yerçekimi ivmesi bileşenini ve hangi doğrultuda olduğunu ifade eder. Standard ve Realizable k-ε modelleri için Prt =0.85 olarak tanımlanır. RNG k-ε modeli için ise Prt 1/ a şeklinde ifade edilir. Isıl genleşme katsayısı  : 1    (3.38)        T  p Şeklinde ifade edilir. İdeal gazlar için ise Gb ifadesi aşağıda gösterilmiştir:   (3.39) Gb  g t i  Prt xi Kaldırma etkilerinden etkilenme derecesi C3 sabiti ile belirlenir. C3 kullanılan yazılım içinde belirtilmemiş, bunun yerine hesaplanmıştır (Henkes ve ark. 1991): v (3.40) C 3  tanh u Yukarıdaki denklemde v , yerçekimi ivmesine paralel, u ise yerçekimi ivmesine dik akış hızının bileşenidir. Bu şekilde, C3 , ana akış yönünün yerçekimi yönüyle aynı hizada olduğu yüzer kayma tabakaları için 1 olacaktır. Yerçekimi vektörüne dik yüzer kayma tabakası için C3 sıfır olacaktır. k-ε Modellerinde Konvektif Isı ve Kütle Transferinin Modellemesi: 23 ANSYS Fluent'te türbülanslı ısı transferi, Reynolds’un türbülanslı momentum transferi benzetimi konseptinden modellenmiştir. “Modellenmiş” enerji denklemi Denklem 3.41’de ifade edilmiştir:     T  E  ui E  p  keff u i Tij    Sh (3.41) t xi x j  x eff j  E ifadesi toplam enerjiyi, keff ise efektif ısı iletimini ve Tij  deviatorik gerilme eff tensörünü ifade eder. Tij  yalnızca yoğunluk bazlı çözücülerde çözülür, basınç tabanlı eff çözücülerde ise çözülmez. Bu çalışmada çözücüler basınç tabanlıdır ve ifade yok sayılmıştır. Standard ve Realizable k-ε modelleri için, keff : c (3.42) pt keff  k  Prt RNG k-ε modeli için ise: k  ac  (3.43) eff p eff Şeklinde ifade edilir. 3.3.2 k-ω Türbülans Modelleri Çalışmada ω modellerinden Standard (Wilcox 1998), Baseline (BSL) (F. R. Menter 1994), Shear-Stress Transport (SST) (F. R. Menter 1994) k-ω modelleri incelenmiştir. Her üç model de k ve ω için benzer transport denklemlerine sahiptir. 24 Standard k-ω Modeli: ANSYS Fluent'taki Standard k-ω modeli, düşük Reynolds sayıları etkileri, sıkıştırılabilirlik ve kayma akışı yayınımı için bir takım modifikasyonlar içeren Wilcox’un k-ω modeline (Wilcox 1998) dayanmaktadır. Wilcox modelinin dezavantajı ise, kayma tabakasının dışında k ve ω değerlerinde görülen hassasiyettir. Bu hassasiyet daha sonra yapılan araştırmalar ve sonucundaki modifikasyonlarla bir miktar aşılmaya çalışılmıştır (Florian R. Menter 2009). Standart k- ω modeli, türbülans kinetik enerjisi (k) ve spesifik yayınım oranı (ω) için oluşturulmuş olan transport denklemlerine dayanan, deneysel bir modeldir. Bu model yıllar içinde pek çok defa modifiye edilmiştir. Bu tezde kullanılan yazılımda var olan ifadeler verilecektir. Türbülans kinetik enerjisi ve spesifik yayınım oranı sırasıyla Denklem 3.44 ve 3.45’teki transport denklemlerinden elde edilmektedir:     k  (3.44) k  kui   k  Gk Yk  Sk t xi x   j  x j        (3.45)   ui     G Y  S t xi x  j  x  j  Gk , ifadesi ortalama hız gradyanlarına bağlı türbülans kinetik enerjisi üretimini, G ise ω üretimini ifade eder. k ve  sırasıyla k ve ω için efektif difüziviteyi ifade eder. Yk ve Y ise türbülansa bağlı olarak sırasıyla k ve ω için yayınımı ifade eder. Sk ve S kaynak terimlerdir. k-ω modelinde efektif difüzivite k ve ω için Denklem 3.46 ve 3.47’de ifade edilmiştir:  (3.46) k    t  k  (3.47)     t  25 Bu denklemlerde  k ve  sırasıyla k ve ω için türbülanslı Prandtl sayılarıdır. Türbülans viskozitesi, t , k ve ω için Denklem 3.48’deki gibi ortak eşitlikle elde edilir: * pk (3.48) t  a  a* katsayısı, türbülans viskozitesini düşürerek düşük Reynolds sayısının düzeltilmesini sağlar. * * *  a0 Re / R  (3.49) a  a t k     1Ret / R k  k (3.50) Ret   R k  6 *  (3.51) a0  i 3 i  0.072 * k-ω modelinde yükek Reynolds sayılarının söz konusu olduğu durumlarda a = a* =1’dir. k için türbülans üretimi, Gk , Denklem 3.52’de olduğu gibi ifade edilir: u (3.52) j Gk   ' 'i j xi Gk 'nın Boussinesq hipotezi ile tutarlı bir şekilde değerlendirilmesi için ifade Denklem 3.52’teki yapıya dönüşür: Gk   S 2 t (3.53) Yukarıdaki ifadede S , k-ε modeli ile aynı şekilde tanımlanan ortalama gerilme tensörü modülüdür. ω için türbülans üretimi, G , Denklem 3.54’te ifade edilmiştir:  (3.54) G  a Gk k a  a0 Ret / R  (3.55) a    a*    1Ret / R  26 k Yukarıdaki denklemde R  2.95.a * şeklinde tanımlanır. Ret ise aynı şekilde Ret   hesaplanır. k’nın yayınımı Denklem 3.56’da gösterilmiştir: Yk   f *k (3.56)   1 (3.57)  k  0 f *  1 680 2 k  k  0 1 400 2 k 1 k  (3.58) k  3 xj x j  *   * i 1 *F M t  (3.59)  4   (3.60) 4 /15 Re* *  t / R i    1 Ret / R      * 1.5 R  8  *  0.09 ω’nın yayınımı Denklem 3.61’de olduğu gibi hesaplanır: Y   f 2 (3.61)   1 70 (3.62) f  180 (3.63) ij jk Ski   3  *  1  u u  (3.64)   i j ij  2     x x  j i  27   *  (3.65)   i 1 i  *F M t   i  F Mt  sıkıştırılabilme fonksiyonudur ve Denklem 3.66’da ifade edilmiştir:  0 M (3.66)   t  M t0 F M t   2 M t M 2 t0 M t  M t0 2 2k (3.67) M t  a2 M t0  0.25 (3.68) a  RT (3.69) Standard k-ω türbülans modelinin model sabitleri ise aşağıda gösterilmiştir: a* 1 , a  0.52 , 1 a0  , 9  *  0.09 , i  0.072 , R  8 , Rk  6 , R  2.95 ,  * 1.5 , M t0 =0.25,  k  2.0 ,   2.0 28 Baseline (BSL) k-ω Modeli: Wilcox’un Standard k-ω modeliyle ilgili temel sorun, serbest akış koşullarındaki güçlü hassasiyetidir. BSL modeli, Menter tarafından k-ω modelinin duvara yakın bölgedeki sağlam ve doğru formülasyonu ile k-ε modelinin cidardan uzak alandaki serbest akışın kararlı çözümünü harmanlamak için geliştirilmiştir (F. R. Menter 1994). Harmanlama fonksiyonu cidara yakın bölgede k-ω modelini aktifleştirmesi için bir, cidardan uzakta ise modifiye edilmiş k-ε modelini aktifleştirmesi için de sıfır olarak tasarlanmıştır. BSL k-ω türbülans modeli için transport denklemleri Denklem 3.70 ve 3.71’de ifade edilmiştir:     k  (3.70) k  kui   k  Gk Yk  Sk t x  i x j  x j        (3.71)   ui     G Y  DS t x i x j  x  j  Gk , Standard k-ω modelinde olduğu gibi türbülans kinetik enerjisi üretimini ifade eder ve aynı yolla hesaplanır, G ise ω üretimini ifade eder. k ve  sırasıyla k ve ω için efektif difüziviteyi ifade eder. Yk ve Y ise türbülansa bağlı olarak sırasıyla k ve ω için yayınımı ifade eder. Sk ve S kaynak terimlerdir. D  ise çapraz difüzyon terimidir. k-ω modelinde efektif difüzivite k ve ω için sırasıyla Denklem 3.72 ve 3.73’te gösterilmiştir:  (3.72)  t k     k  (3.73)  t     Bu denklemlerde  k ve  sırasıyla k ve ω için türbülanslı Prandtl sayılarıdır. Türbülans * pk viskozitesi ise Standard modelde olduğu gibi t  a formülü ile ifade edilir.  29 1 (3.74)  k  F1 / k ,1  1 F1  / k ,2 1 (3.75)   F1 / ,1  1 F1  / ,2 Harmanlama fonksiyonu, F1 , Denklem 3.76’da olduğu gibi ifade edilir: F  tanh  4  (3.76) 1 1       (3.77) k 500 4 k 1  min max  ,  ,    0.09y py 2   D 2  ,2  y    1 1 k   (3.78) D  max 2 ,10 10    ,2  x j x j  G ise Denklem 3.79’da gösterildiği gibi hesaplanır: aa* (3.79) G  Gk vt a  F (3.80)  1a,1  1F1 a,2  2 (3.81) i,1  a,1    * *  ,1   2 (3.82) i,2  a,2    * *  ,2  a terimi denklem 3.55’te ifade edilmiştir.  ’nın değeri 0.41’dir. k’nın yayınımı, Yk , ifadesi, türbülans kinetik enerjisinin yayınımını temsil eder ve Standard k-ω modelinde * olduğu gibi benzer bir şekilde tanımlanır. Aralarındaki fark ise f teriminden kaynaklanır. Bu terim Standard k-ω modelinde parçalı fonnksiyon olarak tanımlanırken BSL k-ω modelinde sabit değer alır ve değeri 1’e eşittir. Bu sebeple Yk aşağıdaki gibi ifade edilebilir: Yk   *k (3.83) 30 ω’nın yayınımı da Standard k-ω modeli ile benzerdir. İkisi arasındaki fark ise i ve f terimlerinden kaynaklanır. BSL modelinde i sabit değilken f ise sabittir. Bu durum Standard model de tam tersidir. Bu sebeple ifade Denklem 3.84’te olduğu gibi değişir: Y 2 (3.84)    i   F (3.85) i ii,1  1F1 i,2 F1  tanh  4  (3.86) 1 BSL modeli Standard k-ω ve Standard k-ε modellerinin harmanlanması ile oluşturulması çapraz difüzyon teriminin, Dw , ortaya çıkmasına sebep olduğu ifade edilmişti. Dw Denklem 3.87’de ifade edilmiştir: 1 k  (3.87) Dw  21 F1   ,2 x j x j BSL k-ε türbülans modelinin model sabitleri aşağıda ifade edilmiştir:  k ,1  2.0 ,  ,1  2.0 ,  k ,2 1.0 , ,2 1.168 , i,1  0.075 , i,2  0.0828 Shear-Stress Transport (SST) k-ω Modeli: SST k-ω modeli, BSL k-ω modelinin tüm geliştirmelerini içerir ve ek olarak, türbülans viskozitesinin hesabında türbülans kayma gerilmesini de hesaba katar (F. R. Menter 1994). BSL k-ω modeli, Wilcox’un Standard k-ω modeli ve k-ε modelinin avantajlarını birleştirir, ancak yine de düzgün yüzeylerden gelen akıştaki ayrılmanın başlangıcını ve miktarını düzgün bir şekilde tahmin edememektedir. Bunun sebebi, her iki modelin de türbülanslı kayma gerilmesinin taşınmasını hesaba katmamasıdır. Bu durum, türbülans 31 viskozitesinin olduğunun üzerinde değerlerde hesaplanmasına yol açar. Uygun taşıma davranışı, türbülans viskozitesinin denklemine eklenen sınırlayıcı bir ifade ile Denklem 3.88 de olduğu gibi elde edilebilir: k 1 (3.88)   t   1 SF  max 2 , a*   a1  Bu denklemde S uzama oranı, F harmanlama fonksiyonudur. 2 * (3.89) a*  a*  a 0 Ret / R k     1Ret / R k  F2  tanh  2  (3.90) 2  k 500  (3.91) 2  max 2 ,   0.09 y pyh2  Bu ifadede y bir sonraki yüzeye olan mesafedir. SST k-ω türbülans modelinin model sabitleri aşağıda ifade edilmiştir:  k ,1 1.176 ,  ,1  2.0 ,  k ,2 1.0 , ,2 1.168 a1  0.31, i,1  0.075 , i,2  0.0828 Diğer sabitler ise Standard k-ω modeli ile aynıdır. 32 3.4. Yakın Duvar Yaklaşımı Türbülanslı akışlar cidarlardan önemli bir biçimde etkilenir. Bunun sebebi cidarda meydana gelen kaymama koşulunun hız alanını etkilemesidir. Ayrıca türbülans da cidarların varlığından dolayı değişir ve bu değişimi çözmek zordur. Kaymama koşulundan dolayı, cidar bölgesinde, viskoz sönümleme teğet hız dalgalanmalarını azaltırken kinematik bloklama normal dalgalanmaları azaltır. Duvarın çok yakınında, viskoz sönümleme teğet hız dalgalanmalarını azaltırken kinematik bloklama normal dalgalanmaları azaltır. Bununla birlikte, duvara yakın bölgenin dış kısmına doğru, türbülans, ortalama hızdaki büyük gradyanlar nedeniyle türbülans kinetik enerjisinin üretilmesiyle hızla artar. Duvar bölgesinin üst katmanında, hız gradyanları büyüdüğü için türbülans kinetik enerjisi üretimi artar ve bu da sonuç olarak türbülansı arttırır. Yakın duvar modellemesi, duvarların ortalama vortisite ve türbülansın ana kaynağı olması nedeniyle, sayısal çözümlerin kalitesini büyük oranda etkiler. Bu sebeple doğru yakın duvar modellemesi yapmak daha makul çözümler elde etmemizi sağlar. Yapılan birçok deney sonucunda, duvarın yakınındaki akış bölgesinin üç katmana bölünebileceğini gösterilmiştir. “Viskoz alt tabaka” olarak adlandırılan en alttaki tabakada akış neredeyse laminerdir ve (moleküler) viskozite momentum ve ısı-kütle transferinde baskın bir rol oynar. Tamamen türbülanslı katman adı verilen dış katmanda türbülans önemli bir rol oynar. Son olarak, viskoz alt katman ile moleküler viskozite ve türbülansın etkilerinin eşit derecede önemli olduğu tam türbülanslı katman arasında bir ara bölge vardır. Aşağıdaki şekilde bu katmanlar gösterilmiştir (ANSYS Fluent Theory Guide 2016). 33 Şekil 3.5. Y+ değerine bağlı olarak sınır tabaka katmanları (ANSYS Fluent Theory Guide, 2016) Geleneksel olarak, duvara yakın bölgenin modellenmesinde iki yaklaşım vardır. Bir yaklaşımda, viskoziteden etkilenen iç bölge (viskoz alt tabaka ve tampon tabaka) çözülmez. Bunun yerine, duvar fonksiyonları olarak adlandırılan yarı ampirik formüller, viskoziteden etkilenen bölgeyi duvar ve tam türbülanslı bölge arasında bağ kurmak için kullanılır. Duvar fonksiyonlarının kullanılması, duvarın etkisini çözüme aktarmak için türbülans modellerini değiştirme ihtiyacını ortadan kaldırır. Başka bir çözüm yolu ise, modellerinin modifiye edilerek, viskoziteden etkilenen bölgede ve viskoz alt tabaka da dahil olmak üzere çözüm almasını sağlamaktır. Bunu başarabilmek için ise cidar etrafında çok sık bir ağ yapısı oluşturmak gerekmektedir. İki modelin farkı şematik olarak Şekil 3.6’da gösterilmiştir. 34 Şekil 3.6 Duvar fonksiyonu ve yakın duvar modeli çözümlerinin gösterimi (ANSYS Fluent Theory Guide 2016) Duvar fonksiyonu yaklaşımındaki en büyük dezavantaj bu yaklaşımın değerlerine karşı aşırı hassas olmasıdır. Genellikle duvar fonksiyonları farklı y+ aralıklarında etkili çalışır ve ağ yapısı oluşturulurken y+ değerlerine dikkat edilmezse elde edilen çözümler büyük ölçüde hatalı olacaktır. Çalışmada yakın duvar yaklaşımı için Güçlendirilmiş Duvar Yaklaşımı (Enhanced Wall Treatment – EWT) kullanılmıştır. EWT yaklaşımı, iki katmanlı model ile geliştirilmiş bir duvar fonksiyonunu birleştiren yakın duvar modellemesidir. Eğer cidardaki ağ yapısı yeterince sıksa iki denklemli modeller için makul duvar modeli EWT’dir (ANSYS Fluent 16.2.3 Theory Guide 2016). ANSYS Fluent’ın EWT yakın duvar modelinde, viskoziteden etkilenen duvara yakın bölge, viskoz alt tabakaya kadar tamamen çözülür. İki katmanlı yaklaşım, güçlendirilmiş duvar yaklaşımın en önemli özelliğidir ve duvarlara yakın hücrelerde hem ε hem de türbülans viskozitesini belirlemek için kullanılır. Bu yaklaşımda, viskoziteden etkilenen bölge ve tamamen türbülanslı bölge olmak üzere ikiye ayrılmıştır. Bu iki bölgenin sınırı, Denklem 3.92’te gösterilen, duvara olan uzaklığa dayalı, türbülanslı Reynolds sayısı ile belirlenir: 35  y k (3.92) Re y   Bu denklemde gösterilen y , hücre merkezinden duvara olan uzaklığı ifade eder. Hücre yakınında birden fazla yüzey varsa y değeri Denklem 3.93’te gösterildiği gibi hesaplanır:   (3.93) y  min r r w  r ww   Bu ifadede r , alan içerisindeki noktanın konum vektörüdür ve r w , duvar sınırının konum vektörüdür. w ilgili tüm duvar sınırlarını ifade eder. Bu yorum, y'nin birden fazla duvar içeren karmaşık hatasız olarak tanımlanmasını sağlar. * Türbülanslı bölgede ( Rey Rey ;Re * y  200 ) k-ε modelleri kullanılır. Yakın duvar * bölgesinde ise ( Rey  Rey ) tek denklemli Wolfstein modeli (Wolfshtein 1969) kullanılır. Türbülans viskozitesi, t , ise Denklem 3.94’te gösterilmiştir: t ,2layer  Cl k (3.94)  Bu ifadede yer alan uzunluk ölçeği, l , ise Denklem 3.95’te gösterildiği gibi hesaplanır: *  Rel  yC 1e y / A  (3.95)  1 İki katmanlı modelde elde edilen türbülans viskozitesi, Denklem 3.96’te gösterildiği gibi dış bölgedeki türbülanslı viskozite ile harmanlanır:  (3.96) t ,enh  t  (1 )t ,2layer Denklem 3.95’te gösterilen harmanlama fonksiyonu,  , 1   Re Re *  (3.97) y y   1 tanh   2   A  Bu denklemde A harmanlama fonsiyonun genişliğini ifade eder: 36 Re (3.98) y A  ar tanh 0.98 * Rey genellikle Rey ’nin %5’i ile %20’si arasında bir değer alır. Viskozitenin de etkilediği alanda  ve l değerleri sırasıyla Denklem 3.99 gösterildiği gibi hesaplanır: k 3/2 (3.99)   l l  yC*  Re / A1e y   (3.100)  l Bu ifadelerde yer alan sabitler ise aşağıda tanımlanmıştır: C*l C 3/4  A  70 A *  2Cl Ansys Fluent yazılımında türbülanslı ve laminer duvar yasası ifadelerini harmanlamak için Denklem 3.101 de gösterilen terim tanımlanmıştır (Kader 1981): u  eu  e1/u (3.101) lam turb Bu denklemde harmanlama fonksiyonu: 4    (3.102) a y    1 by du Harmanlama fonksiyonunda a değeri 0.01 b değeri ise 5’tir. türevi ise Denklem dy 3.103’te gösterilmiştir:  du du.4du (3.103)  e lam  e1/ turb dy dy dy 1 ay' for y  y  (3.104) S  s 1 ay   s for y  ys vw dp  dp (3.105) a    u 3dx  2 u*w  dx 37 2  *  (3.106)  t u   2cpTw Laminer duvar yasası aşağıdaki eşitlik ile Denklem 3.107’de ifade edilmiştir: du (3.107) lam 1 ay dy Bu ifadede yalnızca basınç gradyenlerinin etkisi vardır, ısı transferi ve sıkıştırılabilme özelliklerinin etkisi laminer duvar yasasında ihmal edilebilir (Kader 1981).    (3.108) u   lam  y 1 y   2   Termal duvar fonksiyonları u ’ya benzer bir şekilde elde edilir: Tw T c (3.109)  P pTT   eT  1/r lam  e Tturb q 4 (3.110) a Pr y     1bPr3 y (3.111) T   u  Pr u  * u2  lam  lam   2q       (3.112)    u* 2 Pr 2 2   Tturb  Prt uturb  P u    1uc  u*    2q   Prt   3.5. Isıl Konfor İklimlendirme sistemlerinin temel amacı, “termal çevre ile ilgili memnuniyeti ifade eden zihinsel durum” için uygun şartları sağlamaktır (ANSI/ASHRAE 2004). İlgili standartta da ifade edildiği üzere bu tanım, “zihinsel durum” veya “memnuniyet” yaklaşımlarını bir kenara bırakırsak, ısıl konforun değerlendirilmesinin, fiziksel, fizyolojik, psikolojik ve diğer süreçlerden etkilenen birçok girdiyi içeren bilimsel bir süreç ve aynı zamanda bir mühendislik problemi olduğunu doğru bir şekilde vurgulamaktadır. 38 İnsan zihni, vücudun yüzey sıcaklığından, nem hissinden, deri altı sıcaklıklarından ve vücut sıcaklıklarını düzenlemek için gereken işlemler vasıtası ile bir ısıl konfor ve rahatsızlık sonucuna varmaktadır (Gagge 1937). Isıl konforu etkileyen parametreler kişisel ve çevresel olarak sınıflandırılabilir. Çevresel parametreler olarak ortam sıcaklığı, ortam bağıl nemi, ortam hava hızı ve ortalama ışınım sıcaklığıdır. Kişisel parametreler ise kişinin metabolik aktivite düzeyi ve giyinme durumudur. Çevresel faktörler iç mekanlarda çoğu zaman iklimlendirme sistemleri ile kontrol edilebilmektedir fakat gerek kullanılan sistemin çalışma sistemi ve özellikleri gerek iç ortamın geometrisi nedeniyle incelenen hacmin tamamının konforlu hale getirilmesi her zaman kolay olmamaktadır. Çevresel ve kişisel faktörlerin yanı sıra ısıl konfora “diğer faktörler” de etki edebilir. Bu ikincil faktörler arasında iç mekanın bütünlüğü, görsel uyaranlar, yaş ve dış mekan iklim şartları sayılabilir. Rohles (1973) ve Rohles ve Nevins (1971) tarafından 1600 üniversite öğrencisi üzerinde yapılan çalışmalar, konfor seviyesi, sıcaklık, nem, cinsiyet ve maruz kalma süresi arasındaki ilişkiyi ortaya koymuştur. Bu çalışmalar için duyu ölçeğine ASHRAE ısıl duyarlılık ölçeği tanımlanmıştır: +3 Çok Sıcak +2 Sıcak +1 Ilık 0 Nötr -1 Serin -2 Soğuk -3 Çok Soğuk Isıl konforu sayısal olarak ifade edebilmek ve değerlendirmek için çeşitli çalışmalar yapılmıştır. PMV, insan vücudunun ısıl dengesinden yola çıkarak, büyük bir grup insanın 7 puanlık ısıl duyarlılık ölçeğindeki oylarının ortalama değerini öngören bir endekstir. Vücuttaki iç ısı üretimi ile çevreye olan ısı kaybı eşit olduğunda ısıl denge elde edilir. 39 Ilımlı bir ortamda, insan ısıl duyarlılık sistemi, sıcaklık dengesini korumak için cilt sıcaklığını ve ter salgısını otomatik olarak değiştirmeye çalışacaktır. PMV değerinin hesaplanışı Denklem 3.113’de gösterilmiştir (ISO 2005): PMV   0.036M0.303e  0.028 4 4 M W 3.96E8 f cl t 273  T  273   f h t T  (3.113)  cl r  cl c cl a 3.055.730.007M W  Pa 0.42 M W 58.15  0.0173M 5.87 pa  0.0014M 34Ta  1.0 0.2l (3.114) fcl  cl 1.05 0.1l cl ANSYS Fluent tek başına bu değeri hesaplama kabiliyetine sahip olmadığı için kullanıcı tanımlı fonksiyonlar (User Defined Functions - UDF) hazırlanmıştır. UDF’ler ilgili programın belirlediği limitler dahilinde programa müdahele edilebilmesine ve ek yetenekler kazandırılmasına olanak sağlar. UDF dosyaları C programlama dilinde hazırlanmaktadır. PPD, kendilerini çok sıcak ya da çok soğuk hisseden, ısıl durumundan memnun olmayan insanların yüzdesinin nicel bir tahminini yapan bir endekstir. PPD değerine göre ısıl olarak konforsuz hisseden insanlar, 7 noktalı ısıl duyarlılık ölçeğinde sıcak, çok sıcak, soğuk veya çok soğuk oy verecek olanlardır. PMV ile PPD arasındaki Şekil 3.7’de ifade edilmiştir. PPD Denklem 3.115’te gösterildiği gibi hesaplanır: 0.3353PMV 40.2170PMV 2  (3.115) PPD 10095e  40 Şekil 3.7 PMV ölçeğine göre PPD değerlerinin dağılımı (ASHRAE 2009) PPD sonuçlarının değerlendirilmesi için PMV değerlerinde olduğu gibi UDF hazırlanmış ve kullanılan yazılımla harmanlanmıştır. PMV ölçeğine göre PPD değerlerinin dağılımı. Gagge, enerji depolamanın ihmal edilmediği, geçici rejim için bir ısıl duyarlılık modeli geliştirmiştir (Gagge ve ark. 1971). Gagge’nin oluşturduğu model insan vücudunu iç içe iki silindir olarak kabul eder. İç silindir iç organlar, kaslar ve kemiklerden oluşan vücut içini veya diğer bir ifade ile koru simgelerken, dış silindir ise deri tabakasını simgeler. Gagge modelinde kor ve deri bölmesi arasında ısı ve kütle transferi olmasının yanında, deri tabakasından çevreye de duyulur ve gizli ısı transferi gerçekleşmektedir. Kor ve deri sıcaklığı yüzey boyunca homojendir. Metabolik ısı üretimi ve yapılan mekanik iş, solunum ile yapılan ısı alışverişi kor tabakası ile ilişkilidir, kor tabakasından deri tabakasına iletim ve kan akışı ile enerji alışverişi gerçekleşmektedir. Şekil 3.8’de insan ve çevresi arasında meydana gelen ısı transferinin iki katmanlı modelde gösterimi verilmiştir. 41 Şekil 3.8 İnsan ve çevresinin eş merkezli silindir modeli (ASHRAE 2009) Gagge modelinde kullanılan eşitlikler için Atmaca (2006) ve Arslanoğlu (2015)’nun çalışmalarından faydalanılmıştır. Kor ve deride birim zamanda depolanan enerjiyi ifade eden klasik enerji dengesi denklemleri Denklem 3.116 ve 117’de belirtilmiştir: Scr  M W  Cres E Q (3.116) res cr ,sk S Q  C RE  (3.117) sk cr ,sk sk Anlık olarak depolanan ısıl enerji, iç enerji artışına eşittir: S  (1).m.c dT / d  / A (3.118) cr p,b cr D S (3.119) sk .m.cp,b dTsk / d  / AD 42 Bu ifadede vücut özgül ısısı cp,b = 3490 j/kgK olarak tanımlanmıştır. A ise DuBois D denklemi ile edilen DuBois yüzey alanıdır (DuBois 1916) ve Denklem 3.120’de gösterilmiştir: A  0.202m0.425l0.725 (3.120) D Çizelge 3.1’de her bir vücut parçasının yüzey alanı ve ağırlığı verilmiştir. Bu çalışmada 16 farklı vücut parçasının ısıl davranışı incelenmiştir bu sebeple ısıl dengeyi ifade eden denklemler her bir vücut parçası için düzenlenmiştir: S (3.121) cr (i, )  M W  Cres i,  Eres i,  Qcr ,sk i,  Ssk i,  Qcr ,sk i,  C i,  R i,  E i,  (3.122) sk  Denklem 3.121 ve 3.122’de i vücut parçasını ifade ederken  ise zamanı ifade etmektedir. Bu modelde parçalar arasındaki iletimle olan ısı transferi ve kan dolaşımı ile olan ısı transferi küçük sıcaklık aralıkları arasında meydana geldikleri için ihmal edilmiştir. Denge denklemleri de son olarak Denklem 3.123 ve 3.124’teki forma ulaşmıştır: dT i,  (3.123) Scr i,   1  m  cr i cp,b   / A (i)  d  dTsk i,  (3.124) Ssk i,    mic p,b   / A(i)  d  43 Çizelge 3.1 16 farklı vücut parçası için alan ve ağırlık değerleri (Tanabe ve ark. 2002) Yüzey Alanı Ağırlığı i Vücut Parçası A(i) (m2) m(i) (kg) 1 Sol Ayak 0.056 0.48 2 Sağ Ayak 0.056 0.48 3 Sol diz altı 0.112 3.343 4 Sağ diz altı 0.112 3.343 5 Sol bacak 0.209 7.013 6 Sağ bacak 0.209 7.013 7 Pelvis 0.221 17.57 8 Baş 0.140 4.020 9 Sol el 0.050 0.335 10 Sağ el 0.050 0.335 11 Sol dirsek altı 0.063 1.373 12 Sağ dirsek altı 0.063 1.373 13 Sol kol 0.096 2.163 14 Sağ kol 0.096 2.163 15 Göğüs 0.175 12.40 16 Sırt 0.161 11.03 Toplam 1.87 74 Çizelge 3.2’de 16 farklı vücut parçası için kor sıcaklığı değerleri verilmiştir. Nötr deri sıcaklıkları ise daha önce yapılmış deneysel çalışmalardan elde edilmiştir (Arslanoğlu 2015). dT cr dT ve sk ifadelerinin ileri sonlu farklar açılımı ile ilk sıcaklıkları bilinen vücut d d parçalarının belirli zaman adımlarından yeni deri ve kor sıcaklıkları Denklem 3.125 ve 3.126’dan hesaplanabilmektedir: 44 Scr i,  Ai (3.125) Tcr i, 1 Tcr i,  1  micp,b Scr i,  Ai (3.126) Tsk i, 1 Tcr i,  1  micp,b Çizelge 3.2 16 farklı vücut parçası için kor sıcaklıkları (Tanabe ve ark. 2002) Nötr Kor Sıcaklığı i Vücut Parçası (°C) 1 Sol ayak 35,1 2 Sağ ayak 35,1 3 Sol diz altı 35,6 4 Sağ diz altı 35,6 5 Sol bacak 35,8 6 Sağ bacak 358 7 Pelvis 36,3 8 Baş 36,9 9 Sol el 35,4 10 Sağ el 35,4 11 Sol dirsek altı 35,5 12 Sağ dirsek altı 35,5 13 Sol kol 35,8 14 Sağ kol 35.8 15 Göğüs 36.5 16 Sırt 36.5 Oratalama 35.94 Giysiler, gerek kendi etkileri gerekse vücut ile aralarında oluşturduğu hava boşluğunun etkisi sebebi ile ısı ve kütle transferini azaltmaktadır. Her bir vücut parçası üzerindeki kıyafetin bu duruma farklı bir etkisi vardır bu etki giysi yalıtımı ile ifade edilir ve clo ile gösterilir. 1 clo = 0,155 m2-K/W’ a eşdeğerdir (Standard, 2009). Vücut parçaları üzerindeki giysilerin özellikleri Çizelge 3.3’te verilmiştir. 45 Çizelge 3.3 Vücut parçaları üzerindeki giysilerin özellikleri (Arslanoğlu 2015, Mccullough ve ark. 1989) Buharlaşma Giyimli vücut Kalınlık Isıl direnç Direnci Ref Parçaları (mm) (m2 °C/W) (m2kPa/W) Pelvis 1.270 0.036 0.0040 Sağ ve sol kol, göğüs, 1.270 0.036 0.0040 sırt Sağ ve sol diz altı, sağ 0.787 0.026 0.0041 ve sol bacak, pelvis Sağ ve sol ayak 1.270 0.036 0.0040 Cilt yüzeyinde gerçekleşen duyulur ısı kaybı giysilerden geçerek, insan ile bulunduğu ortam arasında iki yönde de gerçekleşebilir. Bu ısı transferinde cilt ile giysi arasında giysinin yalıtımının ısıl direnci ve giysi ile ortam arasındaki ısıl direncin önemli bir etkisi vardır. Deride gerçekleşen toplam duyulur ısı kaybı 16 farklı vücut parçası için Denklem 3.127’de gösterilmiştir (ASHRAE 2009): C i,  R i,   Tsk i, To i / Rt i (3.127) Yukarıdaki denklemde gösterilen To ifadesi operatif sıcaklıktır ve Ta ortam sıcaklığı ve Tr ışınım sıcaklığının, karşılıklı ısı geçiş katsayılarına göre ağırlıklı ortalaması olarak Denklem 3.128’de gösterilmiştir: T i  h T ih T  / h h  (3.128) o r r c a r c Bu çalışmada ortalama ışınım sıcaklığı, ortam sıcaklığına eşit alınmıştır. Isı taşınım katsayısı aktivite durumu (oturma, ayakta durma, yürüme vb.), ve bulunduğu ortamdaki havanın durgun veya hareketli olmasına bağlı olarak değişir. Kişinin oturma hali için 46 farklı hava hızlarında taşınım katsayısı Denklem 3.129’dan hesaplanabilmektedir (ASHRAE 2009): 8.3V 0.6 0.2 V  4.0 (3.129) h   c  3.1 0 V  0.2 Ofis ortamında oturarak çalışan biri için (dosyalama eylemi) M=70 W/m2 met alınmıştır (ASHRAE 2009). 16 farklı vücut parçası için toplam ısıl direnç ayrı ayrı olarak Denklem 3.120’da olduğu gibi ifade edilmelidir: r i,0 nl  r i,0 r i,0  (3.130) Rt i  Ra i  Ral i, j   R f i, j   r i,nl  j1  r i, j 1 r i, j  Bu denklemde Rf giysilerin iletim direncini, Ral giysi katmanları arasında kalan havanın iletim ve ışınım direncini, Ra dış ortam havasının taşınım ve ışınım direncini ifade eder. nl vücut parçası üzerindeki giysi sayısı, r ise yarıçaptır. 1 (3.131) R  a hc  hr Burada hr = 4.9 W/m 2K’dir. 1 (3.132) R  al hr  k / x Ral ifadesinde hr = 4,9 W/m 2K k = 24 mm, x =1.3 mm’dir (Mccullough ve ark. 1989). Deriden gizli ısı kaybı, Esk , iki şekilde gerçekleşir:  Vücuttan salgılanan terin buharlaşması, Ersw ,  Suyun deriden doğal difüzyonu, Edif , şeklindedir. Esk i,   Ersw  Edif i,  (3.133) 47 Duyulur ısı kaybının 16 farklı vücut parçası için zamana bağlı olarak Denklem 3.134’te gösterilmiştir: Esk i,   wi,   psk ,s i,  p  / R i (3.134) a  e,t Bu ifadede w boyutsuz deri ıslaklığı, p çevre havasının su buharı basıncı, psk ,s deri üzerindeki su buharının doyma basıncı (deri yüzey sıcaklığının doyma basıncı, Re,t ise hava ve giysi katmanının buharlaşma ile ısı geçişi direncidir. Boyutsuz deri ıslaklığı aynı zamanda gerçek ısı kaybının buharlaşma ile olabilecek en fazla ısı kaybına oranı olarak da ifade edilebilir ve Denklem 3.135’te gösterilmiştir: E (3.135) w  sk Emax Salgılanan terin buharlaşması ile oluşan ısı kaybı Ersw , salgılanan ter miktarı ile orantılıdır: E (3.136) rsw  mrswhfg Bu ifadede hfg suyun buharlaşma entalpisi, mrsw ise birim zamandaki terleme miktarıdır. E  1w 0.06E (3.137) dif rsw max Bu ifade ile birlikte boyutsuz deri ıslaklığı Denklem 3.138’de gösterilmiştir: E (3.138) w  wrsw wdif  wrsw  0.061w rsw rsw   0.060.94 Emax Bu denklem son olarak her bir vücut parçası için Denklem 3.139’da olduğu gibi ifade edilebilir: Ersw   (3.139) wi,   0.06 0.94 Emax (i, ) Gizli ısı kaybının hesaplanabilmesi içim 16 farklı vücut parçasının her birinin buharlaşma direncinin bulunması gerekmektedir: r(i,0) nl  r(i,0) r(i,0)  (3.140) Re,t i  Re,a i  Re,al (i, j)  Re, f (i, j)  r(i,nl) j1  r(i, j 1) r(i, j)  48 Bu denklemde Re, f , giysi kumaşlarının buharlaşma direnci, Re,al giysi katmanları arasında kalan durgun havanın buharlaşma direnci ve Re,a ise dış ortam havasının buharlaşma direncidir. 1 (3.141) R e,a  hcLR Bu ifadede LR Lewis oranıdır ve iç ortam koşullarında yaklaşık olarak 16.5 K/kPa olarak alınabilir (ASHRAE 2009). Re,al  a1exp(x / b) (3.142) Burada a = 0.0334 m2kPa/W, b =15 mm’dir (Mccullough ve ark. 1989). Solunum esnasında havaya taşınım ve buharlaşma sonucunda ısı transferi meydana gelir. Solunum ile olan ısı kaybı Denklem 3.143 yardımıyla hesaplanabilir (ASHRAE 2009): C (i, ) E (3.143) res res (i, )  0.0014M 34Ta  0.0173M (5.87 pa ) Vücudun metabolik aktiviteleri sonucunda daima ısı oluşur ve bu ısının, normal vücut sıcaklığını korumak için sürekli olarak dağıtılması ve ayarlanması gerekmektedir. Yetersiz ısı kaybı aşırı ısınmaya (hipertermi) yol açar ve aşırı ısı kaybı vücudun aşırı soğumasına (hipotermi) neden olur. 45°C'den yüksek veya 18°C'den düşük deri sıcaklığı ağrıya neden olur (Hardy ve ark. 1952). Sedenter aktivitelerde, fiziksel aktivitenin çok az olduğu durumlarda, konforlu cilt sıcaklığı 33 - 34°C arasındadır ve artan aktivite ile konfor azalır (Fanger 1967). Buna karşılık, iç sıcaklıklar aktivite ile artar. Beyindeki sıcaklık düzenleme merkezi, rahat bir şekilde dinlendiğinde yaklaşık 36.8°C'dir, yürürken yaklaşık 37.4°C'ye ve koşu sırasında 37.9°C'ye çıkar. 28°C'den düşük bir iç sıcaklık ciddi kalp ritim bozukluğuna ve ölüme neden olabilir. 46°C'den yüksek bir sıcaklık kalıcı beyin hasarına neden olabilir. Bu nedenle, vücut sıcaklığının dikkatli bir şekilde düzenlenmesi, konfor ve sağlık için çok önemlidir (ASHRAE 2009). İnsanın vücut sıcaklığından daha düşük bir sıcaklıktaki ortamda bulunması durumunda insandan çevreye olan ısı transferi artar. Bu durumda deri bölgesinde bulunan damarlar kısılarak bu bölgede kan akışı azalır (vazokonstriksiyon). Böylelikle ısı kaybı daha çok vücudun yüzey kısmında olur ve iç sıcaklık değerleri daha kolay korunur. Bu denetim 49 mekanizmasına bağlı olarak ısıl korunmanın sağlandığı bölge soğuğa karşı vazomotor denetim bölgesidir. Eğer sıcaklık bu bölgedeki korunma mekanizmasına rağmen düşmeye devam ederse kas gerilmesi, titreme ve kendiliğinden hareket etme gibi ek mekanizmalar devreye girmektedir. Vücut kor sıcaklığının 35 °C’nin altında düştüğü durumlarda verimlilikte önemli bir düşüş olur (Arslanoğlu 2015, ASHRAE 2009). Eğer cilt sıcaklığından daha yüksek bir ortamda bulunuluyorsa veya dış ortam ile yeterince ısı transferi sağlanamıyorsa, vücut damarları genişleyerek (vazodilasyon) deriye olan kan akışını artırır, böylece cilt sıcaklığı kor sıcaklığına yaklaşmış olur. Bu denetim mekanizmasının gerçekleştiği bölge, sıcağa karşı vazomotor denetim bölgesi olarak adlandırılır. Vücudun bu önlemine rağmen kor sıcaklığı 37°C değerinin üstüne çıkarsa, vücut ter salgılayarak çevreye olan ısı transferini arttırır. Deneysel verilerle ulaşılmış olan bazı bağıntılar, vazomotor denetimi, terleme ve titreme gibi sıcaklık denetim mekanizmalarının, belirli sıcaklık sinyalleri ile harekete geçmesine olanak tanır. Isıl duyarlılık modelinin çalışmasında çok önemli bir yeri olan bu bağıntılar Denklem 3.144, 3.145, 3.146, 3.147 ve 3.148’de gösterilmiştir:  0 Tcr ,m ( ) T (3.144) cr ,n WSIGcr ( )   Tcr ,m ( )Tcr ,n Tcr ,m ( ) Tcr ,n T (3.145) cr ,n Tcr ,m ( ) Tcr ,m ( ) Tcr ,n CSIG cr ( )    0 Tcr ,m ( ) Tcr ,n  0 Tsk ,m ( ) T (3.146) sk ,n WSIGsk ( )   Tsk ,m ( )Tsk ,n Tsk ,m ( ) Tsk ,n T T (3.147) sk ,n sk ,m ( ) Tsk ,m ( ) Tsk ,n CSIG ( )   sk  0 Tsk ,m ( ) Tsk ,n  0 T (3.148) b,m ( ) Tb,n WSIGb ( )   Tb,m ( )Tb,n Tb,m ( ) Tb,n Bu denklemlerde, WSIGcr kordan gelen sıcak sinyal, CSIGcr kordan gelen soğuk sinyal, WSIGsk deriden gelen sıcak sinyal, CSIGsk deriden gelen soğuk sinyal ve WSIGb vücuttan 50 gelen sıcak sinyaldir ve bu ifadeler daima pozitiftir. Her bir zaman adımı için ağırlıklı kor Tcr ,m , ve cilt sıcaklığı Tsk ,m hesaplanmıştır. Denetim sinyalleri ise bu sıcaklıklar vasıtası ile elde edilmiştir ve denklem 3.149 ve 3.150’de gösterilmiştir: 16 (3.149) Tcr i,  Ai T i1cr ,m ( )  AD 16 (3.150) Tsk i,  Ai Tsk ,m ( )  i1 AD Tb,m , insan vücudunun ortalama sıcaklığıdır ve Denklem 3.151’de ifade edilmiştir: T (3.151) b,m( ) ( )Tsk ,m( )1( )Tcr ,m( ) Daha önce kor veri sıcaklıklarındaki değişimin denetim mekanizmasını uyardığı ve kan akış debisinin değiştiği ifade edilmişti. Sıcaklık sinyallerine göre değişen kan akış debisi Denklem 3.152’de olduğu gibi hesaplanır: mbl ( )  6.3 200WSIGcr ( ) / 10.5CSIGsk ( ) / 3600 (3.152)  Kor tabakasından deri tabakasına iletim ve kan akışı ile olan toplam ısı transferi 16 farklı vücut parçası için Denklem 3.153’te ifade edilmiştir: q   (3.153) cr ,sk i,   K  cp,blmbl ( ) Tcr i, Tsk (i, ) Bu denklemde, K , kor ile deri arasındaki toplam ısı transfer katsayısıdır ve 5.28W/m2K’dir, cp,bl ise kanın özgül ısısıdır ve 4187 J/kgK’dir. Kan debisindeki değişmelerin, deri ve kor bölmelerinin göreceli kütlelerine olan etkisi Denklem 3.154’te olduğu gibi ifade edilir: a( )  0.04180.745 / 3600mbl ()0.585 (3.154) Deri ve kordan gelen sıcaklık sinyalleri ile ter üretimi meydana gelir ve ter miktarı Denklem 3.155 vasıtası ile hesaplanabilir: 51 mrsw( )  4,7.10 5WSIGb()expWSIGsk () /10,7 (3.155) Terlemeye bağlı olarak ortaya çıkan buharlaşma ile ısı kaybı ise Denklem 3.156’da gösterilmiştir: Ersw( )  m (3.156) rsw( )hfg Titreme vasıtası ile gerçekleşen metabolik enerji üretimi vazokonstriksiyona oranla daha etkili bir mekanizmadır. Titreme mekanizmasının gerçekleşmesi için eş zamanlı olarak hem kordan hem de deriden sinyal alınması gerekmektedir: M ( ) 19.4CSIG ( )CSIG ( ) (3.157) shiv sk cr Titremenin gerçekleştiği durumda ise toplam metabolik enerji aktivite sebebi ile oluşan enerji M act ve titreme ile oluşan enerjinin M shiv ’in toplamıdır: (3.158) M ( )  M act M shiv ( ) Bu çalışmada Gagge’in ısıl duyarlılık modeli HAD yazılımı ANSYS Fluent ile harmanlanmıştır ve eş zamanlı çözüm elde edilmiştir. Bu işlemin gerçekleşmesi için C programlama dilinde UDF dosyaları hazırlanmıştır (Bkz. Ek-1). Öncelikle incelenecek olan geometriler BDT (Bilgisayar Destekli Teknik Resim) yazılımı vasıtası ile HAD işlemine uygun olacak şekilde hazırlanmıştır. Ardından bu geometriler vasıtası ile ağ modelleri oluşturulmuştur. Ağ modelleri üzerinde ANSYS Fluent ile sınır koşulları oluşturulmuş ve çözüme hazır hale getirilmiştir. Çözümün başlaması ile beraber termal manken üzerine 16 farklı vücut parçası için başlangıç sınır şartı olarak UDF dosyasından okunan ilk sıcaklıklar atanmıştır. Sonraki zaman adımlarında ise manken etrafında HAD ile elde edilen sıcaklık ve hız değerleri UDF’te tanımlı Gagge modeline aktarılmış ve ardından yeni cilt sıcaklıkları türetilmiştir. Bu işlem her bir zaman adımında 16 farklı vücut parçası için tekrarlanmıştır. Çözüm algoritması Şekil 3.9‘de gösterilmiştir: 52 BDT MODELİ AĞ MODELİ Hayır AĞDAN BAĞIMSIZLIK ELDE EDİLDİ Mİ? Evet Başla 16 farklı vücut parçası için başlangıç deri sıcaklıkları girilir Akış ve sıcaklık için HAD çözümleri yapılır Isıl duyarlılık modeline hız ve sıcaklık değerleri aktarılır 16 farklı vücut parçası için yeni deri sıcaklıkları sonraki zaman adımına uygulanır Hayır Çözüm Yakınsadı mı? Evet Sonlandır Şekil 3.9 Kullanılan HAD yazılımının ve hazırlanan kodu eş zamanlı çalışma ilişkisi. 53 4. BULGULAR Bu bölümde öncelikle HAD çalışmalarında kullanılan yöntemleri doğrulamak için üç farklı ısı transfer mekanizması için teyit amaçlı doğrulama çalışmaları üzerinde araştırma yapılmıştır. Kullanılan yöntemler daha önce yapılmış deneysel çalışmalar ile doğrulandıktan sonra bir ofis odasında alttan ısıtma, radyatör ile ısıtma ve son olarak klima ile ısıtma durumları için hava akışı ve ısıl konfor sonuçları elde edilerek tartışılmıştır. 4.1. Validasyon Çalışmaları Enerji üretimi ve kullanımı, küresel ısınma, hava kirliliği, asit yağmuru, ozon tabakasının incelmesi, ormanların yok oluşu ve radyoaktif maddelerin emisyonu gibi bilinen çevresel problemlerle önemli ölçüde ilişkilidir. Sürdürülebilir bir kalkınma için enerji ile ilgili dolaylı veya dolaysız tüm problemler dikkate alınmalıdır (Dincer 2000). Bu sebeplerden dolayı sürdürülebilir enerji gelişim stratejilerinden biri açık bir şekilde enerji tasarrufudur (Lund 2007). Birleşmiş Milletler Çevre Programı (UNEP United Nations Environment Programme Sustainable Buildings and Climate Initiative 2017) verilerine göre, binalar küresel enerjinin yaklaşık % 40'ını, suyun % 25'ini, kaynakların ise % 40'ını kullanmaktadır. Bu yüzden inşaat sektörü hem enerji kullanımı (Jiang 2011) hem de sera gazı (GHG) emisyonları açısından (Wright ve ark. 2014) dikkat çekmektedir. Ayrıca, binalardaki enerji tüketimi, mevcut teknolojiler kullanılarak % 30 ila % 80 oranında azaltılabilir. Bina enerji verimliliğine yapılan yatırım, yatırım maliyetlerinde kısa sürede getiri sağlayıp, artan maliyetleri dengelemeye yardımcı olur ve doğrudan ve dolaylı tasarruf sağlar. Bina havalandırması binalarda enerji tüketiminde önemli bir yere sahiptir. Ayrıca, çalışanların üretkenliği ve verimliliği, bulundukları ortamun termal koşullarına bağlıdır (Akimoto ve ark. 2010, Shin ichi Tanabe ve ark. 2014). HVAC (Heating Ventilating and Air Conditioning - Isıtma Havalandırma ve iklimlendirme) sistemlerinin temel amacı bina sakinleri için temiz hava ve konforlu 54 koşullar sağlayarak iç mekan havasını düzenlemektir (Niu ve ark. 2007). Havalandırma için farklı hava dağıtım sistemleri mevcuttur ve bu sistemler Awbi (1998) tarafından özetlenip karşılaştırılmıştır. Düşük hızda hava girişlerine rağmen, oda içerindeki hava akışı genellikle mekanik üfleme sistemi nedeniyle türbülanslıdır. Havalandırma, doğal, karıştırmalı ve deplasmanlı havalandırma olarak sınıflandırılabilir. Karıştırmalı ve deplasmanlı havalandırma tipleri yaygın olarak kullanılan havalandırma türleridir. HAD simülasyonları ve ölçümler, deplasmanlı havalandırmanın karıştırmalı havalandırmadan enerji tüketimi açsısından daha verimli olduğunu göstermektedir (Awbi 1998). Geleneksel karıştırmalı havalandırma sistemleri düşük havalandırma verimliliğine ve daha az enerji tasarruflu olmasına rağmen, hala pazarın büyük bir bölümününe sahiptir. Son zamanlarda, yeni bir havalandırma türü olan çarpan hava jetli havalandırma sistemi deplasmanlı havalandırma ile benzer eğilimler göstermektedir, ancak bu yeni sistemin pratikte meydana gelen birçok farklı durum için havalandırma araştırmacıları ve tasarımcıları tarafından bir süre daha incelenmelidir (Karimipanah ve Awbi 2002). Hava akışı ve hız özellikleri, sıcaklık, nem ve kirletici madde konsantrasyonu gibi iç mekan hava parametrelerini düzenlemek ve kontrol etmek için çok önemlidir. Bu parametrelerin hava kalitesi, ısıl konfor, sağlık ve enerji tasarrufu üzerinde de önemli etkisi vardır (Moureh ve Flick 2005). İç ortam tasarımı, hava hızı / sıcaklık dağılımları, bağıl nem, kirletici madde konsantrasyonları ve türbülans seviyelerinin ayrıntılarını gerektiren zorlu bir mühendislik sürecidir. İç hava akışının oldukça karmaşıktır ve genellikle hem basınç gradyanı hem de ısınma etkisi ile kaldırma kuvveti tarafından meydana gelir edilir. Hava akışını kesin olarak tahmin etmenin zorluklarına rağmen, geçmişte hem deneysel ölçümler hem de HAD simülasyonları kullanılmıştır. Çoğu deneysel çalışmada, deney alanını çevreden izole etmek için yapay bir ortam oluşturmak gerekir ve tam ölçekli bir test odası kullanılır. Bu yöntem, kontrol edilebilir akış ve ısıl sınır koşullarına izin verirken, bir yandan da maliyeti yüksek olur ve test süreleri çok uzundur (Horikiri ve ark. 2011). Dolayısıyla, son yıllarda, HAD iç ortam havalandırması simülasyonunda önemli bir rol üstlenmiştir (P. V Nielsen 2004, Pitarma ve ark. 2004). Ancak, oda havalandırması karmaşık, türbülanslı 55 bir sistemdir ve özellikle türbülans modelinin kullanılması gereken simülasyon çalışmaları için güvenilir kıyaslama çalışmaları gerektirir. Bu çalışmada, k-ε ve k-ω türbülans modellerinin geçerliliği için zorlanmış taşınım, doğal taşınım ve iki durumun beraber meydana geldiği birleşik taşınımla ilgili literatürdeki bazı iyi bilinen deneyler kullanılmıştır. Bu deneylerin verileri kullanılarak farklı türbülans modelleri için kıyaslama çalışmaları yapılmıştır. Bu amaçla ANSYS-Fluent 18.0 yazılımı kullanılmıştır. Üç farklı ısı transfer mekanizması için çalışmalar yapılmıştır. 4.1.1. Zorlanmış Taşınım Zorlanmış taşınım, akışkan hareketinin pompa, fan vb. harici bir kaynak tarafından üretildiği bir ısı taşınım türüdür. Bu durum, Şekil 4.1'de gösterilmiş olan, iyi bilinen IEA Annex20 odası ile incelenmiştir. Bu çalışmanın sonuçları, Nielsen’in (1990) Laser Doppler Anemometre (LDA) ölçümlerinden alınan deney sonuçları ile karşılaştırılmıştır. Bütün duvarlar adyabatik olarak kabul edilmiştir. Duvarlar arasında sıcaklık farkı yoktur, bu sebeple kaldırma kuvveti etkisi ihmal edilmiştir. Hava, sol üst köşede bulunan üfleme kanalından 0.455 m/s hızında girer ve sağ alt köşede bulunan emme menfezinden dışarı çıkar. Üfleme menfezinin yüksekliğine (h) göre tanımlanan Reynolds sayısının değeri 5000’dir. Geometrinin ölçüleri Şekil 4.1'de verilmiştir. Ağ yapısı, çözümü etkileyen en önemli parametrelerden biridir. Ağ yapısının eleman sayısı ve kalitesi sadece doğru sonucu elde etmekle kalmaz, aynı zamanda çözüm süresini de etkiler. Duvarlarda daha yoğun bir ağ yapısı elde edilmiş ve dörtgen tip eleman seçilmiştir. Ağdan bağımsızlığı sağlamak için dört farklı eleman sayısında (4000, 18150, 28100, 43100) sonuçlar elde edilmiştir. Neredeyse tüm hız profilleri x=H kesitinde, Şekil 4.2 a, b, c, d,'den görüldüğü gibi birbirleriyle çakışmaktadır ve bunun sonucunda 28100 eleman sayısı seçilmiştir. 56 Şekil 4.1 IEA Annex 20 odasının ölçüleri ve sınır şartları (b) Deney Deney (a) (b) Deney Deney (c) (d) Şekil 4.2 Zorlanmış taşınım için a) x=H, b) x=2H, c) y=0,084m, d) y=2,916m kesitlerinde hız dağılımları Bu çalışmada elde edilen akım çizgileri deneysel veriler (Nielsen 1974) ve literatürdeki bazı diğer sonuçlar ile karşılaştırmalı olarak Şekil 4.3'te verilmiştir. İlk bakışta, k-ε modellerinin bu çalışmada k-ω modellerinden daha iyi performans gösterdiği 57 söylenebilir. BSL k-ω modelinin performansı diğer k-ω modellerine kıyasla nispeten iyi olarak değerlendirilebilir. Lattice Boltzmann Yöntemi (LBY) (Elhadidi ve Khalifa 2013a) ile elde edilmiş sonuçlarda, deneysel verilere benzer şekilde, sağ üst köşede ve sol alt köşede iki girdap vardır. SST k-ω modeli ve bazı diğer modellerde, Realizable k-ε (Dreau ve ark. 2012), Standard k-ω ve Rong ve Nielsen'in çalışmasında (Rong ve Nielsen 2008) sol alt köşedeki girdap deneysel çalışmaya kıyasla daha büyük bulunmuştur. 58 Şekil 4.3 Akım çizgileri (Bu çalışma: ANSYS-Fluent 18.0; Rong ve Nielsen (2008): ANSYS-CFX 11.0; Dreau ve ark. (2012): STAR-CCM+; Pulat and Ersan (2015): ANSYS-CFX 13.0 Non-isothermal) 59 X=H'de, hızın x bileşenin k-ε ve k-ω modelleri için boyutsuz hız profilleri gösterilmiş ve Nielsen'in (Nielsen, 1974) ölçümleriyle karşılaştırılarak, Şekil 4.5 a, b'de gösterilmiştir. k-ε modelleri için, odanın orta bölgesinde (0.3 ˂ y / H ˂ 0.8), en tutarlı modelden en tutarsıza doğru sıralama şu şekildedir: Standard k-ɛ, RNG k-ɛ ve Realizable k-ɛ. Voigt (Voigt 2005) tarafından belirtildiği gibi, bu çekirdek bölgede, hem RNG k-modeli hem de Realizable kε-modeli, açıklanması zor bir şekilde deneysel verilerden küçük de olsa sapma göstermektedir. Genel olarak, hız tahminlerinin bu bölgedeki ölçümlerle tutarlı olduğu söylenebilir. Ancak zeminde ve kısmen tavanda, tahminler ölçümlerden sapmaktadır. Tüm k-ɛ modelleri zeminde düşük hızları, 0 ˂ y / H ˂ 0.3, ve jet içindeki yüksek hızları, 0.9 ˂ y / H ˂ 1.0, çekirdek bölgeye benzer şekilde aynı sıra ile tahmin etmiştir. Bu bulgular aynı zamanda, Şekil 4.3'teki oda içi akış yapısı ile de tutarlıdır Şekil 4.3 d’de gösterilen (Pulat ve Ersan 2015) çalışma izotermal olmayan durum için geçerlidir. Fakat Zuo ve Chen (2010) tarafından belirtildiği gibi izotermal durumda elde edilen akım çizgileri ile izotermal olmayan durumda elde edilen akım çizgileri neredeyse aynıdır. k-ω modellerinde, k-ε modellerine benzer şekilde, hız tahminleri, SST k-ω modeli hariç odanın ortasında yapılan ölçümlerle tutarlıdır. Yine zeminde ve tavanda, tahminler ölçümlerden sapmaktadır. Hem k-ε hem de k-ω modelleri izotropik modellerdir ve bu sapmalar yakın duvar bölgesinde türbülans anizotropisine bağlanabilir. Zeminde, Standard ve SST k-ω modellerinde akışta dönme meydana gelmektedir. Akım çizgilerinin deneysel gösteriminde, sağ üst köşede ve sol alt köşede iki girdap bulunmaktadır. SST modelinin kullanıldığı durumda, Rong ve Nielsen'in (2008) çalışmasına benzer şekilde sol alt köşede oluşan girdap deneysel veriye göre oldukça büyüktür (Rong ve Nielsen 2008). Ayrıca, Şekil 4.5 c’de akışda meydana gelen dönmenin daha iyi anlaşılabilmesi için, SST k-ω modelinde elde edilen hız vektörleriyle birlikte x=H ve x=2H'deki hız profilleri beraber verilmiştir. Ölçümlerle en tutarlı model BSL k-ω modelidir. İkinci en iyi model ise std. k-ω’dır. x = 2H'de, Şekil 4.5 d ve e'de görüldüğü gibi, hemen hemen tüm modellerde benzer sonuçlar elde edilmiştir, ancak tavanda, k-ε tahminleri, k-ω tahminlerinden biraz daha iyidir. Hemen hemen tüm k-ω modelleri, zemindeki SST k-ω modeli sonuçları dışında, benzer sonuçları öngörmüştür. Ölçümlerle en tutarlı tahminler std. k-ε ile BSL k-ω modellerinden ede edilmiştir. 60 Y = 0.084 m'deki zemine yakın bölgede sol duvarın (x / H ˂ 0.3) ve çıkışın yanında ( x / H ˃ 2.8), hem k-ɛ hem de k-ω modellerinin hız tahminleri Şekil 4.4 a,b'den görüldüğü gibi deneysel verilerle tutarlıdır. Odanın ortasında (0.3 ˂ x / H ˂ 2.8) ölçülen hız değerlerinde ise deneysel verilerle farklılıklar vardır ve bu tür farklılıklar türbülanslı hızların dalgalanan kısmına bağlanabilir. En tutarlı modelden en tutarsız modele doğru olan sıralama Standard k-ɛ, RNG k-ɛ ve Realizable k-ɛ şeklindedir. Şekil 4.4 b’de görüldüğü üzere, k-ω modellerinde Standard ve SST k-ω modelleri, odanın ortasından sol duvara kadar olan kısımda (Standard k-ω için x/H ˂ ̴ 1.1 ve SST k-ω için x/H ˂ ̴ 1.6) odanın sol tarafındaki büyük girdap yapısı nedeniyle yüksek hızları daha yüksek tahmin etmişlerdir. Odanın tavana yakın olan kısmı y = 2.916 m'de ölçülen değerlerde, y = 0.084 m'de ölçülen değerlere kıyasla sayısal tahmin daha başarılıdır ve k-ε modellerinden elde edilen sonuçlar Şekil 4.4 c, d'de görüldüğü gibi k-ω modellerinden daha iyidir. Standard ve RNG k-ε modelleri Realizable k-ε modeline göre daha başarılıdır. Tüm k-ω modelleri yakın bir performansa sahiptir. Bazı modellerde odanın sağında oluşan girdap odanın üst kısmındaki akım çizgilerini etkilemediği için genel olarak tüm modeller birbirine yakın sonuçlar vermiştir. Sağ üst köşenin yakınında (2.7 ˂ x / H ˂ 3.0, Şekil 4.4’de görüldüğü gibi), tüm k-ω modelleri, boyutsuz hızı, k-ɛ modellerinden biraz daha iyi tahmin etmişlerdir. Bunun sebebi odanın sağ üst kenarında bulunan girdabın büyüklüğünü k-ɛ modellerinin, k-ω modellerinden daha küçük tahmin etmesidir (Bkz. Şekil 4.3 b, c, d, e, f, g). 61 Deney Deney (a) (b) Deney Deney (c) (d) Şekil 4.4 Boyutsuz hızın y=0.084 m’de x bileşeni for a) k-ε modelleri and b) k-ω modelleri ve Boyutsuz hızın y=2.916 m’de x bileşeni c) k-ε modelleri d) k-ω modelleri 62 Deney Deney (a) (b) (c) Deney Deney (d) (e) Şekil 4.5 Boyutsuz hızın x=H’da x bileşeni for a) k-ε modelleri, b) k-ω modelleri, c) Hız vektörleri ve x=H ve x=2H’da SST k-w modeliyle elde edilmiş U/U0 boyutsuz hız profili, x=2H’da boyutsuz hız profilleri a) k-ε modelleri and b) k-ω modelleri Son olarak, bu çalışmada x = H'deki boyutsuz hızın x-bileşeni tahminleri, Şekil 4.6'da görüldüğü gibi literatürdeki bazı SGS (Sub-Grid-Scale) BGS çalışmalarıyla karşılaştırılmıştır. Standard ve RNG k-ɛ modelleri ve Standard ve BSL k-ω modelleri, Emmerich ve McGrattan'ın (1998) BGS tahminlerine benzer bir performans 63 sergilemektedir (Emmerich ve McGrattan 1998). Tüm modellerin hız profilleri tavandan y / H = 0.8'e kadar olan kısımda birbirine oldukça yakındır. k-ɛ modellerinde, odanın orta yüksekliğinde (0.2 ˂ y / H ˂ 0.8), hız profilleri deneysel sonuçlarla ve Emmerich ve McGrattan'ın BGS modeliyle tutarlıdır (Emmerich ve McGrattan 1998). Ancak, tabanın yakınında (y / H ˂ 0.2), bu çalışmanın hız profilleri Emmerich ve McGrattan 'in BGS çalışmasının profilleri ile tutarlı olmasına rağmen, her iki profil de deneysel değerlerden daha küçüktür. Zemine yakın bölgede, Realizable k-ɛ modeli deneysel verilerle en uyumsuz sonuçları veren k-ɛ modeli olmuştur. k-ω modelleri için ise Şekil 4.6’da görüldüğü gibi, BSL k-ω modelinin hız profili Standard ve RNG k-ɛ modelinin hız profiline oldukça benzerdir. Duvara yakın bölgede Standard k-ω modeli, k-ɛ modellerinden farklı olarak Emmerich ve McGrattan’in (Emmerich ve McGrattan 1998) tahminleri ile değil, Davidson ve Nielsen'in (Davidson ve Nielsen 1996) BGS modeli ile tutarlıdır. SST k-ω modelinin hız profili ise, daha önce açıklandığı gibi diğer tüm profillerle karşılaştırıldığında y / H ˂ 0.8 için tamamen farklıdır. Davidson ve Nielsen'in (Davidson ve Nielsen 1996) BGS modelindeki tutarsızlık Smagorinsky modeline bağlanmıştır, çünkü sonuçlar Smagorinsky sabitine (Cs) çok bağımlıdır ve bu değer muhtemelen hem ağ yapısına hem de hem de akışa bağlıdır (Davidson ve Nielsen 1996). Bu nedenle, Cs değerinin model içinde sabit alınması yerine dinamik olarak hesaplanmasını önermişlerdir. (Bu çalışma) (Bu çalışma) Deney Deney Şekil 4.6 Nielsen (1990)’in, Davidson ve Nielsen (1996)’in SGS BGS tahminleri ve Emmerich ve McGrattan (1998)’ın zorlanmış taşınımda x=H’da boyutsuz hız profillerinin deneysel verilerle karşılaştırılması a) k-ε ve b) k-ω modelleri (Smagorinsky sabiti her BGS modeli için Cs = 0.14 alınmıştır.) 64 4.1.2. Doğal Taşınım Doğal taşınım, akışkanın sıcaklık değişimine bağlı olarak meydana gelen yoğunluk değişimlerinden dolayı kaldırma kuvveti etkilerinden kaynaklanan akışkan hareketinin üretildiği bir başka ısı taşınım türüdür. Doğal taşınım çalışması için literatürde sıklıkla kullanılan bir geometri ele alınmıştır (Betts ve Bokhari 2000). Doğal taşınım sebebi ile meydana gelen hava hareketi tipik bir pencere problemidir ve binalarda temel akış özelliği olarak sınıflandırılmıştır (Zuo ve Chen 2009). Sayısal sonuçları doğrulamak için (Betts ve Bokhari 2000) tarafından yapılan deneysel bir çalışma kullanılmıştır. Şekil 4.7'de doğal taşınım çalışması için kullanılan geometri ve sınır şartları gösterilmiştir. Soğuk ve sıcak duvar sıcaklığı sırasıyla 15.1°C ve 34.7°C'dir. Üst ve alt duvarlar adyabatiktir. Sonuç olarak, dikey duvarlar arasındaki sıcaklık farkı 19.6 K ve genişliğe (W) göre Rayleigh sayısı 0.86 x 106 'dır. Bu koşullar altında, geometrinin merkezindeki akış tam türbülanslı ve sıcaklık ile olan özellik değişimleri küçükse, yerçekimi kuvvetinin etkisini sonuçlarda görmek için Boussinesq yaklaşımı kullanılabilir. Sonuçlar, Şekil 4.7 ve 4.8'de görüldüğü gibi ağdan bağımsızlığı sağlamak için dört farklı eleman sayısında [20x40, 20x80, 20x160, 80x320] elde edilmiştir. Şekil 4.8 ve 4.9'de hemen hemen tüm sıcaklık ve dikey hız profilleri birbirleriyle çakışmaktadır. En tutarlı sıcaklık dağılımının y / H = 0.1'de görüldüğü tespit edilmiştir. Sonuç olarak 40x160 eleman sayısı çözüm için seçilmiştir. Sıcak duvar Soğuk duvar Şekil 4.7 Doğal taşınım için çalışılan geometrinin ölçüleri ve sınır şartları 65 Deney Deney Deney Şekil 4.8 Farklı konumlar için dört farklı ağ yapısında sıcaklık dağılımları y/H=0.1’de Farklı Ağ Yapılarında Boyutsuz y/H=0.5’de Farklı Ağ Yapılarında Boyutsuz y/H=0.9’de Farklı Ağ Yapılarında Hız Profilleri Hız Profilleri Boyutsuz Hız Profilleri Deney y Deney Deney y y Şekil 4.9 Farklı konumlar için dört farklı ağ yapısında hız dağılımları Betts ve Bokhari’nin (2000) deney düzeneği geniş en ve boy oranına sahiptir ve iki boyutlu olarak hız ve sıcaklık profilinin elde edilebilmesi için tasarlanmıştır. Yazarın bilgisi dahilinde, muhtemelen çok ince yapısından dolayı Betts ve Bokhari’nin geometrisinin verilerini içeren herhangi bir kıyaslama çalışmasına rastlanmamıştır. Dolayısıyla, sadece bu çalışma için akım çizgileri Şekil 4.10'da verilmiştir. Tüm modellerde akım çizgileri saat yönünün tersi yönde oluşmuştur ve tüm modellerde hemen 66 hemen akış yapısı aynı olduğundan, çeşitli y / H değerlerinde sıcaklık ve hız profilleri göz önüne alınarak karşılaştırmalar yapılmıştır. Şekil 4.10 Doğal taşınımda farklı türbülans modellerinde elde edilen akım çizgileri Seçilen y / H değerlerinde, Betts ve Bokhari’nin ölçümleriyle karşılaştırılmış olan sıcaklık ve dikey hız profilleri, Şekil 4.11’de verilmiştir. Realizable ve Standard k-ε modellerinin hız ve sıcak profilleri incelendiğinde sonuçların birbiri ile çok benzer hatta çakışık olduğu görülmektedir. Genel olarak, k-ε modelleri incelendiğinde hız ve sıcaklık profillerindeki farklılıklarının çok küçük olduğu söylenebilir. Sol ve sağ duvarların yakınında (x <10 mm ve x ˃ 70 mm), tüm modellerin sıcaklık profilleri, deneysel değerlerle oldukça 67 uyumludur. Akış yapısı tüm modeller için hemen hemen aynı olduğundan, duvarların yakınında sıcak ve soğuk sınır şartlarından dolayı tüm modellerde sıcaklık değerleri tek bir noktada birleşmiştir. Benzer bir yorum hız profilleri için de yapılabilir. Geometrinin orta bölgesinden (10 mm = xmin1) && (x <= xmax1)) /* line 49 */ { if ((y >= ymin1) && (y <= ymax1)) { //if ((z >= zmin1) && (z <= zmax1)) // { thermosensor_temperature1 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity1 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt1++; /* count number */ // } } } if ((x >= xmin2) && (x <= xmax2)) /* line 49 */ { if ((y >= ymin2) && (y <= ymax2)) { // if ((z >= zmin2) && (z <= zmax2)) // { 131 thermosensor_temperature2 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity2 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt2++; /* count number */ // } } } if ((x >= xmin3) && (x <= xmax3)) /* line 49 */ { if ((y >= ymin3) && (y <= ymax3)) { // if ((z >= zmin3) && (z <= zmax3)) // { thermosensor_temperature3 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity3 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt3++; /* count number */ // } } } if ((x >= xmin4) && (x <= xmax4)) /* line 49 */ { if ((y >= ymin4) && (y <= ymax4)) { // if ((z >= zmin4) && (z <= zmax4)) // { thermosensor_temperature4 += C_T(c, t); /* get thermocouple temperature */ 132 thermosensor_velocity4 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt4++; /* count number */ // } } } if ((x >= xmin5) && (x <= xmax5)) /* line 49 */ { if ((y >= ymin5) && (y <= ymax5)) { // if ((z >= zmin5) && (z <= zmax5)) // { thermosensor_temperature5 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity5 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt5++; /* count number */ // } } } if ((x >= xmin6) && (x <= xmax6)) /* line 49 */ { if ((y >= ymin6) && (y <= ymax6)) { // if ((z >= zmin6) && (z <= zmax6)) // { thermosensor_temperature6 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity6 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); 133 nt6++; /* count number */ // } } } if ((x >= xmin7) && (x <= xmax7)) /* line 49 */ { if ((y >= ymin7) && (y <= ymax7)) { // if ((z >= zmin7) && (z <= zmax7)) // { thermosensor_temperature7 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity7 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt7++; /* count number */ // } } } if ((x >= xmin8) && (x <= xmax8)) /* line 49 */ { if ((y >= ymin8) && (y <= ymax8)) { // if ((z >= zmin8) && (z <= zmax8)) // { thermosensor_temperature8 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity8 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt8++; /* count number */ // } 134 } } if ((x >= xmin9) && (x <= xmax9)) /* line 49 */ { if ((y >= ymin9) && (y <= ymax9)) { // if ((z >= zmin9) && (z <= zmax9)) // { thermosensor_temperature9 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity9 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt9++; /* count number */ // } } } if ((x >= xmin10) && (x <= xmax10)) /* line 49 */ { if ((y >= ymin10) && (y <= ymax10)) { // if ((z >= zmin10) && (z <= zmax10)) // { thermosensor_temperature10 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity10 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt10++; /* count number */ // } } } 135 if ((x >= xmin11) && (x <= xmax11)) /* line 49 */ { if ((y >= ymin11) && (y <= ymax11)) { // if ((z >= zmin11) && (z <= zmax11)) // { thermosensor_temperature11 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity11 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt11++; /* count number */ // } } } if ((x >= xmin12) && (x <= xmax12)) /* line 49 */ { if ((y >= ymin12) && (y <= ymax12)) { // if ((z >= zmin12) && (z <= zmax12)) // { thermosensor_temperature12 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity12 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt12++; /* count number */ // } } } if ((x >= xmin13) && (x <= xmax13)) /* line 49 */ 136 { if ((y >= ymin13) && (y <= ymax13)) { // if ((z >= zmin13) && (z <= zmax13)) // { thermosensor_temperature13 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity13 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt13++; /* count number */ // } } } if ((x >= xmin14) && (x <= xmax14)) /* line 49 */ { if ((y >= ymin14) && (y <= ymax14)) { // if ((z >= zmin14) && (z <= zmax14)) // { thermosensor_temperature14 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity14 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt14++; /* count number */ // } } } if ((x >= xmin15) && (x <= xmax15)) /* line 49 */ { if ((y >= ymin15) && (y <= ymax15)) 137 { // if ((z >= zmin15) && (z <= zmax15)) // { thermosensor_temperature15 += C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity15 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt15++; /* count number */ // } } } if ((x >= xmin16) && (x <= xmax16)) /* line 49 */ { if ((y >= ymin16) && (y <= ymax16)) { // if ((z >= zmin16) && (z <= zmax16)) // { thermosensor_temperature16 += (double)C_T(c, t); /* get thermocouple temperature */ thermosensor_velocity16 += (double)sqrt(C_U(c,t)*C_U(c,t) + C_V(c,t)*C_V(c,t) + C_W(c,t)*C_W(c,t)); nt16++; /* count number */ // } } } } end_c_loop(c, t) } 138 thermosensor_temperature1 = thermosensor_temperature1 / nt1; thermosensor_velocity1 = thermosensor_velocity1/ nt1; thermosensor_temperature2 = thermosensor_temperature2 / nt2; thermosensor_velocity2 = thermosensor_velocity2/ nt2; thermosensor_temperature3 = thermosensor_temperature3 / nt3; thermosensor_velocity3 = thermosensor_velocity3/ nt3; thermosensor_temperature4 = thermosensor_temperature4 / nt4; thermosensor_velocity4 = thermosensor_velocity4/ nt4; thermosensor_temperature5 = thermosensor_temperature5 / nt5; thermosensor_velocity5 = thermosensor_velocity5/ nt5; thermosensor_temperature6 = thermosensor_temperature6 / nt6; thermosensor_velocity6 = thermosensor_velocity6/ nt6; thermosensor_temperature7 = thermosensor_temperature7 / nt7; thermosensor_velocity7 = thermosensor_velocity7/ nt7; thermosensor_temperature8 = thermosensor_temperature8 / nt8; thermosensor_velocity8 = thermosensor_velocity8/ nt8; thermosensor_temperature9 = thermosensor_temperature9 / nt9; thermosensor_velocity9 = thermosensor_velocity9/ nt9; 139 thermosensor_temperature10 = thermosensor_temperature10 / nt10; thermosensor_velocity10 = thermosensor_velocity10/ nt10; thermosensor_temperature11 = thermosensor_temperature11 / nt11; thermosensor_velocity11 = thermosensor_velocity11/ nt11; thermosensor_temperature12 = thermosensor_temperature12 / nt12; thermosensor_velocity12 = thermosensor_velocity12/ nt12; thermosensor_temperature13 = thermosensor_temperature13 / nt13; thermosensor_velocity13 = thermosensor_velocity13/ nt13; thermosensor_temperature14 = thermosensor_temperature14 / nt14; thermosensor_velocity14 = thermosensor_velocity14/ nt14; thermosensor_temperature15 = thermosensor_temperature15 / nt15; thermosensor_velocity15 = thermosensor_velocity15/ nt15; thermosensor_temperature16 = thermosensor_temperature16 / nt16; thermosensor_velocity16 = thermosensor_velocity16/ nt16; *P_Tsensor1 = thermosensor_temperature1; *V_Tsensor1 = thermosensor_velocity1; 140 *P_Tsensor2 = thermosensor_temperature2; *V_Tsensor2 = thermosensor_velocity2; *P_Tsensor3 = thermosensor_temperature3; *V_Tsensor3 = thermosensor_velocity3; *P_Tsensor4 = thermosensor_temperature4; *V_Tsensor4 = thermosensor_velocity4; *P_Tsensor5 = thermosensor_temperature5; *V_Tsensor5 = thermosensor_velocity5; *P_Tsensor6 = thermosensor_temperature6; *V_Tsensor6 = thermosensor_velocity6; *P_Tsensor7 = thermosensor_temperature7; *V_Tsensor7 = thermosensor_velocity7; *P_Tsensor8 = thermosensor_temperature8; *V_Tsensor8 = thermosensor_velocity8; *P_Tsensor9 = thermosensor_temperature9; *V_Tsensor9 = thermosensor_velocity9; *P_Tsensor10 = thermosensor_temperature10; *V_Tsensor10 = thermosensor_velocity10; *P_Tsensor11 = thermosensor_temperature11; *V_Tsensor11 = thermosensor_velocity11; *P_Tsensor12 = thermosensor_temperature12; *V_Tsensor12 = thermosensor_velocity12; *P_Tsensor13 = thermosensor_temperature13; *V_Tsensor13 = thermosensor_velocity13; *P_Tsensor14 = thermosensor_temperature14; *V_Tsensor14 = thermosensor_velocity14; *P_Tsensor15 = thermosensor_temperature15; *V_Tsensor15 = thermosensor_velocity15; *P_Tsensor16 = thermosensor_temperature16; *V_Tsensor16 = thermosensor_velocity16; /* 141 printf("thermosensor_temperature1: %g\n", thermosensor_temperature1); printf("thermosensor_velocity1 : %f\n", thermosensor_velocity1); printf("thermosensor_temperature2: %g\n", thermosensor_temperature2); printf("thermosensor_velocity2: %f\n", thermosensor_velocity2); printf("thermosensor_temperature3: %g\n", thermosensor_temperature3); printf("thermosensor_velocity3: %f\n", thermosensor_velocity3); printf("thermosensor_temperature4: %g\n", thermosensor_temperature4); printf("thermosensor_velocity4: %f\n", thermosensor_velocity4); printf("thermosensor_temperature5: %g\n", thermosensor_temperature5); printf("thermosensor_velocity5: %f\n", thermosensor_velocity5); printf("thermosensor_temperature6: %g\n", thermosensor_temperature6); printf("thermosensor_velocity6: %f\n", thermosensor_velocity6); printf("thermosensor_temperature7: %g\n", thermosensor_temperature7); printf("thermosensor_velocity7: %f\n", thermosensor_velocity7); printf("thermosensor_temperature8: %g\n", thermosensor_temperature8); printf("thermosensor_velocity8: %f\n", thermosensor_velocity8); printf("thermosensor_temperature9: %g\n", thermosensor_temperature9); printf("thermosensor_velocity9: %f\n", thermosensor_velocity9); printf("thermosensor_temperature10: %g\n", thermosensor_temperature10); printf("thermosensor_velocity10: %f\n", thermosensor_velocity10); printf("thermosensor_temperature11: %g\n", thermosensor_temperature11); printf("thermosensor_velocity11: %f\n", thermosensor_velocity11); printf("thermosensor_temperature12: %g\n", thermosensor_temperature12); printf("thermosensor_velocity12: %f\n", thermosensor_velocity12); printf("thermosensor_temperature13: %g\n", thermosensor_temperature13); printf("thermosensor_velocity13: %f\n", thermosensor_velocity13); printf("thermosensor_temperature14: %g\n", thermosensor_temperature14); printf("thermosensor_velocity14: %f\n", thermosensor_velocity14); printf("thermosensor_temperature15: %g\n", thermosensor_temperature15); printf("thermosensor_velocity15: %f\n", thermosensor_velocity15); printf("thermosensor_temperature16: %g\n", thermosensor_temperature16); printf("thermosensor_velocity16: %f\n", thermosensor_velocity16); 142 printf("tsk16real: %g\n", TSK[15]); printf("nt16: %g\n", nt5); fflush(stdout); */ } ////////////////////////////////////------------------------////termosensor degerleri calisiyor, dogrulama yapıldı. x,y min max değerlerini dikkatle girmek gerekiyor//-------------------- ---/////////////////////////////////// DEFINE_EXECUTE_AT_END(bolgeler) { int sayac3; int sayac4; int sayac5; float solayak; float sagayak; float soldizalti; float sagdizalti; float solbacak; float sagbacak; float pelvis; float head; float solel; float sagel; float soldirsekalti; 143 float sagdirsekalti; float solkol; float sagkol; float gogus; float sirt; float* P_Temp1 = &temp1; float* P_Temp2 = &temp2; float* P_Temp3 = &temp3; float* P_Temp4 = &temp4; float* P_Temp5 = &temp5; float* P_Temp6 = &temp6; float* P_Temp7 = &temp7; float* P_Temp8 = &temp8; float* P_Temp9 = &temp9; float* P_Temp10 = &temp10; float* P_Temp11 = &temp11; float* P_Temp12 = &temp12; float* P_Temp13 = &temp13; float* P_Temp14 = &temp14; float* P_Temp15 = &temp15; float* P_temp16 = &temp16; face_t f; /* f is a face thread index */ /* if komutu zaman adımı içerisindeki iterasyonlarda */ float flow_time = CURRENT_TIME; float thermosensor_temperature1; float thermosensor_velocity1; float thermosensor_temperature2; float thermosensor_velocity2; 144 float thermosensor_temperature3; float thermosensor_velocity3; float thermosensor_temperature4; float thermosensor_velocity4; float thermosensor_temperature5; float thermosensor_velocity5; float thermosensor_temperature6; float thermosensor_velocity6; float thermosensor_temperature7; float thermosensor_velocity7; float thermosensor_temperature8; float thermosensor_velocity8; float thermosensor_temperature9; float thermosensor_velocity9; float thermosensor_temperature10; float thermosensor_velocity10; float thermosensor_temperature11; float thermosensor_velocity11; float thermosensor_temperature12; float thermosensor_velocity12; float thermosensor_temperature13; float thermosensor_velocity13; float thermosensor_temperature14; float thermosensor_velocity14; float thermosensor_temperature15; float thermosensor_velocity15; float thermosensor_temperature16; float thermosensor_velocity16; thermosensor_temperature1= Tsendor1; thermosensor_velocity1 = Vsendor1; 145 thermosensor_temperature2= Tsendor2; thermosensor_velocity2 = Vsendor2; thermosensor_temperature3= Tsendor3; thermosensor_velocity3 = Vsendor3; thermosensor_temperature4= Tsendor4; thermosensor_velocity4 = Vsendor4; thermosensor_temperature5= Tsendor5; thermosensor_velocity5 = Vsendor5; thermosensor_temperature6= Tsendor6; thermosensor_velocity6 = Vsendor6; thermosensor_temperature7= Tsendor7; thermosensor_velocity7 = Vsendor7; thermosensor_temperature8= Tsendor8; thermosensor_velocity8 = Vsendor8; thermosensor_temperature9= Tsendor9; thermosensor_velocity9 = Vsendor9; thermosensor_temperature10= Tsendor10; thermosensor_velocity10 = Vsendor10; thermosensor_temperature11= Tsendor11; thermosensor_velocity11 = Vsendor11; thermosensor_temperature12= Tsendor12; thermosensor_velocity12 = Vsendor12; thermosensor_temperature13= Tsendor13; thermosensor_velocity13 = Vsendor13; thermosensor_temperature14= Tsendor14; thermosensor_velocity14 = Vsendor14; thermosensor_temperature15= Tsendor15; thermosensor_velocity15 = Vsendor15; thermosensor_temperature16= Tsendor16; thermosensor_velocity16 = Vsendor16; 146 TA1=thermosensor_temperature1-273.15; TA2=thermosensor_temperature2-273.15; TA3=thermosensor_temperature3-273.15; TA4=thermosensor_temperature4-273.15; TA5=thermosensor_temperature5-273.15; TA6=thermosensor_temperature6-273.15; TA7=thermosensor_temperature7-273.15; TA8=thermosensor_temperature8-273.15; TA9=thermosensor_temperature9-273.15; TA10=thermosensor_temperature10-273.15; TA11=thermosensor_temperature11-273.15; TA12=thermosensor_temperature12-273.15; TA13=thermosensor_temperature13-273.15; TA14=thermosensor_temperature14-273.15; TA15=thermosensor_temperature15-273.15; TA16=thermosensor_temperature16-273.15; //TA=24; TR=24; W=0; RH=0.5; 147 Qisinim=80; k=24; a=0.0334; b=15; velo_aver1 = thermosensor_velocity1; velo_aver2 = thermosensor_velocity2; velo_aver3 = thermosensor_velocity3; velo_aver4 = thermosensor_velocity4; velo_aver5 = thermosensor_velocity5; velo_aver6 = thermosensor_velocity6; velo_aver7 = thermosensor_velocity7; velo_aver8 = thermosensor_velocity8; velo_aver9 = thermosensor_velocity9; velo_aver10 = thermosensor_velocity10; velo_aver11 = thermosensor_velocity11; velo_aver12 = thermosensor_velocity12; velo_aver13 = thermosensor_velocity13; velo_aver14 = thermosensor_velocity14; velo_aver15 = thermosensor_velocity15; velo_aver16 = thermosensor_velocity16; //velo_aver = 0.2; if (velo_aver1<= 0.15) { hc1=4 ; } if (0.15<=velo_aver1 && velo_aver1< 100.5) 148 { hc1=14.8*pow(velo_aver1,0.69); } //////////////////////////////////// if (velo_aver2<= 0.15) { hc2=4 ; } if (0.15<=velo_aver2 && velo_aver2< 100.5) { hc2=14.8*pow(velo_aver2,0.69); } //////////////////////////////////// if (velo_aver3<= 0.15) { hc3=4 ; } if (0.15<=velo_aver3 && velo_aver3< 100.5) { hc3=14.8*pow(velo_aver3,0.69); } //////////////////////////////////// if (velo_aver4<= 0.15) { 149 hc4=4 ; } if (0.15<=velo_aver4 && velo_aver4< 100.5) { hc4=14.8*pow(velo_aver4,0.69); } //////////////////////////////////// if (velo_aver5<= 0.15) { hc5=4 ; } if (0.15<=velo_aver5 && velo_aver5< 100.5) { hc5=14.8*pow(velo_aver5,0.69); } //////////////////////////////////// if (velo_aver6<= 0.15) { hc6=4 ; } if (0.15<=velo_aver6 && velo_aver6< 100.5) { hc6=14.8*pow(velo_aver6,0.69); } 150 //////////////////////////////////// if (velo_aver7<= 0.15) { hc7=4 ; } if (0.15<=velo_aver7 && velo_aver7< 100.5) { hc7=14.8*pow(velo_aver7,0.69); } //////////////////////////////////// if (velo_aver8<= 0.15) { hc8=4 ; } if (0.15<=velo_aver8 && velo_aver8< 100.5) { hc8=14.8*pow(velo_aver8,0.69); } //////////////////////////////////// if (velo_aver9<= 0.15) { hc9=4 ; } if (0.15<=velo_aver9 && velo_aver9< 100.5) 151 { hc9=14.8*pow(velo_aver9,0.69); } //////////////////////////////////// if (velo_aver10<= 0.15) { hc10=4 ; } if (0.15<=velo_aver10 && velo_aver10< 100.5) { hc10=14.8*pow(velo_aver10,0.69); } //////////////////////////////////// if (velo_aver11<= 0.15) { hc11=4 ; } if (0.15<=velo_aver11 && velo_aver11< 100.5) { hc11=14.8*pow(velo_aver11,0.69); } //////////////////////////////////// if (velo_aver12<= 0.15) { hc12=4 ; } 152 if (0.15<=velo_aver12 && velo_aver12< 100.5) { hc12=14.8*pow(velo_aver12,0.69); } //////////////////////////////////// if (velo_aver13<= 0.15) { hc13=4 ; } if (0.15<=velo_aver13 && velo_aver13< 100.5) { hc13=14.8*pow(velo_aver13,0.69); } //////////////////////////////////// if (velo_aver14<= 0.15) { hc14=4 ; } if (0.15<=velo_aver14 && velo_aver14< 100.5) { hc14=14.8*pow(velo_aver14,0.69); } //////////////////////////////////// 153 if (velo_aver15<= 0.15) { hc15=4 ; } if (0.15<=velo_aver15 && velo_aver15< 100.5) { hc15=14.8*pow(velo_aver15,0.69); } //////////////////////////////////// if (velo_aver16<= 0.15) { hc16=4 ; } if (0.15<=velo_aver16 && velo_aver16< 100.5) { hc16=14.8*pow(velo_aver16,0.69); } //////////////////////////////////// rt1=1/(hc1+4.9)+1/(4.9+(k/tx))+0.036; rt2=1/(hc2+4.9)+1/(4.9+k/tx)+0.036; rt3=1/(hc3+4.9)+1/(4.9+k/tx)+0.026; rt4=1/(hc4+4.9)+1/(4.9+k/tx)+0.026; rt5=1/(hc5+4.9)+1/(4.9+k/tx)+0.026; rt6=1/(hc6+4.9)+1/(4.9+k/tx)+0.026; rt7=1/(hc7+4.9)+1/(4.9+k/tx)+1/(4.9+k/tx)+0.036+0.026; rt8=1/(hc8+4.9); 154 rt9=1/(hc9+4.9); rt10=1/(hc10+4.9); rt11=1/(hc11+4.9); rt12=1/(hc12+4.9); rt13=1/(hc13+4.9)+1/(4.9+k/tx)+0.036; rt14=1/(hc14+4.9)+1/(4.9+k/tx)+0.036; rt15=1/(hc16+4.9)+1/(4.9+k/tx)+0.036; rt16=1/(hc16+4.9)+1/(4.9+k/tx)+0.036; ret1=1/(hc1*LR)+a*(1-exp(-tx/b))+0.0040; ret2=1/(hc2*LR)+a*(1-exp(-tx/b))+0.0040; ret3=1/(hc3*LR)+a*(1-exp(-tx/b))+0.0041; ret4=1/(hc4*LR)+a*(1-exp(-tx/b))+0.0041; ret5=1/(hc5*LR)+a*(1-exp(-tx/b))+0.0041; ret6=1/(hc6*LR)+a*(1-exp(-tx/b))+0.0041; ret7=1/(hc7*LR)+a*(1-exp(-tx/b))+a*(1-exp(-tx/b))+0.0040+0.0041; ret8=1/(hc8*LR); ret9=1/(hc9*LR); ret10=1/(hc10*LR); ret11=1/(hc11*LR); ret12=1/(hc12*LR); ret13=1/(hc13*LR)+a*(1-exp(-tx/b))+0.0040; ret14=1/(hc14*LR)+a*(1-exp(-tx/b))+0.0040; ret15=1/(hc15*LR)+a*(1-exp(-tx/b))+0.0040; ret16=1/(hc16*LR)+a*(1-exp(-tx/b))+0.0040; RT[0]=rt1; RT[1]=rt2; RT[2]=rt3; RT[3]=rt4; RT[4]=rt5; RT[5]=rt6; RT[6]=rt7; 155 RT[7]=rt8; RT[8]=rt9; RT[9]=rt10; RT[10]=rt11; RT[11]=rt12; RT[12]=rt13; RT[13]=rt14; RT[14]=rt15; RT[15]=rt16; RET[0]=ret1; RET[1]=ret2; RET[2]=ret3; RET[3]=ret4; RET[4]=ret5; RET[5]=ret6; RET[6]=ret7; RET[7]=ret8; RET[8]=ret9; RET[9]=ret10; RET[10]=ret11; RET[11]=ret12; RET[12]=ret13; RET[13]=ret14; RET[14]=ret15; RET[15]=ret16; if (TCRM<= 35.94) { WSIGCR=0.0 ; } 156 else { WSIGCR=TCRM-35.94 ; } /***********/ if (TCRM< 35.94) { CSIGCR=35.94-TCRM ; } else { CSIGCR=0 ; } /***********/ if (TSKM<= 32.35) { WSIGSK=0 ; } else { WSIGSK=TSKM-32.35 ; } /***********/ if (TSKM< 32.35) { CSIGSK=32.35-TSKM; 157 } else { CSIGSK=0 ; } /***********/ MBL=((6.3 + 200*WSIGCR)/(1 + 0.5*CSIGSK))/3600; if (MBL< 0.00014) { MBL=0.000141; } if (MBL>= 0.025) { MBL=0.0249; } if (TBM <= 35.4) { WSIGB=0; } else { WSIGB=TBM-35.4 ; } TO[0]= (hr*TR + hc1*TA1)/(hr + hc1); TO[1]= (hr*TR + hc2*TA2)/(hr + hc2); 158 TO[2]= (hr*TR + hc3*TA3)/(hr + hc3); TO[3]= (hr*TR + hc4*TA4)/(hr + hc4); TO[4]= (hr*TR + hc5*TA5)/(hr + hc5); TO[5]= (hr*TR + hc6*TA6)/(hr + hc6); TO[6]= (hr*TR + hc7*TA7)/(hr + hc7); TO[7]= (hr*TR + hc8*TA8)/(hr + hc8); TO[8]= (hr*TR + hc9*TA9)/(hr + hc9); TO[9]= (hr*TR + hc10*TA10)/(hr + hc10); TO[10]= (hr*TR + hc11*TA11)/(hr + hc11); TO[11]= (hr*TR + hc12*TA12)/(hr + hc12); TO[12]= (hr*TR + hc13*TA13)/(hr + hc13); TO[13]= (hr*TR + hc14*TA14)/(hr + hc14); TO[14]= (hr*TR + hc15*TA15)/(hr + hc15); TO[15]= (hr*TR + hc16*TA16)/(hr + hc16); PA[0]=RH*(exp(18.956-4030.18/(TA1 + 235)))/10; PA[1]=RH*(exp(18.956-4030.18/(TA2 + 235)))/10; PA[2]=RH*(exp(18.956-4030.18/(TA3 + 235)))/10; PA[3]=RH*(exp(18.956-4030.18/(TA4 + 235)))/10; PA[4]=RH*(exp(18.956-4030.18/(TA5 + 235)))/10; PA[5]=RH*(exp(18.956-4030.18/(TA6 + 235)))/10; PA[6]=RH*(exp(18.956-4030.18/(TA7 + 235)))/10; PA[7]=RH*(exp(18.956-4030.18/(TA8 + 235)))/10; PA[8]=RH*(exp(18.956-4030.18/(TA9 + 235)))/10; PA[9]=RH*(exp(18.956-4030.18/(TA10 + 235)))/10; PA[10]=RH*(exp(18.956-4030.18/(TA11 + 235)))/10; PA[11]=RH*(exp(18.956-4030.18/(TA12 + 235)))/10; PA[12]=RH*(exp(18.956-4030.18/(TA13 + 235)))/10; PA[13]=RH*(exp(18.956-4030.18/(TA14 + 235)))/10; PA[14]=RH*(exp(18.956-4030.18/(TA15 + 235)))/10; PA[15]=RH*(exp(18.956-4030.18/(TA16 + 235)))/10; 159 ALFA=(0.0418 + 0.745/(3600*MBL + 0.585)); MRSW = 0.000047*WSIGB*exp(WSIGSK/10.7); ERSW=MRSW*2430000; MSHIV=19.4*(CSIGSK*CSIGCR); M=MACT + MSHIV; NRES[0]=0.0014*M*(34-TA1)+0.0173*M*(5.87-PA[0]); NRES[1]=0.0014*M*(34-TA2)+0.0173*M*(5.87-PA[1]); NRES[2]=0.0014*M*(34-TA3)+0.0173*M*(5.87-PA[2]); NRES[3]=0.0014*M*(34-TA4)+0.0173*M*(5.87-PA[3]); NRES[4]=0.0014*M*(34-TA5)+0.0173*M*(5.87-PA[4]); NRES[5]=0.0014*M*(34-TA6)+0.0173*M*(5.87-PA[5]); NRES[6]=0.0014*M*(34-TA7)+0.0173*M*(5.87-PA[6]); NRES[7]=0.0014*M*(34-TA8)+0.0173*M*(5.87-PA[7]); NRES[8]=0.0014*M*(34-TA9)+0.0173*M*(5.87-PA[8]); NRES[9]=0.0014*M*(34-TA10)+0.0173*M*(5.87-PA[9]); NRES[10]=0.0014*M*(34-TA11)+0.0173*M*(5.87-PA[10]); NRES[11]=0.0014*M*(34-TA12)+0.0173*M*(5.87-PA[11]); NRES[12]=0.0014*M*(34-TA13)+0.0173*M*(5.87-PA[12]); NRES[13]=0.0014*M*(34-TA14)+0.0173*M*(5.87-PA[13]); NRES[14]=0.0014*M*(34-TA15)+0.0173*M*(5.87-PA[14]); NRES[15]=0.0014*M*(34-TA16)+0.0173*M*(5.87-PA[15]); syc=0; for (syc=0;syc<16;syc++) { QCRSK[syc] = (5.28 + 4187*MBL)*(TCR[syc] - TSK[syc]); Ne[syc]=(TSK[syc] - TO[syc])/RT[syc]; PSK[syc]=(exp(18.956-4030.18/(TSK[syc]+235)))/10; EMAX[syc]=(PSK[syc]-PA[syc])/RET[syc]; w[syc]=0.06 + (0.94*ERSW/EMAX[syc]); 160 if (w[syc] >= 1) { w[syc]=1; } WRSW[syc]=ERSW/EMAX[syc]; if (WRSW[syc]>=1) { WRSW[syc]=1; } EDIF[syc]=(1-WRSW[syc])*0.06*EMAX[syc]; ESK[syc]=ERSW + EDIF[syc]; SSK[syc]=QCRSK[syc]-(Ne[syc]+ESK[syc]); SCR[syc]=M-W-(NRES[syc])-QCRSK[syc]; DTCR[syc]=(SCR[syc]*ALAN[syc])/((1-ALFA)*AG[syc]*3490); DTSK[syc]=(SSK[syc]*ALAN[syc])/(ALFA*AG[syc]*3490); TSK[syc]=TSK[syc]+DTSK[syc]; TCR[syc]=TCR[syc]+DTCR[syc]; } TCRM=(TCR[0]*ALAN[0]+TCR[1]*ALAN[1]+TCR[2]*ALAN[2]+TCR[3]*ALAN[3 ]+TCR[4]*ALAN[4]+TCR[5]*ALAN[5]+TCR[6]*ALAN[6]+TCR[7]*ALAN[7]+TCR [8]*ALAN[8]+TCR[9]*ALAN[9]+TCR[10]*ALAN[10]+TCR[11]*ALAN[11]+TCR[1 2]*ALAN[12]+TCR[13]*ALAN[13]+TCR[14]*ALAN[14]+TCR[15]*ALAN[15])/1.8 7; TSKM=(TSK[0]*ALAN[0]+TSK[1]*ALAN[1]+TSK[2]*ALAN[2]+TSK[3]*ALAN[3] +TSK[4]*ALAN[4]+TSK[5]*ALAN[5]+TSK[6]*ALAN[6]+TSK[7]*ALAN[7]+TSK[ 8]*ALAN[8]+TSK[9]*ALAN[9]+TSK[10]*ALAN[10]+TSK[11]*ALAN[11]+TSK[12 ]*ALAN[12]+TSK[13]*ALAN[13]+TSK[14]*ALAN[14]+TSK[15]*ALAN[15])/1.87; TBM=ALFA*TSKM+(1-ALFA)*TCRM; solayak=TSK[0]; sagayak=TSK[1]; soldizalti=TSK[2]; 161 sagdizalti=TSK[3]; solbacak=TSK[4]; sagbacak=TSK[5]; pelvis=TSK[6]; head=TSK[7]; solel=TSK[8]; sagel=TSK[9]; soldirsekalti=TSK[10]; sagdirsekalti=TSK[11]; solkol=TSK[12]; sagkol=TSK[13]; gogus=TSK[14]; sirt=TSK[15]; *P_Temp1 = solayak; *P_Temp2 = sagayak; *P_Temp3 = soldizalti; *P_Temp4 = sagdizalti; *P_Temp5 = solbacak; *P_Temp6 = sagbacak; *P_Temp7 = pelvis; *P_Temp8 = head; *P_Temp9 = solel; *P_Temp10 = sagel; *P_Temp11 = soldirsekalti; *P_Temp12 = sagdirsekalti; *P_Temp13 = solkol; *P_Temp14 = sagkol; *P_Temp15 = gogus; *P_temp16=sirt; zaman=flow_time; 162 printf("TSK[0] = %f \n",TSK[0]); printf("TSK[1] = %f \n",TSK[1]); printf("TSK[2] = %f \n",TSK[2]); printf("TSK[3] = %f \n",TSK[3]); printf("TSK[4] = %f \n",TSK[4]); printf("TSK[5] = %f \n",TSK[5]); printf("TSK[6] = %f \n",TSK[6]); printf("TSK[7] = %f \n",TSK[7]); printf("TSK[8] = %f \n",TSK[8]); printf("TSK[9] = %f \n",TSK[9]); printf("TSK[10] = %f \n",TSK[10]); printf("TSK[11] = %f \n",TSK[11]); printf("TSK[12] = %f \n",TSK[12]); printf("TSK[13] = %f \n",TSK[13]); printf("TSK[14] = %f \n",TSK[14]); printf("TSK[15] = %f \n",TSK[15]); fflush(stdout); /* printf("tsk15: %g\n", sirt); printf("sirt: %g\n", TSK[15]); printf("DTSK[15]: %g\n", DTSK[syc]); printf("SCR[15]: %g\n", SCR[15]); printf("NRES[15]: %g\n", NRES[15]); printf("PA[15]: %g\n", TO[15]); printf("thermosensor_velocity16: %f\n", thermosensor_velocity16); printf("thermosensor_temperature16: %f\n", thermosensor_temperature16); 163 } ////////////////////////////////sol ayak//////////////////////////////////////////////////// DEFINE_PROFILE(SolAyak, t, i) { float solayak; face_t f; /* line 74 */ solayak= temp1; // printf("Sol Ayak Sicakligi: %f\n", solayak+273.15); begin_f_loop(f, t){ F_PROFILE(f, t, i) = (solayak + 273.15); } end_f_loop(f, t) printf("sol ayak: %f\n", solayak+273.15); fflush(stdout); } ////////////////////////////////sag ayak//////////////////////////////////////////////////// DEFINE_PROFILE(SagAyak, t, i) { float sagayak; 164 face_t f; /* line 74 */ sagayak= temp2; // printf("Sag Ayak Sicakligi = %f\n", sagayak+273.15); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sagayak+273.15; } end_f_loop(f, t) //printf("deri sicakligi: %f\n", sagayak+273.15); fflush(stdout); } ////////////////////////////////Sol Diz Alti//////////////////////////////////////////////////// DEFINE_PROFILE(SolDizAlti, t, i) { float soldizalti; face_t f; soldizalti= temp3; // printf("Sol Diz Alti Sicakligi = %f\n", soldizalti+273.15); begin_f_loop(f, t){ F_PROFILE(f, t, i) = soldizalti+273.15; } end_f_loop(f, t) // printf("deri sicakligi: %f\n", soldizalti+273.15); 165 fflush(stdout); } ////////////////////////////////Sag Diz Alti//////////////////////////////////////////////////// DEFINE_PROFILE(SagDizAlti, t, i) { float sagdizalti; face_t f; sagdizalti= temp4; // printf("Sag Diz Alti Sicakligi = %f\n", sagdizalti+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sagdizalti+273.15; } end_f_loop(f, t) fflush(stdout); } ////////////////////////////////Sol Bacak//////////////////////////////////////////////////// DEFINE_PROFILE(SolBacak, t, i) { 166 float solbacak; face_t f; solbacak= temp5; // printf("Sol Bacak Sıcakligi = %f\n", solbacak+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = solbacak+273.15; } end_f_loop(f, t) } ////////////////////////////////Sag Bacak//////////////////////////////////////////////////// DEFINE_PROFILE(SagBacak, t, i) { float sagbacak; face_t f; sagbacak= temp6; // printf("Sag Bacak Sicakligi = %f\n", sagbacak+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sagbacak+273.15; } end_f_loop(f, t) } ////////////////////////////////Pelvis//////////////////////////////////////////////////// 167 DEFINE_PROFILE(PelVis, t, i) { float pelvis; face_t f; pelvis= temp7; // printf("Pelvis Sicakligi = %f\n", pelvis+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = pelvis+273.15; } end_f_loop(f, t) } ////////////////////////////////Head//////////////////////////////////////////////////// DEFINE_PROFILE(Head, t, i) { float head; face_t f; head= temp8; // printf("Baş Sicakligi = %f\n", head+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = head+273.15; } 168 end_f_loop(f, t) } ////////////////////////////////Sol El//////////////////////////////////////////////////// DEFINE_PROFILE(SolEl, t, i) { float solel; face_t f; solel= temp9; // printf("Sol El Sicakligi = %f\n", solel+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = solel+273.15; } end_f_loop(f, t) } ////////////////////////////////Sag El//////////////////////////////////////////////////// DEFINE_PROFILE(SagEl, t, i) { float sagel; face_t f; 169 sagel= temp10; // printf("Sag El Sicakligi = %f\n", sagel+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sagel+273.15; } end_f_loop(f, t) } ////////////////////////////////Sol Dirsek Alti//////////////////////////////////////////////////// DEFINE_PROFILE(SolDirsekAlti, t, i) { float soldirsekalti; face_t f; soldirsekalti= temp11; //printf("Sol Dirsek Alti Sicakligi = %f\n", soldirsekalti+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = soldirsekalti+273.15; } end_f_loop(f, t) } ////////////////////////////////Sag Dirsek Alti//////////////////////////////////////////////////// DEFINE_PROFILE(SagDirsekAlti, t, i) { float sagdirsekalti; 170 face_t f; sagdirsekalti= temp12; // printf("Sag Dirsek Alti Sicakligi = %f\n", sagdirsekalti+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sagdirsekalti+273.15; } end_f_loop(f, t) } ////////////////////////////////Sol Kol//////////////////////////////////////////////////// DEFINE_PROFILE(SolKol, t, i) { float solkol; face_t f; solkol= temp13; //printf("Sol Kol Sicakligi = %f\n", solkol+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = solkol+273.15; } end_f_loop(f, t) } ////////////////////////////////Sag Kol//////////////////////////////////////////////////// 171 DEFINE_PROFILE(SagKol, t, i) { float sagkol; face_t f; sagkol= temp14; //printf("Sag Kol Sicakligi = %f\n", sagkol+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sagkol+273.15; } end_f_loop(f, t) } ////////////////////////////////Gogus//////////////////////////////////////////////////// DEFINE_PROFILE(Gogus, t, i) { float gogus face_t f; gogus= temp15; // printf("Gogus Sıcakligi = %f\n", gogus+273.15); fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = gogus+273.15; } end_f_loop(f, t) } ////////////////////////////////Sirt//////////////////////////////////////////////////// DEFINE_PROFILE(Sirt, t, i) { float sirt; 172 face_t f; sirt= temp16; // printf("Sirt Sicakligi = %f\n", sirt+273.15); // fflush(stdout); begin_f_loop(f, t){ F_PROFILE(f, t, i) = sirt+273.15; } end_f_loop(f, t) } real f_cl=1+0.2*I_cl; real R_cl=0.155*I_cl; DEFINE_PROFILE(inlet_x_velocity, ft, var) { real flow_time; int a; int b; face_t f; /* Face index has its own type */ flow_time = CURRENT_TIME; /* Special Fluent macro */ begin_f_loop(f,ft) /* Special Fluent face loop macro */ { F_PROFILE(f,ft,var) =4.5*sin(1.79*flow_time); } end_f_loop(f,ft) } DEFINE_EXECUTE_AT_END(PMV_PPD) { Domain *domain; Thread *c_thread; cell_t c; real temp; real velo_aver; 173 real h_c; real PD; real PMV; real v_o; real t_cl; domain=Get_Domain(1); real Tu; real ke; /*loops over all cell threads in domain*/ thread_loop_c(c_thread, domain) { /* loops over cells in a cell thread */ begin_c_loop(c, c_thread) { temp = C_T(c, c_thread); velo_aver = sqrt(C_U(c, c_thread)*C_U(c, c_thread)+C_V(c, c_thread)*C_V(c, c_thread)+C_W(c, c_thread)*C_W(c, c_thread)); h_c=12.1*pow(velo_aver, 0.5); t_cl=35.7-0.0275*(Metabolic-WW)-R_cl*((Metabolic-WW)- 3.05*(5.73-0.007*(Metabolic-WW)-Paa) -0.42*((Metabolic-WW)-58.15)-0.0173*Metabolic*(5.87- Paa)-0.0014*Metabolic*(34-(temp-273))); /* PMV is C_UDMI(c,c_thread,0)*/ PMV=(0.303*pow(Euler_num, - 0.036*Metabolic)+0.028)*((Metabolic-WW)-(3.96e-8)*f_cl*(pow(t_cl+273, 4)- pow(t_rad+273, 4))- f_cl*h_c*(t_cl-(temp - 273))-3.05*(5.73- 0.007*(Metabolic-WW)-Paa)-0.42*((Metabolic-WW)-58.15)- 0.0173*Metabolic*(5.87-Paa)-0.0014*Metabolic*(34- (temp - 273))); 174 if (PMV<-3) PMV=-3; if(PMV>3) PMV=3; C_UDMI(c,c_thread,0)=PMV; /* PPD is C_UDMI(c,c_thread,1)*/ C_UDMI(c,c_thread,1)=100-95*pow(Euler_num, - (0.3353*pow(C_UDMI(c,c_thread,0), 4)+0.2179*pow(C_UDMI(c,c_thread,0), 2))); ke=C_K(c, c_thread); Tu=(100/velo_aver)*sqrt((2/3)*ke); if (velo_aver<=0.05) v_o=0.05; else v_o=velo_aver; /* PD is C_UDMI(c,c_thread,2)*/ PD= (34-(temp-273))*pow(v_o-0.05,0.062)*(0.37*v_o*Tu + 3.14); if (PD<0) PD=0; C_UDMI(c,c_thread,2)= PD; } end_c_loop(c, c_thread) } } 175 ÖZGEÇMİŞ Adı Soyadı : Bahadır Erman YÜCE Doğum Yeri ve Tarihi : Erzurum / 05.11.1991 Yabancı Dil : İngilizce Eğitim Durumu Lise : Bursa Atatürk Lisesi (Y.D.A.) Lisans : Selçuk Üniversitesi - Makine Mühendisliği Yüksek Lisans : Uludağ Üniversitesi – Makine Mühendisliği (Enerji) Çalıştığı Kurum/Kurumlar : Bitlis Eren Üniversitesi (2013 - ) Bursa Uludağ Üniversitesi (2013 - ) İletişim (e-posta) : bahadirermanyuce@gmail.com Yayınları : YÜCE, B.E., PULAT, E., “Alttan Isıtma Sisteminin Kullanıldığı Bir Ofis Odasında Isıl Konfor Ve İç Hava Kalitesinin Sayısal Olarak İncelenmesi”, 14. Ulusal Tesisat Mühendisliği Kongresi, 17-20 Nisan, sy.1094-1102, İzmir, 2019. YÜCE, B.E., PULAT, E., “Numerıcal Predıctıon Of Indoor Air Distribution And Thermal Comfort For an Office Room”, Roomvent & Ventilation 2018, 3-5 Haziran, sy.1049- 1054, Espoo, Finland, 2018. YÜCE, B.E., PULAT, E., “Isıtılan Bir Ofis Odasında Hava Akışı Ve Konforun Sayısal Simülasyonu”, IV. Uluslararası Katılımlı Anadolu Enerji Sempozyumu, 18-20 Nisan, sy.1012-1121, Edirne, 2018. YÜCE, B.E., PULAT, E., “Forced, Natural and Mixed Convection Benchmark Studies for Indoor Thermal Environments”, International Communications in Heat and Mass Transfer, Volume 92, sy.1-14, Mart 2018. 176 YÜCE, B.E., PULAT, E., “Bir Ofis Odasındaki Termal Akışın Kış Şartlarında Sayısal Olarak İncelenmesi”, 13. Ulusal Tesisat Mühendisliği Kongresi, 19-22 Nisan, sy.1231- 1240, İzmir, 2017. YÜCE, B.E., PULAT, E., “Numerical Simulation of Thermal Flow and Comfort in an Office Room”, 12th International Conference on Heat Transfer, Fluid Mechanics and Thermodynamics, 11-13 July, sy.468-472, Costa del Sol, Spain, 2017. TAN, F., CANBOLAT, A.S., YÜCE, B.E., TÜRKAN, B., “Experimental and Numerical Investigation of Strawberry Drying”, Uludağ University Journal of The Faculty of Engineering, Vol. 21, No.1, sy.205-218, 2017. YÜCE, B.E., PULAT, E., “Oda Havalandırmasında Isıl Konforun Sayısal Simülasyonu”, 12. Ulusal Tesisat Mühendisliği Kongresi, 8-11 Nisan, sy.2333-2342, İzmir, 2015. 177