# HW5
## Problem 1
### Code
1. **定義常數**:首先定義題目中的常數,如熱流量(Q1, Q2)、傳熱係數(U1, U2)、流量總和(F)和溫度(T)。
```matlab=
% Define constants
Q1 = 800;
Q2 = 1000;
U1 = 2.5;
U2 = 0.2;
F = 10;
T = 100;
```
2. **猜測初始值**:設定變量的初始猜測值。
```matlab=
% Initial guess
x0 = [150, 200, 150, 200, 1, 1, 1, 1, 1, 1, 1, 1];
```
3. **設定變量的上下界**:定義變量的取值範圍。
```matlab=
% Bounds for variables
lb = [100, 100, 100, 100, 0, 0, 0, 0, 0, 0, 0, 0];
ub = [290, 310, 290, 330, 10, 10, 10, 10, 10, 10, 10, 10];
```
4. **優化選項**:設定優化函數的選項,優化算法選擇Sequential Quadratic Programming(SQP)。
```matlab=
% Options for fmincon
options = optimoptions('fmincon', 'Algorithm', 'sqp');
```
5. **優化過程**:使用 `fmincon` 函數來最小化成本函數(cost_function),並同時考慮約束條件(constraints)。初始猜測值為 x0,變量的上下界為 lb 和 ub。
```matlab=
% Solve the optimization problem
[x_opt, fval] = fmincon(@(x) cost_function(x, Q1, Q2, U1, U2), x0, [], [], [], [], lb, ub, @(x) constraints(x, F, T, Q1, Q2), options);
```
6. **顯示結果**
```matlab=
% Display the results
disp('Optimal values:')
disp(x_opt)
disp('Minimum cost:')
disp(fval)
```
7. **成本函數**:計算給定變量值的成本。變量包括 t1, t2, t3, t4(溫度),並使用常數 Q1, Q2, U1, U2 進行計算。
```matlab=
function C = cost_function(x, Q1, Q2, U1, U2)
t1 = x(1);
t2 = x(2);
t3 = x(3);
t4 = x(4);
term1 = (Q1 / (U1 * (2/3 * sqrt((320-t2)*(300-t1)) + (320-t2)*(300-t1)/6)))^0.6;
term2 = (Q2 / (U2 * (2/3 * sqrt((340-t4)*(300-t3)) + (340-t4)*(300-t3)/6)))^0.6;
C = 1200 * (term1 + term2);
end
```
8. **約束條件函數**:定義了等式和不等式約束條件。等式約束條件確保流量和溫度之間的物理關係,且各個變量之間滿足問題描述的平衡條件。
```matlab=
function [c, ceq] = constraints(x, F, T, Q1, Q2)
f1 = x(5);
f2 = x(6);
f3 = x(7);
f4 = x(8);
f5 = x(9);
f6 = x(10);
f7 = x(11);
f8 = x(12);
t1 = x(1);
t2 = x(2);
t3 = x(3);
t4 = x(4);
% Equality constraints
ceq = [
f1 + f2 - F;
f1 + f6 - f3;
f2 + f5 - f4;
f5 + f7 - f3;
f6 + f8 - f4;
T*f1 + t4*f6 - t1*f3;
T*f2 + t2*f5 - t3*f4;
f3*(t2 - t1) - Q1;
f4*(t4 - t3) - Q2;
];
% Inequality constraints
c = [];
end
```
### Result
```
Best optimal values:
200.0000 280.0000 100.0000 200.0000 0.0000 10.0000 10.0000 10.0000 0.0000 10.0000 10.0000 0.0000
Best minimum cost:
1.9778e+03
```
## Problem 2
### Code
1. **定義常數**:定義目標函數的係數。
> 注意!由於 `linprog` 是用來最小化函數的,我們將目標函數係數取負數,這樣最大化問題就轉化為最小化問題。
```matlab=
% Coefficients of the objective function (note: we need to minimize -f)
f = [-2.84, 0.22, 3.33, -1.09, -9.39, -9.51, 0];
```
2. **定義不等式約束**:定義不等式約束的係數矩陣 A 和上界 b,這些約束的形式為 $Ax≤b$。
```matlab=
% Coefficients for inequality constraints (Ax <= b)
A = [
1.1, 0.9, 0.9, 1, 1.1, 0.9, 0;
0.5, 0.35, 0.25, 0.25, 0.5, 0.35, 0;
0.01, 0.15, 0.15, 0.18, 0.01, 0.15, 0
];
b = [2e5; 1e5; 2e4];
```
3. **定義等式約束**:定義了等式約束的係數矩陣 Aeq 和右側向量 beq。這些約束的形式為 $Aeq⋅x=beq$。
```matlab=
% Coefficients for equality constraints (Aeq*x = beq)
Aeq = [
0.4, 0.06, 0.04, 0.05, -0.6, 0.06, 0;
0, 0.1, 0.01, 0.01, 0, -0.9, 0;
-6857.6, 364, 2032, -1145, -6857.6, 364, 21520
];
beq = [0; 0; 2e7];
```
4. **定義變數的下界**:所有變數都不小於0,即 $x≥0$。
```matlab=
% Variable bounds (x >= 0)
lb = zeros(7, 1);
```
5. **解決優化問題**:使用 linprog 函數來求解線性規劃問題。linprog 函數返回最優解 x_opt 和最優值 fval。
```matlab=
% Options for linprog
options = optimoptions('linprog');
% Solve the linear programming problem
[x_opt, fval] = linprog(f, A, b, Aeq, beq, lb, [], options);
```
6. **顯示結果**:`-fval` 為目標函數的最大值。
```matlab=
% Display the results in the desired format
fprintf('Optimal values:\n');
fprintf('x1* = %.4f\n', x_opt(1));
fprintf('x2* = %.4f\n', x_opt(2));
fprintf('x3* = %.4f\n', x_opt(3));
fprintf('x4* = %.4f\n', x_opt(4));
fprintf('x5* = %.4f\n', x_opt(5));
fprintf('x6* = %.4f\n', x_opt(6));
fprintf('x7* = %.4f\n', x_opt(7));
fprintf('Maximum objective function value: %.4f\n', -fval);
```
### Result
```
Optimal values:
x1* = 109090.9091
x2* = 0.0000
x3* = 0.0000
x4* = 0.0000
x5* = 72727.2727
x6* = -0.0000
x7* = 58867.8608
Maximum objective function value: 992727.2727
```
## Problem 3
### Code
1. **定義常數和單位換算**
- $\Delta{P_c}$ 從 $psig$ 轉換為 lb/ft^2,1 psig = 144 lb/ft^2。
- $cp$ 從 c.p. (Pa·s) 轉換為 lb/(ft·min),其中 1 c.p. = 0.001 Pa·s,並且轉換為 lb/(ft·min)。
```matlab=
% Given parameters
Delta_Pc = 20 * 144; % psig -> lb/ft^2
A = 250; % ft^2
cp = 20 / 1488.164 / 60; % c.p.(Pa*s) -> lb/(ft·min)
M = 75; % lb/min
C = 0.01; % lb/lb
beta = 3.2e-8; % (lb/ft)^2
a = 3.643;
b = 2.68; ```
2. **過濾時間函數**:定義了過濾時間函數 $t_f$,`@(xc)` 表示這是一個匿名函數,以 $x_c$ 為輸入變量。
```matlab=
% Filtration time function
filtration_time = @(xc) beta * (Delta_Pc * A^2 * xc) / (mu * M^2 * C) * exp(-a * xc + b);
```
3. **目標函數**:目的是最大化過濾時間 $t_f$。因為 `fmincon` 是最小化函數,所以我們取過濾時間的相反數來最小化。
```matlab=
% Objective function to minimize the negative of filtration time (for maximization)
objective_function = @(xc) -filtration_time(xc);
```
4. **初始猜測值和邊界**:設定了 $x_c$ 的初始猜測值 $x_0$ 和邊界 lb、ub,即 $x_c$ 範圍在0~1之間。
```matlab=
% Initial guess & bounds for xc
x0 = 0.5;
lb = 0;
ub = 1;
```
5. **優化選項和求解**:使用 `fmincon` 函數來求解優化問題,fmincon 返回最優 $x_c$ 和最大過濾時間的相反數,記得要取相反數來得到實際的最大過濾時間 $t_f$。
```matlab=
% Options for fmincon
options = optimoptions('fmincon');
% Solve the optimization problem
[xc_opt, neg_max_tf] = fmincon(objective_function, x0, [], [], [], [], lb, ub, [], options);
optimal_tf = -neg_max_tf;
```
6. **顯示結果**
```matlab=
% Display the results
fprintf('Optimal solid mass fraction xc* = %.4f\n', xc_opt);
fprintf('Maximum filtration time tf = %.4f minutes\n', optimal_tf);
```
### Result
```
Optimal solid mass fraction xc* = 0.2745
Maximum filtration time tf = 673.3294 minutes
```