ODE15i/s sparse Jacobian pattern trick
조회 수: 4 (최근 30일)
이전 댓글 표시
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
채택된 답변
Torsten
2016년 7월 18일
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
추가 답변 (1개)
Damjan Lasic
2017년 10월 19일
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.
참고 항목
카테고리
Help Center 및 File Exchange에서 Geometry and Mesh에 대해 자세히 알아보기
제품
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!