一种非线性函数的曲线拟合方法(函数公式: k = A*(T^a)*exp(E/T) )


上一篇文章说了,函数的曲线拟合我以前没做过,所以是摸着石头过河,不知道所采用的方法是否合理,虽然是完成了拟合,不过我觉得自己采用的拟合方法还是比较原始的,希望做曲线拟合的朋友多多指教。

原始数据如下:
T(K) K
200.00 2.5069E-13
220.00 3.5043E-13
223.00 3.6741E-13
225.00 3.7904E-13
250.00 5.4617E-13
275.00 7.5744E-13
295.00 9.6192E-13
298.00 9.9551E-13
300.00 1.0183E-12
325.00 1.3346E-12
350.00 1.7119E-12
375.00 2.1564E-12
400.00 2.6739E-12
425.00 3.2706E-12
450.00 3.9527E-12
475.00 4.7261E-12
480.00 4.8922E-12
500.00 5.5968E-12
525.00 6.5710E-12
550.00 7.6544E-12
575.00 8.8529E-12
600.00 1.0172E-11
800.00 2.5705E-11
1000.00 5.1733E-11
1250.00 1.0165E-10

目标:拟合成 k = A*(T^a)*exp(E/T) 模式的公式,
其中A、a和E为未知常数,是我们需要通过曲线拟合要求出的数据。

拟合目标中的公式是幂逼近和指数逼近的混合,用Matlab的cftool 工具箱的自定义函数来逼近,效果并不理想,所以我就参考了网上的一些博客和百度知道等资源,采取如下策略:

首先将非线性的拟合公式转化为线性公式,再用求解线性方程组的矩阵方法求出未知常数的值。

具体地说,拟合公式的线性化表达式为: log(k) = log(A) + a*log(T) + E/T 。这里有三个未知常数log(A)、a 和 E,则依次取T,K各三个数据,组成 N 个线性方程组: Cx=b,
其中:x=[log(A), a, E], C=[1, log(T), 1/T], b=log(k) 。
解这些线性方程组,得到所有方程组的解组成的解矩阵 xMat,其大小为 N*3,对解矩阵的每一列求平均,即可得到所求的未知常数值。

根据以上策略,可求得未知常数A、a和E的值如下:

A = 3.8858e-020,a = 3.0595,E = -117.2915

程序源码:

function [A,a,E]= fun_NLFit(T,K)
% 函数 FUN_NLFIT() 根据输入T,K的数据集,求出拟合公式 k = A*(T^a)*exp(E/T)
% 的未知常数 A,a,E 。

logT=log(T);
logK=log(K);
daoT=T.^(-1);
lenT=length(T);
C=ones(3);
xMat=[];
% 为了提高拟合精度,从第一个数据点开始,依次分别取T、K的三个相邻的数据点
% 组成线性方程组,若 T 有 lenT 个元素,则可组成 lenT-2 个方程组
for r=1:lenT-2
C(:,2)=logT(r:r+2);
C(:,3)=daoT(r:r+2);
b=logK(r:r+2);
% C=[1 log(T) 1/T], b=log(k)
x=(C/b)';
xMat=[xMat; x];
% 每解一次方程组,则将解 x 存入解矩阵 xMat
end
% 对解矩阵的每一列求平均,即可得到所求的未知常数值
logA=mean(xMat(:,1));
A=exp(logA);
a=mean(xMat(:,2));
E=mean(xMat(:,3));
% 画出由点集T、K构成的目标曲线
h1=stem(T,K,'bo'); % ‘bo’表示每个点用一个小圆圈表示
set(h1,'MarkerFaceColor','green'); % 小圆圈内的颜色为绿色
set(h1,'LineStyle','none'); % 隐藏基线到点的连线
set(get(h1,'BaseLine'),'LineStyle','none'); % 隐藏基线
hold on; % 保持由点集构成的目标曲线,以便和拟合曲线进行对比

% 根据拟合公式,求出若干的拟合点,画出拟合曲线
t=200:10:1300;
k=A*(t^a)*exp(E/t);
plot(t,k,'r');
% 拟合曲线用红色表示
xlabel('T'); ylabel('K'); title('Nonlinear Curve Fitting');

拟合效果图如下:


  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
实现指数函数曲线拟合需要使用非线性最小二乘法,以下是一个简单的实现: ```c #include <stdio.h> #include <stdlib.h> #include <math.h> #define N 10 // 数据点数目 #define MAX_ITER 100 // 最大迭代次数 #define TOL 1e-6 // 迭代精度 // 指数函数模型 y = a * exp(-b * x) + k double model_func(double x, double a, double b, double k) { return a * exp(-b * x) + k; } // 计算残差 double residuals(double x[], double y[], double a, double b, double k) { double res = 0.0; for (int i = 0; i < N; i++) { res += pow(y[i] - model_func(x[i], a, b, k), 2); } return res; } int main() { double x[N] = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0}; // 数据点的 x 坐标 double y[N] = {2.0, 1.8, 1.5, 1.2, 1.0, 0.8, 0.6, 0.5, 0.4, 0.3}; // 数据点的 y 坐标 double a = 1.0, b = 1.0, k = 1.0; // 初始参数值 double alpha = 0.001; // 步长 double J_curr, J_prev = residuals(x, y, a, b, k); // 当前和上一次的残差 for (int iter = 0; iter < MAX_ITER; iter++) { // 计算梯度 double grad_a = 0.0, grad_b = 0.0, grad_k = 0.0; for (int i = 0; i < N; i++) { double exp_bx = exp(-b * x[i]); grad_a += -2 * (y[i] - a * exp_bx - k) * exp_bx; grad_b += 2 * a * x[i] * (y[i] - a * exp_bx - k) * exp_bx; grad_k += -2 * (y[i] - a * exp_bx - k); } // 更新参数 double a_new = a - alpha * grad_a; double b_new = b - alpha * grad_b; double k_new = k - alpha * grad_k; // 计算新的残差 J_curr = residuals(x, y, a_new, b_new, k_new); // 判断是否收敛 if (fabs(J_curr - J_prev) < TOL) { printf("Converged after %d iterations\n", iter + 1); break; } // 更新参数和残差 a = a_new; b = b_new; k = k_new; J_prev = J_curr; } // 输出拟合结果 printf("a = %f, b = %f, k = %f\n", a, b, k); return 0; } ``` 上述代码中,我们使用了梯度下降法来求解最小二乘问题。首先,我们定义了指数函数模型 `model_func` 和残差函数 `residuals`。然后,我们初始化了参数值和步长,并进行迭代更新。在每次迭代中,我们计算梯度并更新参数,然后计算新的残差。如果新的残差与上一次的残差之差小于设定的迭代精度,我们认为已经收敛,结束迭代并输出拟合结果。 注意:这只是一个简单的实现,实际应用中还需要考虑更多因素,例如初始化参数值的选择、梯度下降的收敛性等。

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值