Transient Heat Conduction in Concrete During Fire – Ghina Khairani – Tugas Besar Komputasi Teknik

A. Project Title

Transient Heat Conduction in Concrete During Fire

B. Author Complete Name

Ghina Khairani ; Student ID – 2306215274 ; Subject – Computational Engineering

C. Affiliation

Departemen Teknik Mesin, Universitas Indonesia, Depok, Jawa Barat.

D. Abstract

Kebakaran dalam kompartemen berdinding beton memanaskan satu sisi dinding, sementara sisi lainnya berbatasan dengan ruang yang masih dihuni. Laporan ini menghitung seberapa cepat panas merambat melintasi dinding beton setebal 15 cm yang salah satu sisinya terpapar kebakaran standar ISO 834, sebagai salah satu acuan waktu penyelamatan dan evakuasi. Masalah dimodelkan sebagai konduksi panas transien satu dimensi, dengan kondisi batas konveksi dan radiasi non-linear pada sisi api. Persamaan diselesaikan secara numerik dengan metode beda hingga: skema FTCS (eksplisit) dan BTCS (implisit) dibandingkan, lalu BTCS dipilih karena stabil tanpa syarat. Sistem persamaan tridiagonal diselesaikan dengan algoritma Thomas, dan non-linearitas radiasi ditangani dengan iterasi Picard. Kode diverifikasi terhadap solusi analitik, neraca energi, uji simetri, dan uji konvergensi, kemudian dibandingkan dengan simulasi ANSYS Transient Thermal. Hasilnya, sisi belakang dinding 15 cm mencapai kenaikan suhu 140 K (kriteria insulasi) pada sekitar 111 menit (kode 110,7 menit; ANSYS sekitar 110,6 menit), dengan selisih suhu kode terhadap ANSYS di bawah 0,4 °C. Waktu ini sangat dipengaruhi tebal dinding dan kondisi batas sisi belakang, serta jauh lebih singkat bila dinding terbakar dari dua sisi. Hasil ini adalah batas insulasi dinding dengan model properti konstan, bukan waktu evakuasi yang utuh.

Kata kunci: konduksi transien, beton, kebakaran ISO 834, BTCS, ANSYS

E. Author Declaration

1. Deep Awareness (of) I

Penulis menyadari bahwa proyek ini bisa selesai karena izin dan bantuan Tuhan Yang Maha Esa sebagai bagian dari pembelajaran serta perkembangan penulis dalam mencapai tujuan kebermanfaatan. Kesadaran ini dijaga sepanjang pengerjaan (CCIT): setiap tahap, dari merumuskan masalah, menurunkan persamaan, membuat kode, sampai menafsirkan grafik, dikerjakan dengan mengingat bahwa hukum fisika yang dipakai adalah ketetapan Sang Pencipta, dan tugas manusia adalah memahaminya untuk menjaga keselamatan sesama.

Hasil ini adalah batas insulasi dinding dengan model properti konstan, bukan waktu evakuasi yang utuh. Penulis memakai asisten AI (Claude) untuk membantu menjelaskan konsep dan menyusun kode, kemudian memeriksa hasilnya sendiri lewat hitungan tangan, uji verifikasi, dan perbandingan dengan ANSYS.

2. Intention of the Project Activity

Mengetahui dan menganalisis metode yang dapat digunakan untuk melihat perpindahan panas pada beton saat kebakaran. Dengan mengetahui seberapa cepat dan seberapa besar panas merambat dari sisi beton di ruang yang terbakar ke sisi beton di ruang lain, diharapkan hasilnya menjadi salah satu acuan waktu penyelamatan dan evakuasi bila terjadi kebakaran dalam kompartemen berdinding beton. Niat ini diarahkan pada kebermanfaatan dan keselamatan jiwa

F. Introduction

1. Latar belakang

Beton banyak dipakai sebagai dinding pemisah antarruang karena tidak mudah terbakar dan lambat menghantarkan panas. Namun “lambat” bukan berarti “tidak sama sekali”. Saat satu ruang terbakar, suhu gas dapat melewati 800–1000 °C dalam setengah jam, sedangkan ruang di baliknya masih berisi orang. Pertanyaan praktisnya: berapa lama dinding itu menahan panas sebelum sisi belakangnya cukup panas untuk menyulut benda mudah terbakar atau melukai orang yang menyentuhnya?

2. Initial Thinking

Penguraian masalah. Kasusnya adalah dinding beton setebal 15 cm yang terbakar di satu sisi. Karena dinding jauh lebih tinggi dan lebar daripada tebalnya, panas hampir seluruhnya mengalir tegak lurus permukaan, sehingga masalah cukup dimodelkan satu dimensi (1D). Yang dicari: suhu sebagai fungsi kedalaman dan waktu, T(x, t), dan terutama waktu ketika sisi belakang naik 140 K.
Pemangku kepentingan. Penghuni ruang sebelah (waktu untuk menyelamatkan diri), petugas pemadam (batas waktu operasi), perancang bangunan (tebal dinding minimum), pemilik gedung (biaya dan keselamatan aset), dan regulator (standar tahan api).
Konteks. Secara fisik, api memanaskan permukaan dengan konveksi dan radiasi. Secara teknis, suhu gas api mengikuti kurva standar ISO 834. Secara sosial, hasilnya berkaitan dengan keselamatan jiwa di bangunan padat.
Akar masalah. Gejalanya adalah “sisi belakang makin panas”. Penyebab mendasarnya: (1) fluks panas dari api sangat besar karena radiasi (berbanding T⁴) pada suhu tinggi; (2) beton punya difusivitas termal rendah (α ≈ 6,5×10⁻⁷ m²/s), sehingga panas merambat lambat tetapi pasti, dengan kedalaman penetrasi sekitar √(αt); (3) tebal dinding dan kondisi di sisi belakang menentukan seberapa cepat panas menumpuk di sana.

