How to send specific conditions to lsqcurvefit?

2 ビュー (過去 30 日間)
Monirul Hasan
Monirul Hasan 2018 年 11 月 26 日
コメント済み: Monirul Hasan 2018 年 11 月 27 日
Hi, I have a experimental data fitting code given below. In this code from the differential equation (line 30), the 'lsqcurvefit' function (line 12) giving B(1) and B(2) estimation after fitting. From the theory, I know B(1)>B(2). I am just wondering, is there a way to send this condition into the function 'lsqcurvefit'?
Thanks
clear all, close all, clc
HighVacData = xlsread('RawData','Sheet1','A2:B1603');% loading raw data from excel file
Timeall = HighVacData(:,1);
Sdataall = HighVacData(:,2);
Time = 0:1:11244.3353200000;
TimeC = Time.';
%parameter estimation (units per ns)
B0 = [2.96e-02,1.45e-3]; %sequence [k_rs]
lb = [0,0];% lower bound estimation
options = optimset('Algorithm', 'trust-region-reflective');
[B,Rsdnrm,Rsd,ExFlg,OptmInfo,Lmda,Jmat] = lsqcurvefit(@ComGenUQFit,B0,Timeall,Sdataall,lb,[],options);%sending data for fitting
F = ComGenUQFit(B,TimeC);% fitted data saving
figure(1)
plot(Timeall,Sdataall,'bo')%raw data
hold on
plot(TimeC,F,'r','linewidth',2)%fit data
grid on
grid minor
xlabel('Time [ns]');
ylabel('Singlet');
legend('Experimental Data','Fit');
set(gca,'xscale','log')
set(gca,'yscale','log')
function S = ComGenUQFit(B,TimeC)% function for fitting
IC_ns = 11000;
function dotn = Mat (TimeC,n)
dotn(1) = -B(1)*n(1)-B(2)*n(1)^2;% equation for singlet Population
end
[TimeC,dotnsolve] = ode45(@Mat,TimeC,IC);% calling ODE solver
S = dotnsolve(:,1);
end

回答 (1 件)

Bruno Luong
Bruno Luong 2018 年 11 月 26 日
編集済み: Bruno Luong 2018 年 11 月 26 日
Change you variable b to
c = [b(1), b(2)-b(1)]
meaning
b = [c(1), c(1)+c(2)]
Make lsqcurvefit works on c, with the constraint c(2)>=0 by seeting lb = [0; 0].
  3 件のコメント
Bruno Luong
Bruno Luong 2018 年 11 月 27 日
Show you code modified.
I don't read xlsx sheet.
Monirul Hasan
Monirul Hasan 2018 年 11 月 27 日
I runned the code with your added lines (from comment 1), unfortunately it showing an error message. Please let me know how to fix it. Cheers.

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeLeast Squares についてさらに検索

製品


リリース

R2017b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by