Метод наименьших квадратов (приближение параболой)
Про метод наименьших квадратов, позволяющий приблизить экспериментальные данные прямой было тут:
http://akostina76.ucoz.ru/blog/2017-09-07-4399
Проблема в том, что данные по теплоотдаче тут:
http://akostina76.ucoz.ru/blog/2017-10-27-4518
… явно не прямая. Теплоотдача сильно растёт с ростом разности температур.
Работать с информаций в виде экспериментальных точек неудобно. Можно конечно, руками нарисовать, что реально это как-то так:
… но математический расчёт явно сделает это точнее.
Прямая задаются функцией y=a*x+b, где константы a и b задают наклон и смещение по оси y/
Следующая по сложности кривая содержит квадрат X (добавляющий искривление функции) и задаётся таким уравнением y = a*x*x+b*x+c.
Как и в случае приближения прямой надо найти те a,b,c при которых кривая (в данном случае парабола) минимально отклоняется от экспериментальных точек:
Отклоняться она может и в плюс ив минус. Чтобы не думать над этим вопросом можно взять квадрат разности. Он-то точно будет положительным, а минимум суммы этих квадратов менее минимальным не станет. Надо только довольно утомительно возвести в квадрат эту разность (временно убрав индекс i для краткости записи):
Полученное выражение это функция F от a,b,c. Все X и Y тут известны. Это экспериментальные данные. А вот минимум этой функции сразу по a,b и с надо найти. Минимум в данном случае это те точки в которых равна нулю первая производная.
Получается система уравнений:
Двойки можно посокращать а слагаемые расставить так чтобы получилась система уравнений в привычном виде:
Осталось засунуть эту систему в Malpe чтобы он её решил. Функция приближения параболой и пример использования (на данных по теплоотдаче):
MNK_x2:=proc(X,Y,x1,x2,n) local N,a11,a12,a13,a21,a22,a23,a31,a32,a33,b1,b2,b3,f1,f2,f3,sys_ur,v,a_r1,b_r1,c_r1,ret;
N:=nops(X):
a11:=sum(X[k]*X[k]*X[k]*X[k],k=1..N):
a12:=sum(X[k]*X[k]*X[k],k=1..N):
a13:=sum(X[k]*X[k],k=1..N):
a21:=a12:
a22:=a13:
a23:=sum(X[k],k=1..N):
a31:=a13:
a32:=a23:
a33:=1:
b1:=sum(X[k]*X[k]*Y[k],k=1..N):
b2:=sum(X[k]*Y[k],k=1..N):
b3:=sum(Y[k],k=1..N):
f1:=a11*a+a12*b+a13*c=b1:
f2:=a21*a+a22*b+a23*c=b2:
f3:=a31*a+a32*b+a33*c=b3:
sys_ur:={f1,f2,f3}:
v:=solve(sys_ur):
subs( v, sys_ur ):
a_r1:=subs(v,a):
b_r1:=subs(v,b):
c_r1:=subs(v,c):
ret:=[a_r1,b_r1,c_r1];
end:
X:=[11.2,11,10.5,10.4,10.1,9.4,9.2,9.1,9,8.9,8.7,8.4,8.3,7.9,7.8,7.7,6.9,6.5,6.1,5.7,5.3,4.4,3.6,3.5]:
Y:=[1.2,1.6,2.4,1.6,1.2,1.2,0.8,1.2,0.8,0.933,1,0.8,0.8,0.8,0.8,0.8,0.4,0.4,0.2,0.3428,0.3,0.4,0.2,0.2]:
abc:=MNK_x2(X,Y):
with(plots):
s_XY:=[seq( [X[i],Y[i]], i=1..nops(X) )]:
g1:=pointplot(s_XY,symbol=CIRCLE,color=BLUE):
g2:=plot(abc[1]*t*t+abc[2]*t+abc[3],t=0..15,color=BLUE,thickness=2):
X2:=[10.4,10.1,9.1,8.7,8.3,7.9,7.8,7.7,6.9,6.5,5.7]:
Y2:=[3,4,2,2.4,1.8,1.2,1.2,1.2,1.2,0.4,0.4]:
s_XY2:=[seq( [X2[i],Y2[i]], i=1..nops(X2) )]:
g3:=pointplot(s_XY2,symbol=BOX,color=RED):
abc2:=MNK_x2(X2,Y2):
g4:=plot(abc2[1]*t*t+abc2[2]*t+abc2[3],t=0..15,color=red):
display(g1,g2,g3,g4);
Получилась вот такая картинка:
Выхода в отрицательные значения конечно быть не может. Но метод выдал то что он мог выдать по засунутым в него данным. Не было в него засунуто данных по малым отклонениям температуры, потому он и «нафантазировал» в этом диапазоне то чего не бывает.
Более менее реальной выглядит первая зависимость (синяя). Это предположительно обычный случай полностью прогретого помещения. А построенный по экспериментальным данным график (точнее параметры параболы a,b,c) позволяет хоть как-то прикинуть что будет происходить при больших отклонениях.
Если при разнице температур в 10 градусов остывание со скоростью 1.31 градус в час:
то при разнице в 20 градусов 5.9 градусов в час:
Как изменится картины при отоплении помещения человеком или курями (т.е при какой наружной температуре можно протопить вечером и не топить всю ночь) вопрос пока отрытый.
|