3. Prinsip Dasar

Perpindahan panas terjadi lewat tiga mode:

Mode	Hukum	Persamaan
Konduksi	Fourier	q″ = −k ∂T/∂x
Konveksi	Pendinginan Newton	q″ = h (T_fluida − T_permukaan)
Radiasi	Stefan–Boltzmann	q″ = εσ (T_gas⁴ − T_s⁴), suhu dalam Kelvin

Pada suhu api yang tinggi, radiasi mendominasi karena bergantung pada pangkat empat suhu. Bilangan tak berdimensi yang dipakai:

Biot Bi = hL/k (perbandingan hambatan permukaan terhadap hambatan di dalam benda)
Didapat Bi ≈ 2,5 hanya dengan h_konveksi yang mana menunjukkan bahwa Bi jauh diatas 0.1 sehingga dapat dikategorikan sebagai thick wall
Fourier Fo = αt/L² (waktu tak berdimensi: seberapa jauh panas sudah merambat dibanding tebal dinding)
Kedalaman penetrasi δ ≈ √(αt). Dengan h = 25 W/m²K, L = 0,15 m, k = 1,5 W/m·K. Artinya suhu di dalam dinding tidak seragam, sehingga model lumped (0D) tidak berlaku dan persamaan konduksi 1D harus diselesaikan

4. Analisis state of the art

Pendekatan yang ada dapat dikelompokkan sebagai berikut.

Tabulasi preskriptif, misalnya Eurocode (EN 1992-1-2) yang memberi tebal minimum dinding untuk kelas ketahanan api tertentu. Cepat dipakai, tetapi tidak memberi riwayat suhu.
Solusi analitik padatan semi-tak-hingga (fungsi erf), bagus untuk verifikasi tetapi tidak dapat menangani radiasi dan suhu api yang berubah.
Metode numerik (beda hingga, elemen hingga) yang fleksibel, termasuk perangkat lunak seperti ANSYS.

Celah yang diisi laporan ini: menunjukkan secara transparan, dari persamaan sampai kode dan pembandingan dengan perangkat lunak, bagaimana waktu insulasi dinding 15 cm dihitung, serta bagaimana ia dipengaruhi tebal, kondisi sisi belakang, parameter h dan ε, dan jumlah sisi yang terbakar.

G. Methods & Procedures

Idealisasi

Asumsi dan alasan yang diterapkan pada proyek ini

No	Asumsi	Alasan / konsekuensi
1	Aliran panas 1D tegak lurus dinding	Rasio tebal terhadap tinggi/lebar kecil; efek tepi diabaikan
2	Sifat beton konstan: k = 1,5 W/m·K, ρ = 2300 kg/m³, c = 1000 J/kg·K (α = 6,52×10⁻⁷ m²/s)	Penyederhanaan; sifat asli berubah terhadap suhu dan kadar air
3	Suhu awal beton seragam 20 °C	Kondisi sebelum kebakaran
4	Suhu gas api mengikuti ISO 834: Tg = 20 + 345 log₁₀(8t + 1), t dalam menit	Kurva kebakaran standar yang umum dipakai untuk uji ketahanan api
5	Sisi api: konveksi h = 25 W/m²K dan radiasi dengan ε = 0,7	Nilai yang lazim dipakai dalam standar rekayasa kebakaran (Eurocode) untuk kurva standar dan beton
6	Sisi belakang: adiabatik (skenario utama), atau konveksi h_b = 9 W/m²K ke udara 20 °C (skenario pembanding)	Adiabatik tidak melepas panas sama sekali, sehingga suhu sisi belakang lebih tinggi dan waktu insulasi lebih pendek, jadi hasilnya konservatif (berada di sisi aman)
7	Kriteria insulasi: sisi belakang naik 140 K dari suhu awal, yaitu mencapai 20 + 140 = 160 °C	kriteria insulasi pada uji ketahanan api standar (ISO 834-1 dan EN 1363-1): dinding dianggap gagal menyekat panas bila rata-rata suhu sisi tak terbakar naik 140 K dari suhu awal

Model Matematika

Energi masuk – Energi keluar = Energi tersimpan (per waktu)

q″(x)·A − q″(x+Δx)·A = ρ·c·A·Δx·∂T/∂t

−∂q″/∂x = ρc ∂T/∂t

Subtitusi h.fourier,

∂/∂x ( k ∂T/∂x ) = ρc ∂T/∂t

Untuk k konstan, k/(ρc) = α, dan tukar ruas

∂T/∂t = α ∂²T/∂x²

Jadi α (difusivitas),

k : seberapa mudah panas dihantarkan material (angka besar = lbh mudah menghantar)

ρc : energi yg harus diserap untuk naik suhu 1K/volume

Persamaan ini kontinu → suhu terdefinisi di setiap titik x & waktu t (tak hingga banyak nilai) sehingga perlu di diskritisasi (numerik)

Kondisi awal (IC)

T(x, 0) = 20 °C di seluruh dinding.

Kondisi batas sisi api (x = 0)

Pada permukaan, panas yang tiba dari gas (konveksi + radiasi) harus sama dengan panas yang dihantarkan masuk ke beton (kekekalan energi di permukaan):

−k ∂T/∂x = h (Tg − T_s) + εσ (Tg⁴ − T_s⁴)

Suku T⁴ membuat kondisi batas ini non-linear.

Kondisi batas sisi belakang (x = L)

−k ∂T/∂x = h_b (T_s − T_ambient)

