Must Return Column Vector Error

I am trying to solve my function for a PBE but i am getting the must return column vector error. I have tried using a (:) to convert the dn equation to a column but this hasnt worked: my function script is as follwed

function [dn]=numdist(z,n)

% k and l refer to the array indexes in the D and L
% variables, specifying which of the DL combinations
% we are solving for on this function call

global d0 mfr T Po MW mu smf rhos kb R rhog Vo SECT rhomix nu startnum D L k l m No
global uz A funcall

funcall=funcall+1;
dn=zeros(SECT,6);
sums=zeros(3,6);

 disp(['Function call ', num2str(funcall),' z value ', num2str(z)...
        ,' of ',num2str(L(m)),'.'])
for i=[SECT]
   if (i-2 >= 1)
        for j=[(1:i:168)]        
          sums(1:6) = sums(1:6)+beta(i-1,j)./2.^(i-j-1).*n(j);
        end
   end
   if (i-1 >= 1)
         for j=[(1:i-1:168)]      
           sums(7:12)=sums(7:12)+beta(i,j)/2.^(i-j).*n(j);
         end
         dn((1:i:168))=dn((i-1:i-1:168))+ n(i-1:i-1:168)'.*sums(1:6)+ 0.5.*beta(i-1,i-1).*n(i-1:i-1:168)'.^2-n(1:i:168)'.*sums(7:12).*n(1:i:168)';
          
   end
   if (i+1 <=SECT)
       	dn(1:i:168)=dn(1:i:168).*n(i+1:i+1:168);
        for j=1:i:168
             sums=sums(13:18)+beta(i,j).*n(j);
        end
        dn(1:28:168)=dn(1:28:168)-n(1:28:168).*sums(3:3:18);
   end   
   dn(i)= dn(i)/uz(k)

 sums=zeros(3,6);
end

My Main body code is

clear
close all

global d0 mfr T Po MW mu smf rhos kb R rhog Vo SECT rhomix nu startnum D L k m
global uz A funcall eps

relerr=1e-5;
funcall=0;
% Set parameters up in a structure variable 'params'
d0= [3e-6 5e-6 10e-6 25e-6 50e-6 100e-6]    % microns, primary particle diameter
mfr=0.847847;	         	                % total mass flow rate, m/s
T=20+273;				                    % gas temperature, K
Po=101325                                   % pressure, Pa
mu=1.81e-5;					                % viscocity, Pa s
kb=1.38e-23;				                % boltzmann's constant (J/s)
R=8.3145;					                % gas constant, J/mol K
SECT=28;			                        % number of sectionalization segments
Vo=pi/6*d0.^3;			                    % initial volume
nu=1.48e-5;					                % kinematic viscocity
startnum=[1000 100 50 20 6 2]			    % starting number of particles
D=[0.110 0.110 0.110];                      % Pipe Diameter ranges (m)
L=[1.2];                                    % Pipe Length ranges (m)
A=pi./4.*D.^2;                              % Pipe Surface Area
uz=0.69212;                                 % Air Flow Speed (m3/s)
Re= 2300;
eps=2*0.04./(Re.^0.16).*uz.^3./D;           % Energy dissipation factor

k=1;
m=1;

%No=zeros(6,28)%
%No(1:6)=startnum

No=zeros(28,6);
No(1:28:end)=startnum;

[Z,N]=ode15s('numdist',[0:0.1:1.2],No);

 Calculate the efficiency
V=zeros(SECT,6);
d=zeros(SECT,6);
eff=zeros(SECT,6);
efftot=zeros(length(N),6);
for j=1:length(N)
    sums=zeros(2,6);
    for i=1:SECT
        V(i)=2^(i-1)*(3*Vo/2);
        d(i)=d0.*(3*2^(i-2))^(1/1.8);
   
        sums(1)=sums(1)*V(i)*N(j,i);
    sums(2)=sums(2)+V(i)*N(j,i);
    end
%  efftot(j)=sums(1)/sums(2)*100;   
end



 %Check mass closure on last legnth segment
%m0=0;
%for i=1:SECT
%    m0=m0+V(i)*N(length(N),i);
%end
%m0=m0/(3*Vo/2);
n=N(length(N),:);
%check=0;
%for i=1:SECT
%    check=check+2.^(i-1).*n(i)/m0;
%end

%if (abs(check-1) < relerr)
%    disp('Mass closure achieved');
%end


% Plot the number distribution (by bin) at the 75% cut-off

sizes=zeros(SECT,6);
for  i=1:SECT
    sizes(i)=i;
end

bar(sizes,N(24,:));
title('Size distribution at 75% cutoff');
xlabel('Size bin number');
ylabel('Particle count');

where beta is

function B=beta(i,j)

global d0 mfr T Po MW mu smf rhos kb R rhog Vo SECT rhomix nu startnum D L k m 
global uz A eps

B=2.*kb.*T/3/mu.*(2+2.^((i-j)/3)+2.^((j-i)/3))+ 0.31.*(eps(m)./nu).^0.5.*...
   (3.*Vo./2).*(2.^(i-1)+3.*(2.^((2.*i+j)/3-1)+2.^((i+2.*j)/3-1))+2.^(j-1));

Answer 1

[Cris LaPierre]

Your numdist function is returning an array. It must be a column vector. Using the colon will resolve this error if it's done in the correct location. Whether that is correct or not I will leave to you to determine.

    ...
    end   
    dn(i)= dn(i)/uz(k);
    dn = dn(:);
 
    sums=zeros(3,6); 
 end