代码之家 › 专栏 › 技术社区 › Ief

Matlab曲线拟合不适用于小值(1e-12),我能做什么?

curve-fitting matlab

Ief · 技术社区 · 7 年前

我安装了曲线拟合工具箱,我正在尝试将扩散数据拟合到特定函数。

该函数是以下形式的误差函数:

y = 3500 - 2500 * erf( ( x-x0 ) / ( 2 * sqrt( D * t )) )

我希望应用程序为我提供合理的 D x0

我知道和 x0约为0.0014 但函数本身无法找到这些解决方案。Matlab输出误差“ 输入必须真实且完整。 “尝试使用曲线拟合应用程序的默认参数时。 D 我需要在方程中设置一个预因子,如本例所示( 1e-12

y = erf( ( x-x0 ) / ( 2 * sqrt( 1e-12 * D * t )) )

通过这种方式,matlab找到了正确的解决方案,但仅适用于1e-10到1e-13之间的预因子。

这是一个很大的问题,因为D的正确解决方案将在1e-3和1e-15之间变化,这取决于我将使用的数据集。x0的值也会有所不同。因此,以这种方式,我无法实现一个通用的解决方案。

下面是我正在使用的示例数据集:

y = [6000 6000 6000 6000 6000 6000 6000 6000 6000 5750 5500 5250 5000 4500 4000 3250 2750 2250 1750 1500 1400 1250 1250 1150];
x = [0:0.0001:0.0023];

当在方程中使用以下固定参数时,得到的直线非常适合数据点。但是matlab找不到它们。

D = 7.1e-11;
t = 900;
x0 = 0.0015;

任何帮助都会很棒!非常感谢你。

%% Fit: 'untitled fit 1'.

   % Input data
    C = [6000 6000 6000 6000 6000 6000 6000 6000 6000 5750 5500 5250 5000 4500 4000 3250 2750 2250 1750 1500 1400 1250 1250 1150]';
    x = [0:0.0001:0.0023]';

    [xData, yData] = prepareCurveData( x, C );

    % Set up fittype and options.
    ft = fittype( '3500-2500*erf((x-x0)/(2*sqrt(1e-10*D*900)))', 'independent', 'x', 'dependent', 'y' );
    opts = fitoptions( 'Method', 'NonlinearLeastSquares' );
    opts.Display = 'Off';
    opts.MaxIter = 4000;
    opts.StartPoint = [0.5 0.0014];
    opts.Upper = [1 Inf];

    % Fit model to data.
    [fitresult, gof] = fit( xData, yData, ft, opts )

    % Plot fit with data.
    figure( 'Name', 'untitled fit 1' );
    h = plot( fitresult, xData, yData );
    legend( h, 'C vs. x', 'untitled fit 1', 'Location', 'NorthEast' );
    % Label axes
    xlabel x
    ylabel C
    grid on

(绘图仅为方便起见)

1 回复 | 直到 7 年前

m7913d 7 年前

你的函数不是为每个实值定义的。假设您只对 D x0 ,只需将下限指定为,即可将搜索范围绑定为正数 [0 0]
由于数值很小,Matlab在数值计算函数导数时存在问题。因此,最好的解决方案是象征性地计算函数的雅可比矩阵。这可以使用符号工具箱完成:
```
syms x x0 D t real;

y = 3500 - 2500 * erf( ( x-x0 ) / ( 2 * sqrt( D * t )) )
J = jacobian(y, [D, x0]);
```
matlabFunction 将其转换为普通的matlab函数或将结果复制并粘贴到脚本中,这样就不必在每次迭代中象征性地计算雅可比矩阵。

您应该设置选项 'SpecifyObjectiveGradient',true

ydata = [6000 6000 6000 6000 6000 6000 6000 6000 6000 5750 5500 5250 5000 4500 4000 3250 2750 2250 1750 1500 1400 1250 1250 1150];
xdata = 0:0.0001:0.0023;
t = 900;

D = 0.1;
x0 = 0.1;

options = optimoptions('lsqcurvefit', 'SpecifyObjectiveGradient',true);
X = lsqcurvefit(@(x, xdata) y(x(1), x(2), t, xdata),[D, x0], xdata, ydata, [0 0], [], options);

D = X(1); % 6.4833e-11
x0 = X(2); % 0.0015

figure
hold on
plot(xdata, ydata);
plot(xdata, y(D, x0, t, xdata));

function [F, J] = y(D, x0, t, x)
  F = 3500 - 2500 * erf( ( x-x0 ) / ( 2 * sqrt( D * t )) );
  J =  [(1250.*t.*exp(-(x - x0).^2/(4.*D.*t)).*(x - x0))/(pi^(1/2).*conj((D.*t).^(3/2))); ...
            (2500.*exp(-(x - x0).^2/(4.*D.*t)))/(pi.^(1/2).*conj((D.*t).^(1/2)))]';
end