pkg load symbolic %ТОЛЬКО ДЛЯ OCTAVE
pkg load optim %ТОЛЬКО ДЛЯ OCTAVE
warning('off','all');
close all
format short
%_______________________________БЛОК КОНСТАНТ
t_vh_im = 160;%Температура входа в ИМ, C
t_vyh_im= 170;%Температура выхода из ИМ,C
t_sred_im = (t_vh_im+t_vyh_im)/2; %средняя температура в ИМ, C
l_pg = 15000; %Длина ПГ, мм. Экстраполяция очень сильно выбивается при длине ПГ, большей 14,5 м
C_mk = 160; %Концентрация ЭДТА в МК, г/кг
f_im = 197; %Площадь поверхности одного модуля, м2
F_im = f_im*10; %Площадь поверхности в расчёте на 10 модулей, м2
G_mr=1000/3.6; %Расход МР, кг/с
C_mr = 1;%Концентрация ЭДТА в МР на выходе из ПГ, г/кг
Kp_edta = 1;%Коэффициент распределения ЭДТА
w_korr = 0.0010; %Оценка скорости коррозии по протоколам ХП БАЭС (по скорости удаления отложений, прнимая 5% корр. железа), г/м2*с
%w_korr = 0.0044 %Оценка Смыкова, г/м2*с
function rho = density(t) % Функция, описывающая плотность воды в зав-ти от температуры.
rho = 0.14395/(0.0112^(1+(1-(t+273.15)/649.727)^0.05107))*10^(-3); % кг/л. Справедливо до 374 °C
end
density(t_sred_im)
%_______________________________расчет опытов ФЭИ (в диапазоне от 170 до 210 градусов)
%1)Расчёт на основе протоколов ХП БАЭС
%function p = wfei(t) %уравнение расчета скорости удаления отложений, г/м2*с
%p = 1.0000e-03+1.0000e-04*t; %На основе протоколов ХП от Белоярской АЭС
%end
% 2) старый расчет опытов ФЭИ (в диапазоне от 170 до 210 градусов)
%function p = wfei(t) %уравнение расчета скорости удаления отложений
%Gfei=33.333; % расход ФЭИ, мм/сек
%Dfei=11; % внутренний диаметр образца ФЭИ, мм
%Lfei=1000; % длина образца ФЭИ, мм
%Sfei=pi*Dfei*Lfei; % площаль поверхности образца, мм^2
%sfei=pi*Dfei^2/4; % внутреннее сечение образца, мм^2
%Tfei=Lfei/Gfei; % время контакта, сек
%t_rra_phei = [170 180 190 200 210]; %Температуры по графику ФЭИ (отчёт, ст. 30)
%c_rra_phei = [50 80 120 160 200]; %Концентрация железа в промывочном растворе по графику ФЭИ
%Vfei=sfei*Lfei; % объем раствора в образце, мм^3
%Mfei=Vfei.*c_rra_phei./1000; % масса удаленного железа, г
%wfei_uzl_0=Mfei./Sfei./Tfei; % уравнение расчета скорости удаления отложений в узловых точках
%wfei_uzl = wfei_uzl_0./0.582*1; %Пересчёт на концентрацию свободной ЭДТА в 1 г/л. 0.582 - средняя концентрация ЭДТА в опыте Смыкова
%wfei_uzl(4) = 0.0206;
%p = polyval(polyfit(t_rra_phei,wfei_uzl,2),t);
%end
%3) Новый расчёт опытов ФЭИ (
function p = wfei(t)
F = 0.0345; %Внутренняя площадь поверхности трубки
m_max = [16.7 15.2 18.1 9.85 12.2]; %масса отложений в трубке, г
m_vym = [14.5 13.5 17.2 9.8 12.1]; %вымытая масса отложений, г
t_pr = 3600*[3 6.833333333 5.66666667 6.5 7.25];
m_max = [484.5 440 524 285 348.3].*F; %масса отложений в трубке, г
m_ost = [10 8.5 10.2 5.4 10.1].*F; %остаточная масса отложений в трубке, г
m_vym = m_max-m_ost;
t_pr = 3600*([3 6.833333333 5.66666667 6.5 7.25]+[20 5 5 5 0]./60+8/60);
C_edta_0 = [1.7 1.5 1.45 1.7 1.8]; %Концентрация ЭДТА в начале опыта
C_edta_kon = [0.3 0 0.1 0 0.05]; %Концентрация ЭДТА в конце опыта
C_edta_sr = (C_edta_0+C_edta_kon)/2; %Средняя концентрация за время опыта
V_kont = 50; %объём контура, л
G = 12.5;% Расход в контуре, л/ч
m_vym_vh = [200 200 248 65.5 102]*1.38*V_kont*10^(-3);
t_trub = 30; %Время прохода раствора через трубку, с
t_kont = 300; %Минимальное время прохода раствора по контуру (с выхода трубки на вход трубки), с
C_sv_vh = {[1.6 1.4 1.1 1.3 1.1 0.5];
[1.4 1.3 1 1.1 0.8 0.6 0.37 0.25 0.1 0];
[0.9 0.78 0.53 0.4 0.23 0.2];
[1.26 1.43 0.4 0.37 0.4 0.13 0.02] ;
[1.58 0.85 0.53 1.05 0.99 0.26 0.2 0.1] };
C_poln_vh = {[1.7 1.7 1.6 1.5 1.7 1.54];
[1.6 1.4 1.61 1.71 1.4 1.37 1.46 1.46 1.52 1.28];
[1.5 1.46 1.34 1.31 1.39 1.4 ];
[1.41 1.61 1.66 1.66 1.66 1.28 1.37];
[1.7 1.7 1.71 1.56 1.4 1.88 1.8 1.65] };
C_sv_vyh = {[1.5 1.15 1 0.9 0.4 0.1] ;
[0.23 0.8 0.55 0.4 0.3 0.25 0.1 0.1 0 0];
[0.96 0.58 0.5 0.1 0.1 0] ;
[1.27 0.44 0.23 0.18 0.1 0.03 0] ;
[1.58 0.64 0.38 0.18 0.13 0.03 0 0] };
C_poln_vyh = {[1.7 1.72 1.92 1.6 1.66 1.54] ;
[1.58 1.3 1.24 1.55 1.55 1.43 1.49 1.28 1.3 1.4];
[1.42 1.34 1.26 1.31 1.26 1.26];
[1.44 1.55 1.67 1.66 1.64 1.28 1.35] ;
[1.8 1.8 1.71 1.71 1.64 1.66 1.6 1.62] ; };
tfei = {[0.583333 0.916666 1.416666 1.66666 2.16666 3];
[0.25 0.833333 1.83333 2.833333 3.33333 3.83333 4.33333 5.66666 6.333333 6.8333333];
[0.333333333333333 1.25 2.25 4 5 5.66666666666667] ;
[0.25 1.5 2.5 3.5 4.5 5.5 6.5] ;
[0.4166666 1.4166666 2.4166666 3.4166666 4.4166666 5.8333333 6.4166666 7.25]}; % Массив времён
C_fe_vh = {[15 156 160 182 200 200] ;
[10 50 87 80 125 160 130 100 140 200];
[16 104 109 180 220 248 ] ;
[25 19.5 17.5 20.6 20 45.2 65.5] ;
[38.4 68 101 52.5 45.8 111 92 102] };
C_fe_vyh = {[33 170 170 191 210 210] ;
[120 134 126 140 137 141 144 150 193 195];
[17.5 108 178 190 240 250];
[38.5 102 118 119 103 120 142] ;
[99.5 124 154 136 126 173 174 178]};
for i=1:length(m_max)
deltaC_smyk = (C_fe_vyh{i}-C_fe_vh{i})*1.38*10^(-3);
mfe_rasch = @(x) m_max(i)*(1-exp(log(1-m_vym(i)/m_max(i))/t_pr(i).*(x-t_trub))).*(x>=t_trub);
dmfe_rasch = @(x) -m_max(i)*log(1-m_vym(i)/m_max(i))/t_pr(i)*exp(log(1-m_vym(i)/m_max(i))/t_pr(i).*(x-t_trub));
Cedta_rasch = @(x) C_edta_0(i)*(1-mfe_rasch(x)/ m_max(i));
wfe_rasch = @(x) dmfe_rasch(x)/F./Cedta_rasch(x);
w_model_0(i) = wfe_rasch(t_kont+15);
end
w_model_0_hp = w_model_0.*1./C_edta_0; %Полученная скорость растворения в начальный момент времени
pH_smyk = [9.35 9.25 9.1 9 9.05];
w_model_0_hp_pH = w_model_0_hp.*(1+(9.6-pH_smyk)*polyfit(pH_smyk,w_model_0_hp,1)(2)); %Учёт изменения скорости растворения отложений от pH
w_model_0_hp_sr = mean(w_model_0_hp_pH); %Средняя скорость растворения магнетита при 185 C
t_rra_phei = [170 180 190 200 210]; %Температуры по графику ФЭИ (отчёт, ст. 30)
c_rra_phei = [50 80 120 160 200]; %Концентрация железа в промывочном растворе по графику ФЭИ
p = w_model_0_hp_sr*polyval(polyfit(t_rra_phei,c_rra_phei,2),t)/100; %Приведённая к температуре скорость (100 - значение концентрации при 185 C)
end
wfei(160:1:170)
%_______________________________Расчёт кривой Смыкова
smyk_temp = [150 160 170 175 180 190]; %Температуры с кривой Смыкова
smyk_t = 0:1:10; %Шаг времени с кривой смыкова
smyk = [0 0 0 0 0 0; %Оцифровка данных Смыкова по программе Graph2digit
1.5 3.6 7.8 11.5 16.7 33.9;
2.6 6.5 14.7 21.5 30.5 56.0;
3.5 9.2 20.8 30.3 42.2 70.8;
4.9 11.9 26.6 38.0 51.6 80.8;
5.8 14.7 32.4 44.8 59.6 86.8;
7.2 17.5 37.0 51.1 66.3 91.4;
8.1 19.9 41.8 56.5 71.9 94.1;
9.3 22.5 46.2 61.4 76.6 96.1;
10.3 24.6 50.0 65.3 80.0 97.2;
11.4 26.8 53.4 69.0 83.2 97.9]
xi = linspace (0, 10, 201); %Индекс времени
yi = linspace (150, 190, 201)'; %Индекс температур
%_______________________________Сетка концентраций для обработки данных Смыкова
C_razl = interp2(smyk_temp,smyk_t,smyk,yi,xi,'spline'); %Концентрация разложившейся ЭДТА по Смыкову
C_svob = 100-C_razl; %Свободная концентрация по Смыкову
for j = 1:length(yi) %Индекс температур
for i = 1:(length(xi)-1) %Индекс времени
dC_sv(i,j) = C_svob(i,j)-C_svob(i+1,j);%Изменение C_sv за delta_t
a(i,j) = dC_sv(i,j)/C_svob(i+1,j); %dC_sv/C_sv_sr pp за delta_t
k_itog(i,j) = a(i,j)/(3600*xi(2)); %секундное изменение С_sv
k_0(j) = k_itog(1,j); %Коэффициент в нулевой момент времени
C_rasch_sred(1,j) = 100;
C_rasch_don(1,j) = 100;
C_rasch_sred(i+1,j) = C_rasch_sred(i,j)*(1-k_0(j)*3600*xi(2)); %Расчётная кривая, полученная последовательным умножением коэффициента
C_rasch_don(i+1,j) = 100*(1-k_0(j)*xi(i+1)*3600/2)/(1+k_0(j)*xi(i+1)*3600/2); %Расчётная кривая по В.Е.
end
end
%_____________________Вычисление поправки
C_razl_sred = 100-C_rasch_sred;
C_razl_don = 100-C_rasch_don;
for j = 1:length(yi) %Индекс температур
delta_C(:,j) = C_razl_don(:,j)-C_razl(:,j);
b(:,j) = LinearRegression ([C_razl(:,j)/100,(C_razl(:,j)/100).^2,(C_razl(:,j)/100).^3], delta_C(:,j)/100); %Аппроксимация (регрессия) кубическим полиномом без свободного члена
end
%_______________________________Обработка кривой удельной загрязнённости
delta = 25; %шаг интерполяции данных по удельной загрязнённости, мм
l_zagr_im = [250 1485 2486 3485 4484 5250 6484 7485 8250 8485 9268 9484 10250 11486 12485 13250 14250 14268 14486]; % координаты точек массива удельной загрязнённости
ud_zagr_im = [73.56 99.98 79.04 78.67 116.8 85.17 103.47 85.87 64.84 122.72 68.50 55.04 53.07 147.39 138.20 150.36 77.9 31.53 118.05]; % соответствующий массив удельных загрязнённостей
ud_zagr_im_ispr = [ud_zagr_im(1:10) 138.20 ud_zagr_im(12:19)];% Исправленный массив уд.загр 11 и 19 индекс
l_zagr_im_ispr = [l_zagr_im(1:10) 9018 l_zagr_im(12:18) 15000];% Исправленный массив длин 11 и 19 индекс
ud_zagr_im_kub_delta = interp1(l_zagr_im_ispr,ud_zagr_im_ispr,0:delta:l_pg,'pchip','extrap'); %Интерполированные значения удельной загрязнённости
%Функция interp1(x,y,xi,'method','extrap') возвращает значения интерполирующей функции в точках xi по узловым точкам x и y по заданному методу (например, pchip - интерполяция кубическим полиномом). С флагом 'extrap' функция также возвращает экстраполированные значения вне узловых точек.
%_______________________________Расчёт удаления отложений с течением времени
d_magn = 0.93 %доля магнетита в отложениях
ud_zagr_im_2023 = (ud_zagr_im_ispr +40)*d_magn; %Удельная загрязнённость ПГ на 2023 год при скорости роста отложений 20 г/м2*год
ud_zagr_im_2023_max = max(ud_zagr_im_2023) %Максимальная удельная загрязнённость
ud_zagr_im_2023_min = min(ud_zagr_im_2023) %Минимальная удельная загрязнённость
ud_zagr_im_sred_2023 = polyval(polyfit(l_zagr_im_ispr,ud_zagr_im_2023,0),0) %Средняя (средняя квадратическая) удельная загрязнённость на 2023 год
delta_zagr = 125; %шаг расчёта ХП ИМ ПГ, мм
delta_t = 120; %шаг по времени, с
for i=1:l_pg/delta_zagr %
l_im_i_sr(i) = [0:delta_zagr:(l_pg+delta_zagr)](i)+delta_zagr/2; %Середины расчётных интервалов
t_im_i_sr(i) = t_vh_im+(t_vyh_im-t_vh_im)*l_im_i_sr(i)/l_pg; %Температуры середин расчётных интервалов, принимая линейное изменение температуры
w_rastv_i_sr(i) = wfei(t_im_i_sr(i)); %г/м^2*с
ud_zagr_im_sr_2023(i) = interp1(l_zagr_im_ispr,ud_zagr_im_2023,l_im_i_sr(i),'pchip','extrap'); %Удельная загрязнённость в середине интервала
ud_zagr_im_sr_2023_t = ud_zagr_im_sr_2023(i); %Начальное значение для цикла
j=1;
while(ud_zagr_im_sr_2023_t>0) %Рачёт полной химической очистки (удельная загрязнённость - до нуля)
ud_zagr_im_sr_2023(i,j) = ud_zagr_im_sr_2023_t-delta_t*wfei(t_im_i_sr(i));
ud_zagr_im_sr_2023_t = ud_zagr_im_sr_2023(i,j);
j =j+1;
end
j=1;
end
for j=1:length(ud_zagr_im_sr_2023(1,:)) %Индекс времени
ud_zagr_sredn_2023_t(j) = polyval(polyfit(1:1:(l_pg+delta_zagr)/delta_zagr-1,ud_zagr_im_sr_2023(:,j).',0),0); %Средняя (средняя квадратическая) удельная загрязнённость с течением времени
t_i(j) = j*delta_t; %Массив времён
end
t_otm_poln = length(ud_zagr_im_sr_2023(1,:))*delta_t/3600 %Расчётное время полной химической очистки, ч
t_otm_ost = find(ud_zagr_sredn_2023_t<20,1,'first')*delta_t/3600 %Расчётное время очистки до средней загрязнённости 20 г/м^2, ч
C_razl_pg = interp2(smyk_temp,smyk_t,smyk,t_im_i_sr,([0 t_i]/3600)','spline')/100; %Массив зависимости доли разложенной коцентрации от времени для ПГ
% _______________________________Расчёт количества необходимой моющей композиции для полной очистки
for j=1:length(ud_zagr_im_sr_2023(1,:)) %Индекс времени
for i=1:length(ud_zagr_im_sr_2023(:,1)) %Индекс длины
t_ud_otloj(i) = find(ud_zagr_im_sr_2023(i,:)<0,1,'first')*delta_t; %Полное время удаления отложений на участке delta_zagr
t_korr(i) = t_otm_poln*3600-t_ud_otloj(i); %Полное время коррозии металла участка delta_zagr после удаления отложений с него, с
if (delta_t*j<t_ud_otloj(i)) %Удаление отложений
G_otloj_i(i,j) = w_rastv_i_sr(i)*F_im*delta_zagr/l_pg; %Расход удаляемого железа на отрезке ПГ delta_zagr, г/с
G_mk_otloj(i,j) = (326.3+36)/55.8/C_mk*G_otloj_i(i,j)*density(t_im_i_sr(i)); %Стехиометрический расход МК на растворение окислов на отрезке ПГ delta_zagr, кг/с
G_mk_korr(i,j) = 0;
else %Коррозия
G_korr_i(i,j) = w_korr*F_im*delta_zagr/l_pg; %Расход удаляемого железа на коррозию, г/с
G_mk_korr(i,j) = (326.3+36)/55.8/C_mk*G_korr_i(i,j)*density(t_im_i_sr(i)); %Стехиометрический расход МК на коррозию, кг/с
G_mk_otloj(i,j) = 0;
end
G_otl_kum(j) = sum(G_otloj_i(:))*10^(-3)*delta_t;%Массив накопленных отложений
C_term(i,j) = (C_razl_pg(2,i)-C_razl_pg(1,i))*delta_zagr/l_pg; %Доля разлож. ЭДТА, отнесённая на участок ПГ
C_term_kum(j) = sum(C_term(:))*0.98; %Накопленная доля разлож. ЭДТА (коэфф. 0.98 учитывает удаление продуктов с выпаром)
k_term(i,j) = interp1(yi,b(1,:),t_im_i_sr(i),'*nearest')*C_term_kum(j)+interp1(yi,b(2,:),t_im_i_sr(i),'*nearest')*(C_term_kum(j)).^2+interp1(yi,b(3,:),t_im_i_sr(i),'*nearest')*(C_term_kum(j)).^3;% Поправка на термолиз
C_term_ispr(i,j) = C_term(i,j)*(1-k_term(i,j));
C_term_kum_ispr(j) = sum(C_term_ispr(:))*0.98;
G_mk_term(i,j) = G_mr*C_mr/C_mk*C_term_ispr(i,j)*(1-k_term(i,j));%Расход МК на термолиз реагента, отнесённый на один участок ПГ, кг/с
end
G_vip_0 = 1; %Первое приближение
deltaC = 1;%Первое приближение
while (abs(deltaC))>0.01 %Расчёт выпара и расхода МК методом последовательных приближений
% _______________________________расчет контура ДТПН (НЕ ЗАВЕДЕНО НА РАСЧЁТ СКОРОСТИ КОРРОЗИИ)
G5vh = G_mr; % расход МР, кг/с
h5vh = 719.43; % энтальпия МР после ИМ ПГ, кДж/кг
h5vhv = 675.57; % энтальпия воды при параметрах насыщения, кДж/кг
h5vhp = 2757.4; % энтальпия пара при параметрах насыщения, кДж/кг
h7 = 126.31; % энтальпия ХОВ, кДж/кг
G_mk_vip(j) = G_vip_0*Kp_edta/C_mk*density(t_sred_im); %Расход МК на компенсацию выпара, кг/с
G_mk_i(j) = sum(G_mk_otloj(:,j))+sum(G_mk_korr(:,j))+sum(G_mk_term(:,j))+G_mk_vip(j);%Полный расход МК в момент времени j, кг/с
G9(j)=G_mk_i(j); % заданный расход МК, кг/с
G5vhv(j) = G5vh*(h5vhp - h5vh)/(h5vhp - h5vhv);% расход воды в Д при параметрах насыщения, кг/с
G5vhp(j) = G5vh*(h5vh - h5vhv)/(h5vhp - h5vhv);% расход пара в Д при параметрах насыщения, кг/с
G7(j) = (G5vh*h5vh - G5vh*h5vhv - G9(j)*h5vhp + G9(j)*h5vhv)/(h5vhp-h7);% расход ХОВ, кг/с
G5vihv(j) = (G5vh*h5vh - G5vh*h5vhv - G9(j)*h5vhp + G9(j)*h5vhv)/(h5vhp - h5vhv); % расход ХОВ с конденсатом в Д, кг/с
G5vih(j) = G5vihv(j)+G5vhv(j); % расход МР из Д с учетом выпара, кг/с
Ckonvih(j) = C_mr*(1-C_term(length(ud_zagr_im_sr_2023(:,1)),j)*(1-k_term(length(ud_zagr_im_sr_2023(:,1)),j)))*G5vhp(j)/(G5vhp(j)+G7(j)); % концентрация в конденсате и выпаре, г/кг
G_mk_vip(j) = G_mk_vip(j)*Ckonvih(j);
G6(j) = (G5vh*h5vh - G5vh*h5vhv + G9(j)*h5vhv - G9(j)*h7)/(h5vhp-h7);% расход выпара, кг/с
deltaC = G6(j)-G_vip_0;
G_vip_0 = G6(j); %Расход ЭДТА в выпаре, г/с
C_vh_i(j) = (C_mr*G_mr+C_mk*G_mk_i(j))/(G_mr+G_mk_i(j)); %Концентрация в точке смешения, г/кг
G_mk_term_popr(:,j) = G_mk_term(:,j).*(C_vh_i(j)+C_mr)/2/1.35; %Расход МК на термолиз с поправкой (линейной) на среднюю концентрацию МР в ПГ (отнесено на все участки ПГ)
end
end
%________________________________Расчёт итоговых расходов
G_dob_0 = G_mk_i(1)%Расход добавки в начальный момент времени, кг/с.
G_otloj_sum = sum(G_otloj_i(:))*10^(-3)*delta_t %Суммарное количество удаляемых отложений за время delta_t, кг
G_fe_otloj_sum = G_otloj_sum*167.5/231.5 %Количество железа (после восстановления гидразином - преимущественно двухвалентного)
G_fe_korr_sum = sum(G_korr_i(:))*10^(-3)*delta_t %Количество железа (после восстановления гидразином - преимущественно двухвалентного)
G_mk_otloj_sum = sum(G_mk_otloj(:))*delta_t %Суммарное количество МК на удаление отложений, кг
G_edta_otloj_sum = G_mk_otloj_sum *C_mk*10^(-3) %расход ЭДТА на удаление отложений за время ХП, кг
G_mk_korr_sum = sum(G_mk_korr(:))*delta_t %Суммарное количество МК на удаление отложений, кг
G_mk_term_sum = sum(G_mk_term_popr(:))*delta_t %Суммарный расход ЭДТА на термолиз, кг
G_mk_vip_sum = sum(G_mk_vip)*delta_t %Суммарный расход ЭДТА на выпар, кг
G_mk_sum = G_mk_otloj_sum+G_mk_korr_sum+G_mk_term_sum+G_mk_vip_sum %Общий расход МК на ХП
G_edta_sum = G_mk_sum*C_mk*10^(-3) %Общий расход ЭДТА на ХП, кг
%________________________________Оценка погрешности
delta_pogr = 0.15;% Относительная погрешность определения концентрации железа фотоколориметрическим методом с сульфосалициловой кислотой. Пренебрежём погрешностью с разбавлением пробы.
%_______________________________ГРАФИКИ
figure(1)
t_rra_phei = [170 180 190 200 210]; %Температуры по графику ФЭИ (отчёт, ст. 30)
c_rra_phei = [50 80 120 160 200]; %Концентрация железа в промывочном растворе по графику ФЭИ
plot(t_rra_phei,c_rra_phei,'-o',150:1:210,polyval(polyfit(t_rra_phei,c_rra_phei,2),150:1:210))
title('Обработка результатов ФЭИ','fontsize', 18)
%Функция polyfit(x,y,n) возвращает массив коэффициентов полинома n степени [an a(n-1)... a1 a0], аппроксимирующий массив данных по МНК. Polyval(polyfit(x,y,n),a) возвращает значение полинома в точке a.
axis ([150, 210, 0, 200])
grid on
figure(2)
plot(l_zagr_im_ispr,ud_zagr_im_ispr,'--o',0:delta:l_pg,ud_zagr_im_kub_delta,'--',0:delta/10:l_pg,interp1(0:delta:l_pg,ud_zagr_im_kub_delta,0:delta/10:l_pg,'spline')) %сплайн интерполяция исправленного массива
title('Сплайн интерполяция массива','fontsize', 18)
legend({"Массив данных", "Интерполяция кубическим полиномом","Интерполяция сплайном"},'location','northwest')% Интерполяция сплайном кубического полинома не приводит к уточнению значений между узлами
axis ([0, l_pg, 30, 175])
grid on
xticks(0:500:l_pg)
set(gca,'xminortick','on','xminorgrid','on')
yticks(0:2:200)
figure(3)
plot((1:length(ud_zagr_im_sr_2023(1,:))).*delta_t/3600,ud_zagr_sredn_2023_t,[0 t_otm_poln], [20 20],'--')
title('Изменение средей удельной загрязнённости ПГ с течением времени','fontsize', 18)
grid on
axis ([0, t_otm_poln, 0, ud_zagr_sredn_2023_t(1)+20])
figure(4) %Поверхность уровня для удаления отложений с течением времени
shading interp
surf((1:length(ud_zagr_im_sr_2023(1,:))).*delta_t./3600,l_im_i_sr,ud_zagr_im_sr_2023, "edgecolor", "none")
grid off
colormap('cubehelix')
h = colorbar('southoutside');
set( h, 'XDir', 'reverse' );
set( h, 'xtick', 0:20:200);
title('График поверхности удельной загрязнённости в зав-ти от времени и длины трубки','fontsize', 18)
axis ([0, t_otm_poln, 0, l_pg, 0, ud_zagr_im_2023_max+20])
view(2)
for i = 150:5:190
figure(5)
plot(xi, interp2(smyk_temp,smyk_t,smyk,i,xi,'cubic'))
title('Интерполированные Смыковские кривые в промежуточных точках','fontsize', 18)
grid on
xticks(0:0.5:10)
set(gca,'xminortick','on','xminorgrid','on')
yticks(0:2:100)
hold on
end
for i = 10:10:120
figure(6) %Дельта
plot(C_term_kum,k_term(i,:),'--') %График для кубической регрессии
title('Поправка по скорости термолиза ЭДТА для 160-170','fontsize', 18)
grid on
text(C_term_kum(length(C_term_kum)),k_term(i,length(C_term_kum)),num2str(round(yi(i)*10)/10))
set(gca,'xminortick','on','xminorgrid','on')
hold on
end
figure(7) %Дельта
plot((1:length(ud_zagr_im_sr_2023(1,:))).*delta_t./3600,C_term_kum_ispr,'--') %График для кубической регрессии
title('Доля разложения с течением ХП','fontsize', 18)
grid on
set(gca,'xminortick','on','xminorgrid','on')
hold on
figure(8) %Отложения
plot(160:1:170,wfei(160:1:170),'o--') %График для скорости по ФЭИ в зависимости от температуры
title('Доля разложения с течением ХП','fontsize', 18)
grid on
set(gca,'xminortick','on','xminorgrid','on')
hold on
