Implementasi proyek

Membangun fungsi proyek, inti tervektorisasi, dan figur hasil analisis.

Daftar isi

Spesifikasi yang telah ditetapkan kini diterjemahkan menjadi kode yang berjalan. Topik ini membahas arsitektur proyek yang memisahkan fungsi perhitungan dari skrip alur, serta linierisasi model orde satu menjadi regresi linier yang diselesaikan dengan operator backslash. Fungsi penyesuai fit_termal dibangun untuk merangkum pencarian kandidat $T_{\text{aset}}$ , skrip utama analisis_termal merangkai seluruh alur sampai figur hasil, dan verifikasi vektorisasi diterapkan pada perhitungan error per titik. Setelah membaca topik ini, seluruh keluaran yang dijanjikan spesifikasi telah dihasilkan dengan angka yang nyata.

Arsitektur dua berkas proyek#

Proyek ini memisahkan perhitungan dari pengaturan alur pada dua berkas yang berperan berbeda. Berkas fit_termal.m adalah fungsi yang menerima data pengukuran lalu mengembalikan struktur berisi parameter model beserta mutunya. Berkas analisis_termal.m adalah skrip yang memuat berkas data, memanggil fungsi itu, memeriksa hasilnya, dan menyimpan figur. Pemisahan ini mengikuti prinsip yang sama dengan pengelolaan berkas proyek yang rapi: berkas berubah dengan ritme yang berbeda. Fungsi yang dapat diuji terpisah juga lebih mudah dipercaya daripada satu skrip besar yang mencampur semuanya.

Fungsi penyesuai sengaja tidak membaca berkas datanya sendiri. Masukan berupa pasangan vektor waktu dan suhu membuat fungsi dapat diuji dengan data apa pun, termasuk data contoh yang dipakai ketika argumen tidak diberikan. Sebaliknya, skrip yang membaca berkas tetap menjadi penghubung antara data dan perhitungan, sehingga penggantian sumber data tidak pernah menyentuh rumus. Dengan pembagian ini, setiap berkas pendek dan setiap tanggung jawabnya tunggal.

Linierisasi model menjadi regresi linier#

Model orde satu bukan model linier pada parameternya, tetapi dapat dilinierkan dengan satu langkah aljabar. Mengambil logaritma alami dari selisih suhu terhadap asimtot menghasilkan

\ln(T_{\text{aset}} - T) = \ln(T_{\text{aset}} - T_{\text{awal}}) - \frac{t}{\tau}

sehingga pasangan $(t,\ \ln(T_{\text{aset}} - T))$ terletak pada satu garis lurus dengan kemiringan $-1/\tau$ dan intersep $\ln(T_{\text{aset}} - T_{\text{awal}})$ . Kuadrat terkecil untuk garis itu dikerjakan operator backslash pada matriks rancangan $[t\ \ 1]$ , tepat seperti regresi sensor pada modul aljabar linear.

Langkah aljabar itu membawa satu masalah: nilai $T_{\text{aset}}$ justru belum diketahui padahal dibutuhkan untuk menghitung logaritmanya. Pendekatan yang jujur dilakukan dua tahap. Tahap pertama menetapkan beberapa kandidat nilai $T_{\text{aset}}$ yang tetap, tahap kedua menjalankan regresi linier untuk tiap kandidat lalu menghitung error RMS modelnya terhadap data, dan kandidat dengan RMS terkecil dipilih. Cara ini sederhana, deterministik, dan tidak menyembunyikan asumsi apa pun di balik pemecah nonlinear. Blok berikut mengerjakan satu kandidat lebih dahulu agar mekanismenya terlihat sebelum dibungkus fungsi.

matlab Linierisasi untuk satu kandidat T_aset

% muat data pengukuran lalu linierisasi untuk satu kandidat T_aset
M = csvread('data-suhu.csv');
t = M(:, 1);
T = M(:, 2);

T_kandidat = 85;                    % nilai asimtot di atas suhu tertinggi
y = log(T_kandidat - T);            % menjadi linear terhadap t