Apabila dianggap adiabatik atau tidak diperhitungkan konveksinya, maka h_b = 0

Bidang simetri (dinding terbakar dua sisi)

fluks nol di tengah, ∂T/∂x = 0, sehingga cukup separuh dinding yang dihitung

Metode Numerik

Kenapa numerik?

Karena kalau analitik/eksak kondisinya konstan atau ideal, sedangkan case ini dikategorikan cukup kompleks. Radiasi (T⁴) dan suhu api yang berubah (ISO 834) membuat persamaan tidak punya rumus tertutup. Tetap nyoba make analitik dengan BC sederhana di uji 1

Uji 1 Metode : Analitik (perbandingan/Verifikasi)

kondisi : dinding sgt tebal, properti konstan, suhu permukaan jadi T_s

T(x, t) = T_s + (T₀ − T_s)·erf( x / (2√(αt)) )

Metode : Numerik

FDM (Finite Diff Method) → dinding dipotong dengan jarak Δx ; waktu dipotong dengan langkah Δt ; → persamaan aljabar di tiap titik

Grid : titik x_i = i·Δx ; waktu t_n = n·Δt

FDM ada beberapa skema : FTCS (Forward Time, Central Space) & BTCS (Backward Time, Central Space)

FTCS

∂T/∂t = α ∂²T/∂x²

( T(i,n+1) − T(i,n) ) / Δt = α . ( T(i+1,n) − 2T(i,n) + T(i−1,n) ) / Δx²

dikali Δt pada kedua ruas

( T(i,n+1) − T(i,n) ) = α . Δt/ Δx² . ( T(i+1,n) − 2T(i,n) + T(i−1,n) )

penyederhanaan α . Δt/ Δx² = r

( T(i,n+1) − T(i,n) ) = r . ( T(i+1,n) − 2T(i,n) + T(i−1,n) )

T(i,n+1) = (1 − 2r)·T(i,n) + r·T(i+1,n) + r·T(i−1,n)

Suhu baru adalah rata-rata berbobot dari titik itu dan dua tetangganya

Pengerjaan Manual (Numerik)

α = 6,5×10⁻⁷ ; Δx = 5.09 mm ; Δt = 10 s ; r = 0,25

Pola yang terlihat: panas merambat satu titik per langkah

* Syarat stabilitas

Syarat stabilitas :|1 − 4r| ≤ 1 sama dengan r ≤ 0,5

Syarat dengan kondisi batas konveks i: r·(1 + Bi_Δ) ≤ 0,5

stabilitas pada percobaan r>5 menghasilkan hasil zigzag/negatif

Perbandingan dengan analitik (pada t = 40s)

Titik	Numerik	Analitik	Selisih
1	60,625	58,360	+2,265
2	34,375	32,584	+1,791
3	23,125	22,712	+0,413
4	20,000	20,374	−0,374

Pada selisih titik 4 dapat dilihat perbedaannya, dimana numerik dibatasi 20°C, sedangkan analitik mengandaikan dinding tak hingga.

BTCS

∂T/∂t = α ∂²T/∂x²

( T(i,n+1) − T(i,n) ) / Δt = α . ( T(i+1,n+1) − 2T(i,n+1) + T(i−1,n+1) ) / Δx²

dikali Δt pada kedua ruas

( T(i,n+1) − T(i,n) ) = α . Δt/ Δx² . ( T(i+1,n+1) − 2T(i,n+1) + T(i−1,n+1) )

penyederhanaan α . Δt/ Δx² = r

( T(i,n+1) − T(i,n) ) = r . ( T(i+1,n+1) − 2T(i,n+1) + T(i−1,n+1) )

T(i,n) = (1 − 2r)·T(i,n+1) + r·T(i+1,n+1) + r·T(i−1,n+1)

Pengerjaan Numerik secara manual – dapat juga diubah jadi matriks

Sel setengah di tepi

Titik di permukaan hanya punya beton di satu sisi, sehingga selnya setengah (lebar Δx/2). Neraca energi di titik 0

ρc(Δx/2)·dT₀/dt = h(Tg − T₀) + k(T₁ − T₀)/Δx

Massa yang menyimpan energi hanya separuh, sehingga kapasitas panasnya separuh tetapi energi yang di dapat sama sehingga suhunya lebih besar 2 kali lipat (koef nya jadi 2r, dengan Bi_Δ = hΔx/k)

• FTCS: T₀(n+1) = [1 − 2r(1 + Bi_Δ)]·T₀(n) + 2r·T₁(n) + 2r·Bi_Δ·Tg

• BTCS: [1 + 2r(1 + Bi_Δ)]·T₀(n+1) − 2r·T₁(n+1) = T₀(n) + 2r·Bi_Δ·Tg

Radiasi

εσ(Tg⁴ − Ts⁴) = (Tg − Ts)(Tg + Ts)(Tg² + Ts²)εσ

h_rad = εσ(Tg² + Ts²)(Tg + Ts)

h efektif = h + h_rad

Masalahnya h_rad bergantung pada Ts yang sedang dicari. Solusinya iterasi Picard: tebak Ts (awalnya suhu permukaan langkah sebelumnya, yang cukup dekat karena langkahnya kecil), hitung h_rad, selesaikan sistem, lalu perbarui Ts; ulangi sampai Ts tidak berubah

Titik sisi belakang (N)

Adiabatik dan bandingannya apabila titik belakang dianggap konveksi

Tabel Ringkasan

Titik	FTCS	BTCS
Dalam	T(i,n+1) = (1−2r)T(i,n) + r·T(i+1,n) + r·T(i−1,n)	−r·T(i−1,n+1) + (1+2r)·T(i,n+1) − r·T(i+1,n+1) = T(i,n)
0, api	[1−2r(1+Bi_Δ)]T₀ + 2r·T₁ + 2r·Bi_Δ·Tg	[1+2r(1+Bi_Δ)]T₀(n+1) − 2r·T₁(n+1) = T₀(n) + 2r·Bi_Δ·Tg
N, adiabatik	(1−2r)T_N + 2r·T(N−1)	(1+2r)T_N(n+1) − 2r·T(N−1)(n+1) = T_N(n)
N, konveksi	[1−2r(1+Bi_b)]T_N + 2r·T(N−1) + 2r·Bi_b·T_amb	[1+2r(1+Bi_b)]T_N(n+1) − 2r·T(N−1)(n+1) = T_N(n) + 2r·Bi_b·T_amb

4. Instruksi set (proses, langkah demi langkah)

Adapun instruksi setnya :

Tetapkan sifat material, geometri (L, Δx, Δt), dan kondisi awal 20 °C di semua titik.
Untuk setiap langkah waktu t = nΔt, hitung suhu gas ISO 834 (rumus memakai menit).
Tebak suhu permukaan (pakai nilai langkah sebelumnya), hitung h efektif = h + h_rad, lalu susun matriks BTCS tridiagonal.
Selesaikan dengan algoritma Thomas: sapuan maju (eliminasi variabel satu per satu sehingga baris terakhir hanya memuat satu bilangan tak diketahui), lalu substitusi mundur. Biaya komputasinya sebanding dengan jumlah titik.
Bandingkan suhu permukaan baru dengan taksiran. Jika berbeda, ulangi dari langkah 3 (iterasi Picard); jika sama, lanjut.
Ulangi dari langkah 2 sampai waktu akhir simulasi (120 menit untuk profil; 240 menit untuk mencari waktu insulasi).
Catat waktu pertama saat T sisi belakang ≥ 160 °C (interpolasi linear antarlangkah), lalu buat grafik.

Perangkat: Python 3 (NumPy, modul math), VS Code, ANSYS Student 2026 R1 (Transient Thermal), dengan Δx = 2,5 mm dan Δt = 5 s sebagai pengaturan utama.

Verifikasi dan validasi.

(1) solusi analitik, termasuk BC konveksi pada sel setengah (Uji 1);

(2) konvergensi orde dua;

(3) neraca energi (Uji 2);

(4) uji simetri setengah domain (Uji 4);

(5) perbandingan dengan ANSYS

Iterasi dan keberlanjutan.

Pengerjaan dilakukan bertahap dari kasus sederhana (properti konstan, suhu permukaan konstan) ke realistis (konveksi + radiasi + ISO 834), dan tiap lompatan diuji dulu supaya penyebab kesalahan mudah dilacak.

Kode dibuat modular (kondisi batas bisa diganti tanpa mengubah solver) sehingga bisa dikembangkan, misalnya ke sifat bergantung suhu atau tebal dan kondisi batas lain.

H. Results & Discussion

Verifikasi

Uji 1: Solusi Analitik (konveksi saja, Tgas = 500 °C, h = 25, t = 600 s).

Δx (mm)	Δt (s)	Error/Galat maks (°C)	Error/Galat perm (°C)
5	4.0	0.4027	0.3416
2.5	1.0	0.1002	0.00847

Galat turun sekitar 4 kali lipat saat Δx dibagi dua, sesuai orde dua. Uji pendahuluan FTCS terhadap solusi erf pada t = 40 s menunjukkan tren yang sama : semakin halus grid, semakin kecil galat.

Uji 2: Neraca Energi

(ISO 834 + radiasi, 120 menit): energi masuk dari api = energi tersimpan = 152,8422 MJ/m², galat relatif 4,7×10⁻¹³. Uji ini memeriksa konsistensi kode, bukan kecocokan dengan kenyataan.

Uji 4: Perbandingan tiga skema

(kasus 5 titik, r = 0,25, Δx = 5,1 mm, Δt = 10 s; contoh pada Δt kasar)

T(i,n)	FTCS	BTCS	Analitik
T(1,10)	40,00	33,73	32,58
T(2,10)	20,00	22,35	20,37
T(3,10)	20,00	20,39	20,00
T(1,40)	60,63	55,79	58,36
T(2,40)	34,38	32,75	32,58
T(3,40)	23,13	23,69	22,71

Hasil Utama

Pada dinding 15cm dengan kondisi satu sisi terdapat kebakaran dan sisi lainnya dianggap adiabatik, h = 25, ε = 0,7

x (mm)	T_30mnt	T_60mnt	T_90mnt	T_120mnt
0	719	874,9	954,5	1007,8
20	420,1	624,1	735,1	810,5
40	220,9	420,6	545,5	635,1
60	108,1	268,7	390,8	486,3
80	53,7	164,3	272,1	366,9
100	31,2	98,8	187,7	277,9
150	20,8	44,1	105,3	186,2
Tgas ISO 834	842	945	1006	1049

Pada saat sisi belakang sampai di T = 140k, apabila sisi belakang adiabatik maka temperature tercapai di t = 110,7 menit, apabila sisi belakang dianggap konveksi (h = 9) maka temperature tercapai di t = 122,4 menit

Grafik

suhu permukaan beton selalu di bawah suhu gas, karena panas mengalir ke permukaan hanya selama gas lebih panas dari permukaan, lalu harus dihantarkan ke dalam beton yang konduktivitasnya rendah. Selisihnya mengecil seiring waktu

garis curam di dekat permukaan dan datar di dekat 150 mm menunjukkan gradien tajam (perilaku thick). Ujung kiri tiap garis sama dengan suhu permukaan pada Grafik 1. Suhu naik seiring waktu karena panas terus masuk dan menumpuk.

