-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathcostfun_RT_Bev_Pem.m
63 lines (38 loc) · 1.35 KB
/
costfun_RT_Bev_Pem.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
%--------------------------------------------------------------------------
%Costfunction used to optimize
%--------------------------------------------------------------------------
function J = costfun_RT_Bev_Pem(pars,data)
epsilon = pars(1);
gamma0 = pars(2);
lambda = pars(3);
s = pars(4);
vpre = pars(5);
% Calculate prediction
ydata = data.ydata;
ypred = analyticSol_RT_Bev_Pem_eval_withPre(data.xdata*24, epsilon,gamma0,lambda,s,vpre);
% Relative error
err = (ypred - ydata)./ydata;
% Ignore 0 points
err(ydata==0) = (ypred(ydata==0) - ydata(ydata==0))/0.5;
% Weigh small data point more than the rest
weights = ones(size(err));% mean(ydata)./(ydata);%1./(ydata).^2;%
weights(abs(ypred - ydata)<0.5) = 0.0001;
weights(end) = 1; % Even if the last datapoint is small it should be accounted for
if ydata(end)<0.3
weights(end) = 2;
if err(end) < -0.5
weights(end) = 2000;
end
end
if any(abs(ypred(ydata>2.5) - ydata(ydata>2.5))./ydata(ydata>2.5)>0.3)
inds1 = find(ydata>2.5);
inds = inds1(find(abs(ypred(ydata>2.5) - ydata(ydata>2.5))./ydata(ydata>2.5)>0.3));
err(inds) = abs(ypred(inds) - ydata(inds));
end
if ypred(end)<0.1
if weights(end)<20
weights(end) = 200;
end
end
J = err'*(weights.*err);
end