A = [t, ones(size(t))];             % matriks rancangan
c = A \ y;                          % kuadrat terkecil

fprintf('kemiringan : %.6f per detik\n', c(1));
fprintf('intersep    : %.6f\n', c(2));
fprintf('tau         : %.4f detik\n', -1 / c(1));
fprintf('T_awal      : %.4f derajat\n', T_kandidat - exp(c(2)));

text Keluaran

kemiringan : -0.033226 per detik
intersep    : 4.087612
tau         : 30.0967 detik
T_awal      : 25.4026 derajat

Kemiringan −0.033226 per detik menghasilkan konstanta waktu $\tau = -1/\text{kemiringan} = 30.0967$ detik, dan intersep 4.087612 adalah $\ln(85 - T_{\text{awal}})$ sehingga $T_{\text{awal}} = 85 - e^{4.087612} = 25.4026$ derajat. Satu kandidat telah menghasilkan model yang utuh, tetapi kecocokannya terhadap data belum diukur dan kandidat lain belum diuji, sehingga belum ada dasar untuk memilih.

Fungsi penyesuai fit_termal#

Fungsi fit_termal merangkum kedua tahap menjadi satu panggilan. Daftar kandidat dibawa sebagai argumen ketiga yang bersifat pilihan, data contoh dipakai bila argumen tidak diberikan, dan seluruh hasil dikemas menjadi struktur agar pemanggil menerima parameter beserta catatan kandidatnya sekaligus.

matlab fit_termal.m

function hasil = fit_termal(t, T, kandidat)
% FIT_TERMAL menyesuaikan model pemanasan orde satu pada data suhu.
%   hasil = FIT_TERMAL(t, T) mencari parameter model
%     T(t) = T_aset - (T_aset - T_awal) * exp(-t/tau)
%   dengan menguji daftar kandidat T_aset tetap lalu regresi linear
%   pada ln(T_aset - T) terhadap t memakai operator backslash.
%   Kandidat dengan error RMS terkecil dipilih.
%   hasil = FIT_TERMAL(t, T, kandidat) memakai daftar kandidat sendiri.
%   Struktur hasil memuat bidang T_aset, T_awal, tau, rms, T_model,
%   kandidat, tau_kandidat, dan rms_kandidat.
%   Tanpa argumen, data contoh pemanasan logam dipakai.

if nargin < 1
  t = (0:10:120)';
  T = [25.00; 42.21; 54.05; 63.10; 69.25; 73.45; 77.00; ...
       79.40; 80.65; 82.15; 82.75; 83.60; 83.80];
end
if nargin < 3
  kandidat = [85 90 95 100];
end

t = t(:);
T = T(:);

% buang kandidat yang tidak melewati suhu terukur tertinggi,
% karena ln(kandidat - T) hanya terdefinisi untuk nilai positif
kandidat = kandidat(kandidat > max(T));

A = [t, ones(size(t))];           % matriks rancangan regresi linear
n = numel(kandidat);
tau_kandidat = zeros(n, 1);
rms_kandidat = zeros(n, 1);
T_awal_kandidat = zeros(n, 1);

for k = 1:n
  y = log(kandidat(k) - T);       % linierisasi logaritma
  c = A \ y;                      % kuadrat terkecil: [kemiringan; intersep]
  tau_kandidat(k) = -1 / c(1);
  T_awal_kandidat(k) = kandidat(k) - exp(c(2));
  T_model = kandidat(k) - exp(c(2)) * exp(-t / tau_kandidat(k));
  rms_kandidat(k) = sqrt(mean((T_model - T).^2));
end

[rms, terbaik] = min(rms_kandidat);

hasil.T_aset = kandidat(terbaik);
hasil.T_awal = T_awal_kandidat(terbaik);
hasil.tau = tau_kandidat(terbaik);
hasil.rms = rms;
hasil.T_model = hasil.T_aset - ...
    (hasil.T_aset - hasil.T_awal) * exp(-t / hasil.tau);
hasil.kandidat = kandidat;
hasil.tau_kandidat = tau_kandidat;
hasil.rms_kandidat = rms_kandidat;
end

Tiga keputusan pada fungsi itu patut dicatat. Matriks rancangan dibangun sekali di luar perulangan karena isinya sama untuk semua kandidat, sedangkan hasil tiap kandidat ditampung array berpraalokasi zeros yang diisi per indeks. Saringan kandidat > max(T) membuang kandidat yang melanggar syarat positif sebelum logaritma dihitung, sehingga penjaga pada callout sebelumnya ditegakkan oleh kode dan bukan oleh kedisiplinan pemanggil. Pemilihan kandidat memakai min yang mengembalikan nilai sekaligus posisinya, sehingga parameter terbaik diambil dari baris yang sama pada tabel kandidat. Pemanggilan berikut membandingkan seluruh kandidat pada data proyek.

matlab Membandingkan seluruh kandidat T_aset

% bandingkan seluruh kandidat T_aset lalu pilih yang terbaik
M = csvread('data-suhu.csv');
t = M(:, 1);
T = M(:, 2);
hasil = fit_termal(t, T);

fprintf('T_aset   tau     RMS\n');
for k = 1:numel(hasil.kandidat)
  fprintf('%5.0f  %6.2f  %7.4f\n', hasil.kandidat(k), ...
          hasil.tau_kandidat(k), hasil.rms_kandidat(k));
end
fprintf('terpilih : T_aset %.2f derajat, T_awal %.2f derajat, tau %.4f detik\n', ...
        hasil.T_aset, hasil.T_awal, hasil.tau);

text Keluaran

T_aset   tau     RMS
   85   30.10   0.1888
   90   50.15   4.3949
   95   65.69   5.3525
  100   79.93   5.8300
terpilih : T_aset 85.00 derajat, T_awal 25.40 derajat, tau 30.0967 detik

Kandidat 85 unggul dengan RMS 0.1888 derajat, sedangkan tiga kandidat lain melenceng dengan RMS di atas 4 derajat. Sebabnya terlihat dari kolom tau: ketika kandidat terlalu tinggi, kurva logaritmanya melengkung sehingga garis lurus hasil regresi tidak mewakili data. Konstanta waktu hasil regresi pun melebar sampai puluhan detik tanpa makna fisis. Data pengukuran sendiri menolak kandidat yang salah dengan RMS yang besar, sehingga pemilihan tidak bergantung pada dugaan awal. Parameter final proyek adalah $T_{\text{aset}} = 85{,}00$ derajat, $T_{\text{awal}} = 25{,}40$ derajat, dan $\tau = 30{,}0967$ detik.

Skrip utama analisis_termal#

Skrip utama merangkai seluruh alur dalam empat langkah bernomor: memuat data, menyesuaikan model, memprediksi waktu ambang dengan dua cara, lalu menyimpan figur. Seluruh angka yang dipakai lebih dari sekali dikumpulkan sebagai konstanta bernama pada bagian atas, yaitu ambang dan berkas_data, sehingga perubahan kebijakan proyek cukup dilakukan pada satu tempat.

matlab analisis_termal.m

% analisis_termal.m: analisis lengkap respon termal alat pemanas
%
% Memuat data pengukuran, menyesuaikan model pemanasan orde satu,
% memeriksa mutu lewat error RMS dan interpolasi, lalu menyimpan
% figur data beserta model. Skrip ini memanggil fit_termal.m dan
% membutuhkan data-suhu.csv pada folder yang sama.

% konstanta proyek dikumpulkan di satu tempat
ambang = 80;                       % derajat
berkas_data = 'data-suhu.csv';

% 1. muat data pengukuran
M = csvread(berkas_data);
t = M(:, 1);
T = M(:, 2);
fprintf('data termuat : %d titik pengukuran\n', numel(t));

% 2. sesuaikan model pemanasan orde satu
hasil = fit_termal(t, T);
fprintf('T_aset : %.2f derajat\n', hasil.T_aset);
fprintf('T_awal : %.2f derajat\n', hasil.T_awal);
fprintf('tau    : %.2f detik\n', hasil.tau);
fprintf('RMS    : %.4f derajat\n', hasil.rms);