Dalam rentang h dan ε yang realistis, waktu insulasi bergeser hanya sedikit. Terlihat radiasi menyumbang sekitar 88% perpindahan panas di permukaan dan Bi besar, sehingga hambatan di dalam beton yang dominan. Hanya saat radiasi dimatikan (ε = 0, tidak realistis untuk beton) waktunya bergeser jauh.

Verifikasi/Validasi terhadap ANSYS Transient Thermal

Simulasi ANSYS (Student 2026 R1, Transient Thermal, beban konveksi dan radiasi pada sisi api sesuai ISO 834, 7200 s) memberikan hasil berikut dibanding kode (h = 25, ε = 0,7):

Waktu	T_x0, coding	T_x0, ansys	Selisih	T_xL, coding	T_xL, ansys	Selisih
30 mnt (1800 s)	719,03	718,64	−0,39	20,76	20,85	0,09
60 mnt (3600 s)	874,88	874,75	−0,13	44,09	44,48	0,39
90 mnt (5400 s)	954,52	954,46	−0,07	105,3	105,69	0,39
120 mnt (7200 s)	1007,81	1007,8	−0,01	186,25	186,53	0,28

Selisih maksimum sekitar 0,39 °C di muka api maupun di muka belakang. Dari data keluaran ANSYS, suhu minimum (sisi belakang) melewati 160 °C antara 6631 s (159,9 °C) dan 6651 s (160,83 °C), yaitu sekitar 110,6 menit, sangat dekat dengan 110,7 menit dari kode. Kedua solver yang dikembangkan terpisah menghasilkan jawaban yang hampir identik, sehingga implementasi kode dapat dipercaya.

I. Conclusion, Closing Remarks, Recommendations

Dapat disimpulkan beberapa hal dari langkah hinggal hasil yang di dapat :

Waktu sisi belakang dinding beton 15 cm naik 140 K adalah sekitar 111 menit (kode 110,7; ANSYS 110,6), dengan selisih suhu kode terhadap ANSYS di bawah 0,4 °C.
BTCS dipilih karena stabil tanpa syarat, sedangkan FTCS meledak untuk r > 0,5. Radiasi ditangani dengan h efektif dan iterasi Picard, dan sel setengah (koefisien 2r) di tepi terbukti benar lewat analitik dan neraca energi.
Waktu insulasi sangat bergantung pada tebal (57, 111, 181 menit untuk 10, 15, 20 cm), kondisi batas sisi belakang (adiabatik vs konveksi bergeser sekitar 8–14%), dan jumlah sisi terbakar (dua sisi: sekitar 37 menit di bidang tengah), tetapi hampir tidak bergantung pada nilai h dan ε yang realistis (sekitar ±3%).
Hasil ini adalah batas insulasi pada model properti konstan, bukan waktu evakuasi yang utuh.

Laporan dan pengujian ini masih terbatas dan dapat ditingkatkan lagi dengan menghilangkan batas yang ditetapkan disini seperti sifat kebakaran pada beton yang akan bergantung dan mempengaruhi kandungan airnya, terjadinya retak atau kekuatan struktur yang menurun, serta dilakukannya validasi langsung dengan data uji kebakaran atau tabel eurocode.

J. Acknowledgments

Terimakasih kepada rekan-rekan satu kelas yang telah berusaha belajar dan berproses bersama, kepada Prof DAI yang telah membersamai proses penyelesaian laporan ini, serta rekan dan orang tua yang mendukung selama perjalanan laporan ini.

K. (References) Literature Cited

Incropera, F. P., DeWitt, D. P., Bergman, T. L., Lavine, A. S. Fundamentals of Heat and Mass Transfer. Wiley.
Chapra, S. C., Canale, R. P. Numerical Methods for Engineers. McGraw-Hill.
Buchanan, A. H., Abu, A. K. Structural Design for Fire Safety. Wiley.
Kodur, V. (2014). Properties of Concrete at Elevated Temperatures. ISRN Civil Engineering, 2014.
ISO 834-1: Fire-resistance tests, elements of building construction, general requirements.
EN 1991-1-2 : Actions on structures exposed to fire.
EN 1992-1-2 : Design of concrete structures, structural fire design.
ANSYS, Inc. ANSYS Mechanical, Transient Thermal Analysis Guide, rilis 2026 R1.

L. Appendices

Lampiran 1 – Daftar Simbol

Aturan notasi: i selalu nomor titik (lokasi), n selalu nomor langkah waktu. Suhu di titik i pada langkah n ditulis T(i,n), dan waktunya t = n·Δt. Mulai dari n = 0 untuk t = 0.

Simbol	Arti	Satuan
T	suhu	°C (radiasi ; K, K = °C + 273,15)
t	waktu	s
x	posisi diukur dari permukaan api	m
L	tebal dinding	m
k	konduktivitas termal	W/m·K
ρ	densitas	kg/m³
c	kalor jenis	J/kg·K
α = k/(ρc)	difusivitas termal	m²/s
q″	fluks panas (laju panas per luas)	W/m²
h	koefisien konveksi	W/m²K
Tg	suhu gas api	°C
ε	emisivitas permukaan	–
σ	konstanta Stefan–Boltzmann, 5,67×10⁻⁸	W/m²K⁴
Δx, Δt	jarak antar titik, langkah waktu	m, s
r = αΔt/Δx²	angka difusi per langkah	–
Bi = hL/k	bilangan Biot (seluruh dinding)	–
Bi_Δ = hΔx/k	bilangan Biot sel (satu sel grid)	–
Fo = αt/L²	bilangan Fourier (waktu tak berdimensi)	–

