ODE15i/s sparse Jacobian pattern trick
3 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Hi,
I am using ODE15i/ODE15s to solve a large-stiff-sparse Differential Algebraic Equation (DAE) system (method of characteristics for a hyperbolic PDE-ODE system) in the order of thousands of equations. Now, without providing the sparsity pattern for the Jacobian ('JPattern',{dFdy, dFdyp}) this can take a very long time to solve. However, deriving the sparsity matrix for this system, particularly for dFdy is tricky. Instead, I have generated a sparsity pattern, based on the numerical Jacobian from the odenumjac function and fed this to odeset. This seems to work well, significantly reducing the solver time (often by an order of magnitude) and giving results consistent with those without providing the sparsity pattern.
Firstly, I was interested in whether this is a sound approach for large systems where the Jacobian sparsity pattern is possibly too difficult to derive analytically?
Secondly, my understanding is that implicit solvers, like ODE15i need to evaluate the Jacobian frequently. Has Matlab missed a trick, or am I missing something? If the above approach is indeed sound, then why does ODE15i/s not generate a sparsity pattern using the method above (if the user has not provided one) and use this to accelerate subsequent calculations?
Thanks
0 commentaires
Réponse acceptée
Torsten
le 18 Juil 2016
Your approach is sound if the sparsity pattern does not change with time.
I guess this possible time-dependence is the reason why MATLAB does not generate it automatically at the beginning of the integration.
Best wishes
Torsten.
0 commentaires
Plus de réponses (1)
Damjan Lasic
le 19 Oct 2017
Wanted to add some remarks even if this is an old question.
Indeed MATLAB has to recompute the entire Jacobian, because one can never be sure when and if some of the components may change. This has to do both with time dependence, but in both the sense that the functions may have some time-dependant terms as well as that the values of some functions may be for example zero at some times, which may give a value of 0 at some particular point in the Jacobian, then the value goes >0 as the function value changes.
Your trick is good and i think it can work wonders in many cases, just be careful when using - consider the above effects. I think a safe(r) approach would be to loop through various randomly generated input time and y values, so you can really be sure that the 0's in the Jacobian matrix are due to actual non-dependence rather than due to current input conditions.
Voir également
Catégories
En savoir plus sur Geometry and Mesh dans Help Center et File Exchange
Produits
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!