% 3. prediksi waktu mencapai ambang, dua cara yang saling memeriksa
t_rumus = hasil.tau * log((hasil.T_aset - hasil.T_awal) / (hasil.T_aset - ambang));
t_interp = interp1(T, t, ambang, 'linear');
fprintf('waktu ke %.0f derajat (rumus)      : %.2f detik\n', ambang, t_rumus);
fprintf('waktu ke %.0f derajat (interpolasi): %.2f detik\n', ambang, t_interp);

% 4. figur data pengukuran dan kurva model
t_halus = linspace(min(t), max(t), 200);
T_halus = hasil.T_aset - (hasil.T_aset - hasil.T_awal) * exp(-t_halus / hasil.tau);
plot(t, T, 'bo', t_halus, T_halus, 'b-', ...
     [min(t) max(t)], [ambang ambang], 'k--');
xlabel('waktu (detik)');
ylabel('suhu (derajat Celsius)');
title('Respon termal alat pemanas: data dan model');
legend('data pengukuran', 'model terfit', ...
       sprintf('ambang %.0f derajat', ambang), 'location', 'southeast');
print('-dsvg', 'respon-termal.svg');

text Keluaran

data termuat : 13 titik pengukuran
T_aset : 85.00 derajat
T_awal : 25.40 derajat
tau    : 30.10 detik
RMS    : 0.1888 derajat
waktu ke 80 derajat (rumus)      : 74.58 detik
waktu ke 80 derajat (interpolasi): 74.80 detik

Kedua kriteria mutu spesifikasi terpenuhi. Error RMS 0.1888 derajat berada jauh di bawah ambang 0,5 derajat. Prediksi waktu mencapai 80 derajat lewat rumus model, yaitu 74.58 detik, berselisih hanya 0.22 detik dengan interpolasi terbalik pada data, yaitu 74.80 detik. Interpolasi terbalik dilakukan dengan menukar kedua sumbu pada interp1 menjadi interp1(T, t, ambang) karena suhu naik searah, sehingga setiap nilai suhu di dalam rentang data dimiliki tepat satu waktu pengukuran. Kesepakatan dua cara yang berbeda asumsi ini menjadi bukti yang lebih kuat daripada satu angka tunggal, karena kekeliruan pada salah satunya akan langsung terlihat sebagai selisih.

Figur hasil skrip memperlihatkan penanda bulat data pengukuran serta kurva model terfit yang mendatar di dekat 85 derajat. Garis ambang 80 derajat memotong kurva pada waktu sekitar 75 detik, sesuai dua prediksi di atas.

Data pengukuran suhu, kurva model terfit, dan garis ambang 80 derajat

Verifikasi vektorisasi#

Error per titik pada proyek ini dihitung dengan satu ekspresi tervektorisasi, yaitu hasil.T_model - T, yang mengurangi seluruh elemen sekaligus. Blok berikut membandingkannya dengan versi perulangan yang mengisi array berpraalokasi elemen demi elemen, lalu mencetak tabel per titik beserta pemeriksaan akhirnya.

matlab Error per titik: tervektorisasi dan perulangan

% muat ulang data dan hasil fit agar blok ini berdiri sendiri
M = csvread('data-suhu.csv');
t = M(:, 1);
T = M(:, 2);
hasil = fit_termal(t, T);

% error per titik dengan operasi tervektorisasi
selisih = hasil.T_model - T;

% pembanding: perulangan dengan praalokasi
selisih_ulang = zeros(size(T));
for k = 1:numel(T)
  selisih_ulang(k) = hasil.T_model(k) - T(k);
end

fprintf('  t   data   model  selisih\n');
for k = 1:numel(t)
  fprintf('%3d  %6.2f  %6.2f  %8.4f\n', t(k), T(k), hasil.T_model(k), selisih(k));
end
fprintf('selisih maksimum kedua cara : %.1e\n', max(abs(selisih - selisih_ulang)));
fprintf('RMS dari selisih            : %.4f\n', sqrt(mean(selisih.^2)));

text Keluaran

  t   data   model  selisih
  0   25.00   25.40    0.4026
 10   42.21   42.25    0.0409
 20   54.05   54.34    0.2861
 30   63.10   63.00   -0.0952
 40   69.25   69.22   -0.0271
 50   73.45   73.68    0.2331
 60   77.00   76.88   -0.1176
 70   79.40   79.18   -0.2228
 80   80.65   80.82    0.1733
 90   82.15   82.00   -0.1459
100   82.75   82.85    0.1010
110   83.60   83.46   -0.1414
120   83.80   83.89    0.0943
selisih maksimum kedua cara : 0.0e+00
RMS dari selisih            : 0.1888

Selisih maksimum kedua cara bernilai 0.0e+00, artinya hasil tervektorisasi identik bit demi bit dengan perulangannya, dan RMS yang dihitung ulang dari selisih bernilai sama dengan bidang hasil.rms. Selisih terbesar, 0.4026 derajat, terjadi pada titik pertama karena $T_{\text{awal}}$ ditaksir regresi dan tidak dipaksa sama dengan pengukuran pertama. Transformasi logaritma menimbang titik-titik akhir yang selisihnya kecil lebih berat daripada titik awal, sehingga titik awal paling longgar. Perilaku semacam ini wajar pada linierisasi logaritma dan menjadi salah satu alasan error RMS dihitung pada skala suhu asli, bukan pada skala logaritmanya.

Kesalahan umum#

Membiarkan array membesar di dalam perulangan. Penyebabnya mengisi selisih_ulang(k) tanpa menyiapkan zeros(size(T)) lebih dahulu. Akibatnya setiap penambahan elemen menyalin ulang seluruh array sehingga perhitungan melambat pada data panjang. Perbaikannya selalu mempraalokasi array hasil sebelum perulangan diisi per indeks.
Menyebar angka ajaib pada banyak baris. Penyebabnya nilai seperti 80 dan 85 diketik ulang di setiap tempat dipakainya. Akibatnya perubahan ambang harus dilakukan di banyak baris dan mudah terlewat satu. Perbaikannya mengumpulkan konstanta bernama pada satu tempat, seperti ambang dan berkas_data di bagian atas skrip utama.
Menerima hasil fit tanpa membandingkan kandidat. Penyebabnya satu kandidat langsung dipakai karena angkanya masuk akal. Akibatnya kandidat yang salah menghasilkan parameter yang tampak meyakinkan padahal RMS-nya buruk. Perbaikannya selalu mencetak RMS seluruh kandidat lalu memilih yang terkecil, seperti tabel pada fungsi penyesuai.

Ringkasan#

Linierisasi $\ln(T_{\text{aset}} - T)$ mengubah model orde satu menjadi regresi linier berkemiringan $-1/\tau$ yang diselesaikan operator backslash pada matriks rancangan $[t\ \ 1]$ .
Nilai $T_{\text{aset}}$ yang belum diketahui ditangani dua tahap: beberapa kandidat tetap diuji dengan regresi yang sama, lalu kandidat dengan error RMS terkecil dipilih.
fit_termal mengembalikan struktur parameter beserta tabel kandidatnya, dan analisis_termal merangkai pemuatan data, penyesuaian, prediksi waktu ambang dua cara, serta figur SVG.
Hasil proyek: $T_{\text{aset}} = 85{,}00$ derajat, $T_{\text{awal}} = 25{,}40$ derajat, $\tau = 30{,}10$ detik, RMS 0.1888 derajat, dan waktu mencapai 80 derajat sebesar 74.58 detik menurut rumus serta 74.80 detik menurut interpolasi.

Referensi#

mldivide, \ — MathWorks — Rujukan operator backslash yang mengerjakan kuadrat terkecil pada sistem yang dilinearkan.
interp1 — MathWorks — Rujukan interpolasi satu dimensi yang dipakai untuk interpolasi terbalik dari suhu ke waktu.
print — MathWorks — Rujukan penyimpanan figur ke berkas gambar, termasuk format vektor SVG.