Помимо моделей, образованных линейной комбинацией функций, рассмотренных ранее, кривые зарядки-разрядки могут аппроксимироваться также более сложными выражениями, которые приближают экспериментальные данные в смысле метода наименьших квадратов. Одна из функций, реализующих такой подход - lsqcurvefit. Рассмотрим ее применение на примере двойной экспоенциальной модели (double exponential fit) - первой из перечня в статье.
Как и ранее, эксперимантальные данные, подготовленные заранее, считываются из файла:
load V_SoC
Следующих шагом необходимо описать саму модель в виде явной функции
demodel=@(k,x) k(1)+k(2)*x+k(3)*(1-exp(-k(4)*x))+k(5)*(1-exp(-k(6)./(1-x)));
Здесь demodel - имя функции, знак @ показывает, что это - функция, первый аргумент (k) - вектор неизвестных параметров, которые надо определенить, второй (x) - независимая переменная. Работая с умножением/делением - не забывать ставить точки, обеспечивающие выполение поэлементных операций.
Указываются нчачальные значения:
k0=[10 1 1 1 -10 1];
Важно: для функций с достаточно сложным поведением (как здесь, где есть быстрый рост/падение на краях интеревала и плавное изменение в середине), процесс регресии может существенно зависеть от начального приближения (процесс может свалиться в локальный, а не глобальный максимум), поэтому, если кривая приближается плохо, стоит поменять начальные значения, подбирая их так, чтобы исход улучшался (можно контролировать графически, перезапуская счет).
Процедура нелинейного фита:
ind=find(SOC1dn<1);
k = lsqcurvefit(demodel,k0,SOC1dn(ind),Vd1(ind));
Первая из строк выше нужна так как аргумент одной из составляющих используемой функции обращается в бесконечность при х=1 и эту точку надо исключить.
Построение получишихся графиков:
figure
plot(SOC1dn,Vd1,'o','MarkerSize',3)
hold on
plot(SOC1dn(ind),demodel(k,SOC1dn(ind)),'-','LineWidth',1)
l=legend({'Data','Double exponential model'},'Location','SouthEast');
set(l,'FontSize',16,'FontName','Times');
xlim([0 1])
ylim([3 4.2])
xlabel('SoC')
ylabel('{\it U}, V')
set(gca,'FontSize',16,'FontName','Times')
Результат:
