# Finding best parametric function estimation for ODE of first order

3 views (last 30 days)
Lyudmil Yovkov on 18 Jul 2021
Edited: Lyudmil Yovkov on 31 Jul 2021
Hello to every one.
Until now I've solved differential equations numerically under given formulae for their coefficients and given initial / boundary conditions. But recently I've been dealt with the inverse problem, namely: find the best fit to the solution of the ODE regarding some unknown paramateres. The hard part here is due to the fact that we need to determine not a single value but the values of a whole parametric function.
Consider the following initial value problem for ODE of first order: Solving the equation by separation of the variables we obtain the exact solution Let's consider however this problem statement: given the Cauchy's problem whose solution is known over the grid Find the estimated values of the function fitting the solution , i.e. establish that they're eventually located close enough to the curve .
I read a lot about the topic but I didn't find what I searched for. Could you give me some short explanations or guidelines on how to perform this idea in general and later using MATLAB toolboxes?
Kind regards,
Lyudmil Yovkov, PhD

Torsten on 18 Jul 2021
Edited: Torsten on 18 Jul 2021
Search for a book on optimal control of ODEs and DAEs, e.g.
Matthias Gerdts:
Optimal control of ODEs and DAEs
De Gruyter 2011
Finally, you will have to use Matlab's lsqcurvefit or fmincon to determine the coefficients c_i on the grid x_i that minimize
sum_{i=0}^{i=n} (u_i - u_i(c_i))^2
where u_i(c_i) are the values for u obtained from the ODE by using the c_i as coefficients in front of du/dx.
Torsten on 21 Jul 2021
The structure of the program will be the same. In f, you will have to call the pdesolver instead of ode15s, and the computation time will be longer.

Alan Stevens on 18 Jul 2021
First plot the points, ui vs xi to see what sort of curve it might be. If it looks like it could be a polynomial (as in the case of u = x^2 + 1) then lookup help on polyfit.
Lyudmil Yovkov on 18 Jul 2021
@Alan Stevens, I gave quite a simple example in order to illustrate what I want to achieve. The real problem I've been dealing with is more complicated and nonlinear. In its formulation u(x) is a temperature and c(u) is a heat capacity. The heat capacity is usually a polynomial function of the temperature: For my goals it will be enough accurate to fit the function with a quadratic one. I've never before performed this inverse problem and this is why I decided to beg for help.

Lyudmil Yovkov on 28 Jul 2021
Edited: Lyudmil Yovkov on 28 Jul 2021
Hello to everyone,
I'm writing again for an advice. First let me introduce into what I've been working on. Given the Cauchy problem where are known constants and is a known function of time. The functions have meaning of heat capacity and depend on the temperature T. For problem (1) an experimental data set as well is provided.
To solve (1) at known capacities and validate the results against the experimental data is a trivial procedure. I encounter difficulties when solving the inverse problem:
(2) Find the best fit for at known physical parameters and functions I'm searching the capacity of the form , i.e. I'm trying to estimate the values of the best fit. There are some constraints imposed upon the values .
So thank to @Torsten, I successfully solved a modal problem of the type (2) but for test solutions prepared beforehand. Even at volatile constraints upon the values of the estimated parameters I obtain great results. I applied the ideas @Torsten suggested above to my real physics problem (2). It turned out that at small changes in the constraints or even at their full removal from the lsqcurvefit() procedure the results for also change significantly. Can you tell me why? Is the problem itself sensitive to the initial guess and the constraints, for instance?
I'm performing the experiments for pure water whose physical parameters are known. In particular, the function for pure water is known and I'm comparing my results with the analytical formula.
Best regards!
Lyudmil Yovkov on 31 Jul 2021
Can anybody advise me how to overcome the diffuculties when invoking lsqcurvefit() with different parameter constraints and explain to me the reason for the inaccurate numerical results? I'm attaching some plots.
Test 1, test 2, test 3 are performed on a test example with known function I'm comparing with. Here . Fig. 1 Fig. 2 Fig. 3
Here are some of the results for the real physics problem I'm dealing with (Test 1 and test 2 below). Small changes, for instance, in the input constraints lead to significant displacement. The estimated -curve is not close to the exact one. Fig. 4:
lb = [4.1, -2.3e-03, 4.8e-05, -4.8e-07, 1.7e-09],
ub = [4.3, -2.1e-03, 5.0e-05, -4.6e-07, 1.9e-09] Fig. 5:
lb = [3.0, -3.0e-03, 3.5e-05, -5.0e-07, 1.0e-09],
ub = [5.0, -1.5e-03, 5.0e-05, -3.5e-07, 2.5e-09]
Thanks if you find a little time to respond to me and help me find the source of the errors.
Best regards!