# Vector output of a function inside a for loop

2 views (last 30 days)

Show older comments

Hello,

So sorry for such a long code. I have a loop and a nested function whose input and output variables use the index of this loop. This index is n in the code.

When I debug, for n=1 u and p are in 1*101 and u are all 1 and p are all 0.

For n=2, u and p are 2*101, but the first rows are zero in u and p.

For n=3, u and p are 3*101, but the first and second rows are all zero in u and p.

So, Matlab forgets about the previous index. ii=101 in the code which is the second dimension of u and p, where first dimension is 3.

See attached please.

global ii

nn=3.;

ii=101;

b=1.;

dx=1./(ii-1.);

dt=0.001;

Nt=20001;

rmass=1.; rmomi=2.;

t=0.;

i_vect=(1:ii); x_vect=(i_vect-1)*dx;

dc=0.001; dp=0.001;

dckept=dc; dpkept=dp;

% Preallocation-----------------------------------------------------------%

...

%-------------------------------------------------------------------------%

%Initial conditions for H, A, HD, AD, cb, Anm1, ADnm1---------------------%

for n=1:nnsay

H(n)=1./nn; HD(n)=0.;

A(n)=0.;

abcd=(n-11.)^2.;

AD(n)=0.05*exp(-abcd);

A(nn)=0.; AD(nn)=0.;

cb(n)=1./nn;

Anm1=0.; ADnm1=0.;

if(n~=1)

Anm1=A(n-1);

ADnm1=AD(n-1);

end

end

%-------------------------------------------------------------------------%

for nt=1:Nt

% At the previous time step:

for n=1:nn

for i=1:ii

x=(i-1.)*dx;

ub(n,i)=(cb(n)-HD(n)*(x-b/2.)-(ADnm1-AD(n))/2.*(x-b/2.)^2)...

/(H(n)+(Anm1-A(n))*(x-b/2.));

pb(n,i)=0.;

end

end

for n=1:nn

Hbar(n)=H(n);

HDbar(n)=HD(n);

Abar(n)=A(n);

ADbar(n)=AD(n);

c(n)=cb(n);

end

flag3=0;

A(nn)=0.;

AD(nn)=0.;

for lmn=1:8

mmm=1;

pend=pb(1,ii); % At the first time step, pend=0;

flag2=10;

while (flag2==10)

m=1;

while (flag3==0)

for n=1:nn

kkk=1;

k = 1;

nm=n;

[this_u,this_p]=uqp(ii,nm, dx,dt,ub,c,HD,b,AD,H,A);

u=[u; this_u];

p=[p; this_p];

st(k)=p(n,ii)-pend;

k=k+1;

if(k==2)

c(n)=c(n)+dc;

end

if(k==3)

den=st(1)-st(2);

c(n)=(c(n)*st(1)-(c(n)-dc)*st(2))/den;

end

if (k==4)

kkk=kkk+1;

if (kkk==2)

flag3=1;

end

end

dc=dckept;

end

end

ss=0.;

for n=1:nn

ss=ss+c(n);

end

rt(m)=ss-1.;

m=m+1;

if (m==2)

pend=pend+dp;

end

if (m==3)

deno=rt(1)-rt(2);

pend=(pend*rt(1)-(pend-dp)*rt(2))/deno;

end

if (m==4)

mmm=mmm+1;

if (mmm==2)

m=1;

else

flag2=11;

end

end

dp=dpkept;

end

end

for n=1:nn

cb(n)=c(n);

for i=1:ii

ub(n,i)=u(n,i);

pb(n,i)=p(n,i);

end

end

t = t + dt;

end

which is the main code. The nested function is:

function [u,p]=uqp(ii,nm,dx,dt,ub,c,HD,b,AD,H,A)

n=nm;

Anm1=0.;

ADnm1=0.;

ub(n,:)

if(n~=1)

Anm1=A(n-1);

ADnm1=AD(n-1);

end

for i=1:ii

x=(i-1.)*dx;

u(n,i)=(c(n)-HD(n)*(x-b/2.)-(ADnm1-AD(n))/2.*(x-b/2.)^2)/(H(n)+(Anm1-A(n))*(x-b/2.));

end

p(n,1)=0.5*(1.-ub(n,1)*u(n,1));

q(n,1)=0.;

for i=2:ii

q(n,i)=q(n,i-1)-dx*(u(n,i-1)-ub(n,i-1))/dt;

p(n,i)=0.5*(1.-ub(n,i)*u(n,i))+q(n,i);

end

return

end

Any idea will highly be appreciated! In attached, main.m is main code which uses secant method double times. uqp.m and HA.m are nested functions.

##### 2 Comments

John BG
on 25 Dec 2016

why do i have the impression that such long string of input variables, perhaps one fell off and by not passing the variable MATLAB has to assume or ignore?

would you please be so kind to supply a chunk of working code with variables. It hasn't got to be actual code, or the actual values in variables, but it's highly appreciated by readers to have some kind of code with variables, you only supplied part of the algorithm.

When asking, a bit too much of information is usually better than a few too less, don't you agree?

### Answers (1)

Image Analyst
on 25 Dec 2016

Edited: Image Analyst
on 25 Dec 2016

You need to pass p in if you want it to retain all values, or else have the function return just the newest row, and use a temporary p called "this_p" and append it to the "master" p like this:

[u, this_p] = uqp(ii, nm, dx, dt, ub, c, HD, b, AD, H, A);

p = [p; this_p]; % Append new row to the bottom of the existing p.

##### 4 Comments

Image Analyst
on 26 Dec 2016

### See Also

### Categories

### Community Treasure Hunt

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

Start Hunting!