Lampiran 2 – Kode program

"""
solver_btcs.py
Solver 1D konduksi transien pada dinding beton: BTCS + algoritma Thomas + sel setengah di kedua tepi.

  Sisi kiri  (x = 0): api  -> konveksi + radiasi, suhu gas mengikuti kurva ISO 834
  Sisi kanan (x = L): adiabatik / konveksi / api (api di dua sisi dipakai untuk uji simetri)

Radiasi ditangani dengan trik "h efektif":
  eps*sigma*(Tg^4 - Ts^4) = h_rad*(Tg - Ts),  h_rad = eps*sigma*(Tg^2 + Ts^2)*(Tg + Ts)   (suhu dalam Kelvin)
Jadi radiasi terlihat seperti konveksi dengan koefisien yang bergantung pada Ts. Karena Ts belum
diketahui, kita ITERASI (Picard): tebak Ts -> hitung h_rad -> selesaikan sistem -> perbarui Ts -> ulangi.

Hanya butuh numpy dan modul math bawaan Python.
"""
import time
import numpy as np
from math import erfc, exp, sqrt

# ----------------------------- Properti (konstan) -----------------------------
k = 1.5            # W/m.K   konduktivitas beton
rho = 2300.0       # kg/m3
c = 1000.0         # J/kg.K
alpha = k / (rho * c)
sigma = 5.67e-8    # W/m2K4
T_AWAL = 20.0      # C


def iso834(t_detik):
    """Suhu gas api ISO 834 (C). Rumus memakai waktu dalam menit."""
    return 20.0 + 345.0 * np.log10(8.0 * t_detik / 60.0 + 1.0)


class Sisi:
    """Kondisi batas konveksi + radiasi di satu sisi. h=0 dan eps=0 berarti adiabatik."""

    def __init__(self, h=0.0, eps=0.0, T_fluida=lambda t: 20.0):
        self.h = h
        self.eps = eps
        self.T_fluida = T_fluida

    def h_efektif(self, Tf, Ts):
        """h konveksi + h radiasi terlinearisasi, W/m2K. Tf dan Ts dalam C."""
        if self.eps == 0.0:
            return self.h
        TfK, TsK = Tf + 273.15, Ts + 273.15
        return self.h + self.eps * sigma * (TfK**2 + TsK**2) * (TfK + TsK)

    def fluks_masuk(self, Tf, Ts):
        """Fluks panas netto yang MASUK ke beton, W/m2 (positif = masuk)."""
        q = self.h * (Tf - Ts)
        if self.eps != 0.0:
            q += self.eps * sigma * ((Tf + 273.15) ** 4 - (Ts + 273.15) ** 4)
        return q

    @property
    def nonlinear(self):
        return self.eps != 0.0


def adiabatik():
    return Sisi()


def konveksi(h, T_udara=20.0):
    return Sisi(h=h, eps=0.0, T_fluida=lambda t: T_udara)


def api(h=25.0, eps=0.7):
    return Sisi(h=h, eps=eps, T_fluida=iso834)


def thomas(a, b, cc, d):
    """Selesaikan a[i]*x[i-1] + b[i]*x[i] + cc[i]*x[i+1] = d[i] (a[0]=0, cc[-1]=0)."""
    n = len(d)
    c_baru = np.zeros(n)
    d_baru = np.zeros(n)
    c_baru[0] = cc[0] / b[0]
    d_baru[0] = d[0] / b[0]
    for i in range(1, n):                                  # sapuan maju
        penyebut = b[i] - a[i] * c_baru[i - 1]
        c_baru[i] = cc[i] / penyebut
        d_baru[i] = (d[i] - a[i] * d_baru[i - 1]) / penyebut
    x = np.zeros(n)
    x[-1] = d_baru[-1]
    for i in range(n - 2, -1, -1):                         # substitusi mundur
        x[i] = d_baru[i] - c_baru[i] * x[i + 1]
    return x


def selesaikan(L, dx, dt, t_akhir, kiri, kanan, catat_menit=(), tol=1e-9):
    """
    Jalankan BTCS dari t = 0 sampai t_akhir (detik).
    Titik: i = 0..N (N+1 titik), titik 0 dan N adalah sel setengah.
    Mengembalikan dict berisi deret waktu, profil di waktu yang dicatat, dan neraca energi.
    """
    N = int(round(L / dx))
    M = N + 1
    r = alpha * dt / dx**2
    n_langkah = int(round(t_akhir / dt))
    T = np.full(M, T_AWAL)

    waktu = np.arange(n_langkah + 1) * dt
    T_kiri = np.zeros(n_langkah + 1); T_kiri[0] = T[0]
    T_kanan = np.zeros(n_langkah + 1); T_kanan[0] = T[-1]
    profil = {}
    Q_kiri = 0.0
    Q_kanan = 0.0
    nonlinear = kiri.nonlinear or kanan.nonlinear

    for n in range(1, n_langkah + 1):
        t = n * dt
        Tf_L = kiri.T_fluida(t)
        Tf_R = kanan.T_fluida(t)
        T_iter = T.copy()
        for _ in range(100):                               # iterasi Picard
            Bi_L = kiri.h_efektif(Tf_L, T_iter[0]) * dx / k
            Bi_R = kanan.h_efektif(Tf_R, T_iter[-1]) * dx / k
            a = np.full(M, -r)
            b = np.full(M, 1.0 + 2.0 * r)
            cc = np.full(M, -r)
            d = T.copy()                                   # ruas kanan: suhu LAMA
            a[0] = 0.0
            b[0] = 1.0 + 2.0 * r * (1.0 + Bi_L)            # baris titik 0 (sel setengah, api)
            cc[0] = -2.0 * r
            d[0] = T[0] + 2.0 * r * Bi_L * Tf_L
            a[-1] = -2.0 * r                               # baris titik N (sel setengah, belakang)
            b[-1] = 1.0 + 2.0 * r * (1.0 + Bi_R)
            cc[-1] = 0.0
            d[-1] = T[-1] + 2.0 * r * Bi_R * Tf_R
            T_baru = thomas(a, b, cc, d)
            selisih = max(abs(T_baru[0] - T_iter[0]), abs(T_baru[-1] - T_iter[-1]))
            T_iter = T_baru
            if (not nonlinear) or selisih < tol:
                break
        T = T_iter
        Q_kiri += dt * kiri.fluks_masuk(Tf_L, T[0])
        Q_kanan += dt * kanan.fluks_masuk(Tf_R, T[-1])
        T_kiri[n] = T[0]
        T_kanan[n] = T[-1]
        for m in catat_menit:
            if abs(t - m * 60.0) < 1e-9:
                profil[m] = T.copy()

    # neraca energi: energi tersimpan (sel setengah di ujung dihitung separuh)
    bobot = np.ones(M); bobot[0] = 0.5; bobot[-1] = 0.5
    E_simpan = rho * c * dx * np.sum(bobot * (T - T_AWAL))
    return dict(x=np.arange(M) * dx, waktu=waktu, T_kiri=T_kiri, T_kanan=T_kanan, T_akhir=T,
                profil=profil, Q_masuk=Q_kiri + Q_kanan, E_simpan=E_simpan, r=r)


