issue with vectorization of ode

Hi All,

This is a follow up to my previous question here.

I'm trying run the code below.

global mat1 mat2 k
k = 2;
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]';
y0 = [1 1 1 1 1 1 1 1]';
z0 = [x0; y0];
tspan = 0:0.01:0.1;

f0 = fun(0, z0);
joptions = struct('diffvar', 2, 'vectvars', [], 'thresh', 1e-8, 'fac', []);
J = odenumjac(@fun,{0 z0}, f0, joptions); 
sparsity_pattern = sparse(J~=0.);
options = odeset('Stats', 'on', 'Vectorized', 'on'); %, 'JPattern', sparsity_pattern);

ttic = tic();
[t, sol]  =  ode15s(@(t,z) fun(t,z), tspan , z0, options);
ttoc = toc(ttic)
fprintf('runtime %f seconds ...\n', ttoc)
plot(t, sol)

function dz = fun(t,z)
fprintf('size of z %d %d...\n', size(z))
dz = zeros(size(z), 'like', z);
f(:,:) = fun_f(z(1:10), z(11:end));
g(:,:) = fun_g(z(11:end));
dz(:,:) = [f;g];
size(dz)
end

function f = fun_f(x, y)
global mat1 mat2 k
f = zeros(size(x), 'like', x);
fprintf('size of f %d %d...\n', size(f))
% size(x)
% size(y)
f(1,:) = 0;
f(2:9,:) = mat1*x + mat2*x -k*y;
f(10,:) = 2*(x(end-1) - x(end));
end

function g = fun_g(y)
global k
g = zeros(size(y), 'like', y);
fprintf('size of g %d %d...\n', size(g))
g(:,:) = k*y;
end

The size of z in (18,1) for few iterations and later it changes to (18,18) becasue of which the following error occurs.

size of z 18 1...
size of f 10 1...
size of g 8 1...

ans =

    18     1

size of z 18 18...
size of f 1 10...
Error using  * 
Incorrect dimensions for matrix multiplication. Check that the number of columns in the first matrix matches the number of rows in the second matrix. To perform elementwise
multiplication, use '.*'.

Error in cse_2_20_21>fun_f (line 60)
f(2:9,:) = mat1*x + mat2*x -k*y;

I'm not sure what's going wrong here. Could someone please have a look? The code works fine when when 'Vectorized' is set 'off'.

[Walter Roberson]

f(:,:) = fun_f(z(1:10), z(11:end));
g(:,:) = fun_g(z(11:end));

Should be

f = fun_f(z(1:10,:), z(11:end,:));
g = fun_g(z(11:end,:));

[Deepa Maheshvare]

Thanks a lot!