What's going wrong in the vectorization of my ode function?

3 ビュー (過去 30 日間)
Deepa Maheshvare
Deepa Maheshvare 2021 年 2 月 19 日
コメント済み: Deepa Maheshvare 2021 年 2 月 19 日
Hi All,
I'm trying to setting `vectorization` = 'on' in the ode settings and I clearly don't understand what's going wrong
global mat1 mat2
mat1=[
1 -2 1 0 0 0 0 0 0 0;
0 1 -2 1 0 0 0 0 0 0;
0 0 1 -2 1 0 0 0 0 0;
0 0 0 1 -2 1 0 0 0 0;
0 0 0 0 1 -2 1 0 0 0;
0 0 0 0 0 1 -2 1 0 0;
0 0 0 0 0 0 1 -2 1 0;
0 0 0 0 0 0 0 1 -2 1;
];
mat2 = [
1 -1 0 0 0 0 0 0 0 0;
0 1 -1 0 0 0 0 0 0 0;
0 0 1 -1 0 0 0 0 0 0;
0 0 0 1 -1 0 0 0 0 0;
0 0 0 0 1 -1 0 0 0 0;
0 0 0 0 0 1 -1 0 0 0;
0 0 0 0 0 0 1 -1 0 0;
0 0 0 0 0 0 0 1 -1 0;
];
x0 = [1 0 0 0 0 0 0 0 0 0]';
tspan = 0:0.01:5;
f0 = fun(0, x0);
joptions = struct('diffvar', 2, 'vectvars', [], 'thresh', 1e-8, 'fac', []);
J = odenumjac(@fun,{0 x0}, f0, joptions);
sparsity_pattern = sparse(J~=0.);
options = odeset('Stats', 'on', 'JPattern', sparsity_pattern, 'Vectorized', 'on');
ttic = tic();
[t, sol] = ode15s(@(t,x) fun(t,x), tspan , x0, options);
ttoc = toc(ttic)
fprintf('runtime %f seconds ...\n', ttoc)
plot(t, sol)
function f = fun(t,x)
global mat1 mat2
% f = zeros('like', x)
% size(f)
f(1,:) = 0;
f(2:9,:) = mat1*x + mat2*x;
f(10,:) = 2*(x(end-1) - x(end));
% df = [f(1, :); f(2:9, :); f(10, :)];
end
I've tried to vectorize f in fun and also set vectvars=2 to vectorize the jacobian calculation. However, there is a problem again when vectvars = 2 in joptions and/or 'Vectorized', 'on' in options defined for ode15s.
Unable to perform assignment because the size of the left side is 8-by-1 and the size of the right side is 8-by-10.
Error in cse_11_5_20>fun (line 44)
f(2:9,:) = mat1*x + mat2*x;
Error in odenumjac (line 143)
Fdel = feval(F,Fargs_expanded{:});
Error in cse_11_5_20 (line 31)
J = odenumjac(@fun,{0 x0}, f0, joptions);
Suggestions on how to fix this error will be of great help. Could someone please have a look ?

採用された回答

Walter Roberson
Walter Roberson 2021 年 2 月 19 日
Put back the
% f = zeros('like', x);
  1 件のコメント
Deepa Maheshvare
Deepa Maheshvare 2021 年 2 月 19 日
f = zeros(size(x), 'like', x) ;
This helped! thank you . f = zeros('like', x); didn't work for me.

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

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeOrdinary Differential Equations についてさらに検索

製品


リリース

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by