def waktu_naik(waktu, T_seri, naik=140.0):
    """Waktu pertama (menit) saat T_seri >= T_AWAL + naik, dengan interpolasi linear. None jika tak tercapai."""
    target = T_AWAL + naik
    idx = np.where(T_seri >= target)[0]
    if len(idx) == 0:
        return None
    j = idx[0]
    if j == 0:
        return 0.0
    t1, t2, y1, y2 = waktu[j - 1], waktu[j], T_seri[j - 1], T_seri[j]
    return (t1 + (target - y1) / (y2 - y1) * (t2 - t1)) / 60.0


def analitik_konveksi(x, t, h, T_inf, T_i):
    """Solusi padatan semi-infinit dengan BC konveksi (Incropera, bagian padatan semi-infinit)."""
    eta = x / (2.0 * sqrt(alpha * t))
    beta = h * sqrt(alpha * t) / k
    return T_i + (T_inf - T_i) * (erfc(eta) - exp(2.0 * eta * beta + beta**2) * erfc(eta + beta))


def nilai_di(x_arr, T_arr, x_mm):
    """Ambil suhu di kedalaman tertentu (mm) dari titik grid terdekat."""
    i = int(round((x_mm / 1000.0) / (x_arr[1] - x_arr[0])))
    return T_arr[i]


