I want to fit experimental data by solving 3 coupled differential equations.The ODEs are rate equations modeling the dynamics of hot carriers in semiconductors after excitation. The function use fot fitting data is x2 and the parameters are ntrmax, trec and tc. I was inspired by comments from Peng Bai when writting the following code but It does not work and I can't figure out why yet.
P[t_]= 4.9560 * 10^18 Exp[-1732.87 t^2];
Nphoton = 2.1102 * 10^17;
tr = 0.2;
data = {{{0., 0., -0.000367822, 0., 0.000435918}, {0.1334, 0.01008,0.0000753379, -0.000200166, -0.000520431},{0.267,0.02302,-0.000766664,-0.000244647,-0.000209062},{0.4004,0.03347, 0.00207, 0.00133,-0.0000533769},{0.5338,0.06225,0.00772, 0.00654, 0.00186},{0.6672, 0.07971, 0.01499, 0.01052,0.00362}, {0.8006, 0.07766, 0.0196, 0.01717, 0.00522},{0.934,0.07577, 0.03668, 0.02662, 0.00975}, {1.0674, 0.07153, 0.03945,0.02604,0.01389},{1.2008,0.0715,0.03943,0.02398,0.01189}, {1.3344, 0.06899,0.03748,0.02269,0.01147},{1.4678,0.06541, 0.03571, 0.02017, 0.01113},{1.6012,0.06087,0.03342, 0.01928, 0.01104},{1.7346, 0.06103, 0.03238, 0.01886, 0.01029}, {1.868, 0.06061, 0.02904, 0.0179, 0.00922},{2.0014,0.05896, 0.02939, 0.01806, 0.00771}, {2.1348, 0.05629, 0.0285,0.01657, 0.00849},{2.2684,0.05387, 0.02817,0.01579,0.00782}, {2.4018, 0.05418, 0.02622, 0.0157, 0.00693}, {2.5352,0.05342,0.02558,0.01506,0.00633},{2.6686, 0.05169, 0.0235,0.01417, 0.00626}, {2.802, 0.04953, 0.02494,0.01348,0.00568},{2.9354,0.04977,0.02392, 0.01312, 0.00528},{3.0688, 0.04788, 0.02347, 0.01241, 0.00573},{3.2022, 0.04753, 0.02294,0.01152,0.00609},{3.3358, 0.04713, 0.02148, 0.01212,0.00488},{3.4692, 0.04553, 0.0213, 0.01136, 0.00504},{3.6026,0.04461, 0.02146, 0.01079,0.00504},{3.736, 0.04348, 0.02108,0.01025, 0.00455},{3.8694, 0.04279,0.01969, 0.01021,0.00419},{4.0028, 0.04337, 0.01864, 0.01021,0.00522},{4.1362,0.04321, 0.01913, 0.01003, 0.00466},{4.2696,0.04168, 0.01842,0.00914, 0.00448},{4.4032, 0.0395, 0.01758, 0.00952,0.00442},{4.5366,0.03814, 0.01689, 0.00843, 0.00399},{4.67,0.0403, 0.01647, 0.00854,0.00391}, {4.8034, 0.03894, 0.01654,0.0083, 0.00364},{4.9368,0.03865, 0.01654, 0.0091,0.00397},{5.0702, 0.03732, 0.0168, 0.00896, 0.00233},{5.2036, 0.03672, 0.01696, 0.00812, 0.00364},{5.337, 0.0359, 0.01592, 0.00825, 0.00344},{5.4706, 0.03541, 0.01621, 0.00823, 0.00322},{5.604, 0.03381, 0.0157, 0.0083, 0.00397},{5.7374, 0.03472, 0.01514, 0.00747, 0.00275},{5.8708, 0.03456, 0.01525,0.00754, 0.00368},{6.0042, 0.03414, 0.01454,0.00705,0.00342},{6.1376, 0.03198, 0.01472, 0.00656, 0.00266},{6.271,0.03209, 0.01421, 0.00705, 0.00326},{6.4044, 0.03245, 0.01355,0.0079,0.00315},{6.538, 0.0312, 0.01386, 0.00696, 0.00295},{6.6714, 0.03103, 0.01377, 0.00623, 0.00213},{6.8048, 0.03118, 0.01304, 0.00643, 0.00264},{6.9382, 0.03016, 0.01135, 0.00718, 0.00262},{7.0716, 0.02962,0.01189, 0.00718, 0.00159},{7.205, 0.02998, 0.01189, 0.00585, 0.00244},{7.3384, 0.02956, 0.01171, 0.00607, 0.00226},{7.472, 0.02885, 0.01246,0.00621, 0.00233},{7.6054, 0.02869, 0.01235, 0.00616,0.00239}, {7.7388, 0.02882, 0.01162, 0.00614, 0.00235},{7.8722, 0.02831,0.01233,0.00583, 0.00255}, {8.0056, 0.02664, 0.01111,0.00634,0.0027}, {8.139, 0.02747, 0.01038, 0.00518,0.00295}, {8.2724, 0.0262, 0.0106, 0.00589,0.00264}, {8.4058,0.02662, 0.01009, 0.00547, 0.00226}, {8.5394,0.0256, 0.00947, 0.00527, 0.00206}, {8.6728, 0.02524, 0.0102,0.00467,0.00206}, {8.8062, 0.02415, 0.00971, 0.00458, 0.00197}, {8.9396,0.02302, 0.00945, 0.00434, 0.0019}, {9.073, 0.02382, 0.00969, 0.005, 0.0021}, {9.2064, 0.02346, 0.00989, 0.00527,0.00175}, {9.3398, 0.02289, 0.00967, 0.00483, 0.00179}, {9.4732,0.02226, 0.01007, 0.00527, 0.00177}, {9.6068, 0.02291, 0.00925, 0.00451, 0.00133}, {9.7402, 0.02264, 0.00985, 0.0052,0.00204}, {9.8722, 0.02362, 0.00909, 0.00405, 0.00104}, {10.007,0.0226, 0.00954, 0.00451, 0.0013}, {10.1404, 0.02137, 0.00907,0.0044, 0.00199}, {10.2738, 0.02186, 0.00958, 0.00451, 0.00195}}}
ti = data[[1, All, 1]];
ci = Nphoton*data[[1, All, 2 ;; 5]];
pfun = ParametricNDSolveValue[{x1'[t] == P[t] - x1[t]/tr,x2'[t] == x1[t]/tr - x2[t] (1 - (x3[t]/ntrmax))/tp - x2[t]/trec, x3'[t] == x2[t] (1 - (x3[t]/ntrmax))/tp - x3[t]/trec, x1[0] == P[0], x2[0] == 0, x3[0] == 0}, {x1, x2, x3}, {t, 0, 20}, {ntrmax, trec, tp}];(*three dependent variables*)
f[ntrmax_?NumericQ, trec_?NumericQ, tp_?NumericQ] := Sum[Total[(ci[[All, i]] - Map[pfun[ntrmax, trec, tp][[i]], ti])^2], {i, 1, 4}] // Quiet;
fit = NMinimize[f[ntrmax, trec, tp], {ntrmax, trec, tp}];
params = fit // Last
f
yourpfun
only has 3 interpolating functions and you are taking parts{i, 1, 4}
which gives a string of"Part 4 of {<<1>>} does not exist"
errors. Changing that 4 to 3 seems to get rid of those errors, but you have another error which may be that you aren't correctly getting a numeric result frompfun
. Can you track down and correct those? You also have extremely large and extremely small coefficients and only a single digit of precision fortr
, all of which may make for very uncertain results. Without trying to do the fit, can you get good results fromf
? $\endgroup$data
has five columns. What do these columns correspond to among the parameterst, x1, x2, x3
? $\endgroup$pfun
tried to simply the problem by normalizing thedata ci = data[[1,All, 2;;5]]
andpfun = ParametricNDSolveValue[{x1'[t] == P[t]/P[0] - x1[t]/tr,x2'[t] == x1[t]/tr - x2[t] (1 - (x3[t]/ntrmax))/tp - x2[t]/trec, x3'[t] == x2[t] (1 - (x3[t]/ntrmax))/tp - x3[t]/trec, x1[0] == 1, x2[0] == 0, x3[0] == 0}, {x1, x2, x3}, {t, 0, 20}, {ntrmax, trec, tp}];
Computed parameters do not converge to experimental data. $\endgroup$