# ================================ UJI-UJI ================================
if __name__ == "__main__":
    mulai = time.time()

    # ---- UJI 1: verifikasi BC konveksi pada sel setengah dengan solusi analitik ----
    print("UJI 1: konveksi saja (Tgas = 500 C, h = 25) vs solusi analitik semi-infinit, t = 600 s")
    print(f"{'dx (mm)':>8} {'dt (s)':>7} {'error maks (C)':>15} {'error permukaan (C)':>20}")
    kiri_uji = konveksi(25.0, 500.0)
    for dx_u, dt_u in [(0.005, 4.0), (0.0025, 1.0)]:
        hasil = selesaikan(0.30, dx_u, dt_u, 600.0, kiri_uji, adiabatik())
        T_an = np.array([analitik_konveksi(xi, 600.0, 25.0, 500.0, T_AWAL) for xi in hasil["x"]])
        err = np.abs(hasil["T_akhir"] - T_an)
        print(f"{dx_u*1000:8.2f} {dt_u:7.1f} {err.max():15.4f} {err[0]:20.4f}")

    # ---- UJI 2: neraca energi (ISO 834, radiasi, satu sisi terbakar) ----
    L = 0.15
    dx, dt = 0.0025, 5.0
    hasil2 = selesaikan(L, dx, dt, 120 * 60.0, api(), adiabatik())
    galat = abs(hasil2["Q_masuk"] - hasil2["E_simpan"]) / hasil2["Q_masuk"]
    print("\nUJI 2: neraca energi, dinding 15 cm, ISO 834 + radiasi, sisi belakang adiabatik, 120 menit")
    print(f"  Energi masuk dari api : {hasil2['Q_masuk']/1e6:10.4f} MJ/m2")
    print(f"  Energi tersimpan      : {hasil2['E_simpan']/1e6:10.4f} MJ/m2")
    print(f"  Galat relatif         : {galat:10.2e}")

    # ---- UJI 3: hasil utama beton satu sisi terbakar ----
    print("\nUJI 3: profil suhu dinding 15 cm (adiabatik di belakang), ISO 834 + radiasi")
    menit = (30, 60, 90, 120)
    hasil3 = selesaikan(L, dx, dt, 120 * 60.0, api(), adiabatik(), catat_menit=menit)
    kedalaman = (0, 20, 40, 60, 80, 100, 150)
    print(f"{'kedalaman (mm)':>15} " + " ".join(f"{str(m)+' mnt':>9}" for m in menit))
    for xm in kedalaman:
        print(f"{xm:15d} " + " ".join(f"{nilai_di(hasil3['x'], hasil3['profil'][m], xm):9.1f}" for m in menit))
    print("  Suhu gas ISO 834: " + ", ".join(f"{m} mnt = {iso834(m*60.0):.0f} C" for m in menit))
    print("  Kedalaman penetrasi sqrt(alpha*t): " + ", ".join(f"{m} mnt = {sqrt(alpha*m*60.0)*1000:.0f} mm" for m in menit))

    print("\n  Waktu sisi belakang naik 140 K (simulasi sampai 240 menit):")
    print(f"  {'tebal':>7} {'adiabatik (mnt)':>16} {'konveksi h=9 (mnt)':>20}")
    for tebal in (0.10, 0.15, 0.20):
        h_ad = selesaikan(tebal, dx, dt, 240 * 60.0, api(), adiabatik())
        h_kv = selesaikan(tebal, dx, dt, 240 * 60.0, api(), konveksi(9.0, 20.0))
        w_ad = waktu_naik(h_ad["waktu"], h_ad["T_kanan"])
        w_kv = waktu_naik(h_kv["waktu"], h_kv["T_kanan"])
        f = lambda w: f"{w:.1f}" if w is not None else "tidak tercapai"
        print(f"  {tebal*100:5.0f} cm {f(w_ad):>16} {f(w_kv):>20}")

    # ---- UJI 4: dua sisi terbakar, setengah domain vs domain penuh ----
    print("\nUJI 4: dua sisi terbakar simetris, dinding 15 cm: domain penuh vs setengah domain (adiabatik di tengah)")
    penuh = selesaikan(L, dx, dt, 120 * 60.0, api(), api())
    setengah = selesaikan(L / 2, dx, dt, 120 * 60.0, api(), adiabatik())
    tengah_penuh = penuh["T_akhir"][len(penuh["T_akhir"]) // 2]
    # deret waktu: permukaan api dan titik tengah
    mid_penuh = []
    # hitung ulang deret titik tengah domain penuh lewat run kedua yang mencatat T_kanan setengah domain
    beda_permukaan = np.max(np.abs(penuh["T_kiri"] - setengah["T_kiri"]))
    beda_tengah = abs(tengah_penuh - setengah["T_kanan"][-1])
    print(f"  Beda maksimum suhu permukaan api sepanjang waktu : {beda_permukaan:.2e} C")
    print(f"  Beda suhu tengah dinding pada 120 menit          : {beda_tengah:.2e} C")
    print(f"  Suhu tengah dinding pada 120 menit (setengah domain): {setengah['T_kanan'][-1]:.1f} C")
    w_tengah = waktu_naik(setengah["waktu"], setengah["T_kanan"])
    print(f"  Waktu tengah dinding naik 140 K (setengah domain): "
          f"{'tidak tercapai dalam 120 menit' if w_tengah is None else f'{w_tengah:.1f} menit'}")

    print(f"\nWaktu komputasi total: {time.time() - mulai:.1f} detik")

HASIL

UJI 1: konveksi saja (Tgas = 500 C, h = 25) vs solusi analitik semi-infinit, t = 600 s

dx (mm) dt (s) error maks (C) error permukaan (C)

5.00 4.0 0.4027 0.3416

2.50 1.0 0.1002 0.0847

UJI 2: neraca energi, dinding 15 cm, ISO 834 + radiasi, sisi belakang adiabatik, 120 menit

Energi masuk dari api : 152.0572 MJ/m2

Energi tersimpan : 152.0572 MJ/m2

Galat relatif : 4.92e-13

UJI 3: profil suhu dinding 15 cm (adiabatik di belakang), ISO 834 + radiasi

kedalaman (mm) 30 mnt 60 mnt 90 mnt 120 mnt

0 710.5 871.7 952.6 1006.5

20 410.8 619.8 732.5 808.6

40 213.7 416.1 542.5 632.8

60 103.8 264.6 387.8 483.9

80 51.5 161.2 269.3 364.4

100 30.3 96.5 185.2 275.4

150 20.7 43.0 103.3 183.9

Suhu gas ISO 834: 30 mnt = 842 C, 60 mnt = 945 C, 90 mnt = 1006 C, 120 mnt = 1049 C

Kedalaman penetrasi sqrt(alpha*t): 30 mnt = 34 mm, 60 mnt = 48 mm, 90 mnt = 59 mm, 120 mnt = 69 mm

Waktu sisi belakang naik 140 K (simulasi sampai 240 menit):

tebal adiabatik (mnt) konveksi h=9 (mnt)

10 cm 58.0 62.3

15 cm 111.5 123.3

20 cm 182.0 206.9

UJI 4: dua sisi terbakar simetris, dinding 15 cm: domain penuh vs setengah domain (adiabatik di tengah)

Beda maksimum suhu permukaan api sepanjang waktu : 1.93e-11 C

Beda suhu tengah dinding pada 120 menit : 7.05e-12 C

Suhu tengah dinding pada 120 menit (setengah domain): 724.3 C

Waktu tengah dinding naik 140 K (setengah domain): 37.6 menit

Waktu komputasi total: 10.5 detik

Lampiran 3 – Keluaran Ansys

Parameter ANSYS

Item	Nilai
Geometri	balok 150 × 10 × 10 mm (150 mm = arah tebal); sisi samping tanpa beban (adiabatik)
Material	k = 1,5 W/m·K; ρ = 2300 kg/m³; c = 1000 J/kg·K
Mesh	Element Size 2,5 mm
Suhu awal	20 °C
Waktu	Step End Time 7200 s; Auto Time Stepping aktif
Sisi api	Convection h = 25 W/m²·°C dan Radiation (Surface to Ambient) ε = 0,7, keduanya dengan suhu lingkungan fungsi =20+345log10(8time/60+1) (time dalam detik)
Sisi belakang	tanpa beban (adiabatik)
Hasil	Temperature (Maximum = muka api, Minimum = muka belakang)

Lampiran 4

ccitonline.com

Quick Links