How to adapt my for loop to get it working with parfor? (ODE solver)

1 view (last 30 days)
Sanne
Sanne on 15 Mar 2015
Commented: Sanne on 16 Mar 2015
Hi guys,
I finally got as far to changing my for loop code in the ODE solver into parfor code. However, something is off which keeps me from getting my old for-loop results. Im really stuck. Could someone please be so kind to look at it? I dont know what else to do besides put the skimmed down code here. Please ask for more/less details!
What youre looking at is there are 66+48+48 ODEs per SACn, for which I reserve blocks of 200 ODEs. My thought was the following: per time step each SACn block of 200 ODEs is uncoupled. Only in the NEXT time step they are coupled at the line gGABA /gtemp = ...Codraai(...). 'Codraai' is a lookup matrix to check to which other ODE of another SACn block each ODE is coupled.
For version:
function [dy t] = encisoODE(t,y,Sgt,Syngg,SACnum,Codraai)
ODEs=SACnum*200;
dy = zeros(ODEs,size(y,2));
gGABA=zeros(ODEs,size(y,2));
...........some constants..........
gGABA=gGABA+gGABA_rest;
Syg = interp1q(Sgt,Syngg,t); % INTERP1 = interpolate for noninteger values of t
Sg = repmat(Syg,size(y,2),1); % template voor alle mogelijkheden y voor dezelfde t
for SACn=200*(0:SACnum-1)
i=19+SACn:1:66+SACn;
dy(i-18+100,:) = alpha.*(1-y(i-18+100,:)).*(1./(1+exp(-((y(i,:) -th1)./k1))))-beta.*y(i-18+100,:);
dy(i-18+150,:) = alpha.*(1-y(i-18+150,:)).*(1./(1+exp(-((y(i-18+100,:)-th2)./k2))))-beta.*y(i-18+150,:);
i=1+SACn:1:66+SACn;
gGABA(i,:)=gGABA_rest+(gGABA_bound-gGABA_rest).*...
(y(Codraai(i,1)-18+150,:)+y(Codraai(i,2)-18+150,:)+y(Codraai(i,3)-18+150,:)+y(Codraai(i,4)-18+150,:)+...
y(Codraai(i,5)-18+150,:)+y(Codraai(i,6)-18+150,:)+y(Codraai(i,7)-18+150,:)+y(Codraai(i,8)-18+150,:));
i=1+SACn:3:18+SACn;
dy(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i)'.*(Egluta-y(i,:))+gGABA(i,:).*(EGABAp-y(i,:))+gK.*(EK-y(i,:))+delta.*((y(i+1,:)-y(i,:))));
i=2+SACn:3:18+SACn;
dy(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i)'.*(Egluta-y(i,:))+gGABA(i,:).*(EGABAp-y(i,:))+gK.*(EK-y(i,:))+delta.*((y(i-1,:)-y(i,:))+(y(i+1,:)-y(i,:))));
i=3+SACn:3:18+SACn;
dy(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i)'.*(Egluta-y(i,:))+gGABA(i,:).*(EGABAp-y(i,:))+gK.*(EK-y(i,:))+delta.*((y(i-1,:)-y(i,:))+(y(i+17,:)-y(i,:))));
i=20+SACn:3:36+SACn;
dy(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i)'.*(Egluta-y(i,:))+gGABA(i,:).*(EGABAd-y(i,:))+gK.*(EK-y(i,:))+delta.*((y(i-17,:)-y(i,:))+(y(i+19+(i-(SACn+20))/3*2,:)-y(i,:))));
i=39+SACn:5:66+SACn;
dy(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i)'.*(Egluta-y(i,:))+gGABA(i,:).*(EGABAd-y(i,:))+gK.*(EK-y(i,:))+delta.*((y(i-19-(i-(SACn+39))/5*2,:)-y(i,:))));
end
end
parfor version:
function [dy t ggabavt] = encisoODE(t,y,Sgt,Syngg,SACnum,L,Co,Codraai)
ODEs=SACnum*200;
dy = zeros(ODEs,size(y,2));
gGABA=zeros(ODEs,size(y,2));
..........some constants...............
gGABA=gGABA+gGABA_rest;
Syg = interp1q(Sgt,Syngg,t); % INTERP1 = interpolate for noninteger values of t
Sg = repmat(Syg,size(y,2),1); % template voor alle mogelijkheden y voor dezelfde t
%parfor temporaries
dsgt1=zeros(SACnum,48,size(y,2));
dsgt2=zeros(SACnum,48,size(y,2));
sgt1=zeros(SACnum,48,size(y,2));
sgt2=zeros(SACnum,48,size(y,2));
gt=zeros(SACnum,66,size(y,2));
dyt3=zeros(SACnum,66,size(y,2));
yt3=zeros(SACnum,66,size(y,2));
parfor SACn=0:SACnum-1
dsgtemp1=zeros(48,size(y,2));
dsgtemp2=zeros(48,size(y,2));
sgtemp1=zeros(48,size(y,2));
sgtemp2=zeros(48,size(y,2));
gtemp=zeros(66,size(y,2));
dytemp3=zeros(66,size(y,2));
ytemp3=zeros(66,size(y,2));
i=19:1:66;
dsgtemp1(i-18,:) = alpha.*(1-sgtemp1(i-18,:)).*(1./(1+exp(-((ytemp3(i,:) -th1)./k1))))-beta.*sgtemp1(i-18,:);
dsgtemp2(i-18,:) = alpha.*(1-sgtemp2(i-18,:)).*(1./(1+exp(-((sgtemp1(i-18,:)-th2)./k2))))-beta.*sgtemp2(i-18,:);
dsgt1(SACn+1,:,:)=dsgtemp1;
dsgt2(SACn+1,:,:)=dsgtemp2;
sgt1(SACn+1,:,:)=sgtemp1;
sgt2(SACn+1,:,:)=sgtemp2;
i=1:1:66;
gtemp(i,:)=gGABA_rest+(gGABA_bound-gGABA_rest).*...
(y(Codraai(i+200*SACn,1)-18+150,:)+y(Codraai(i+200*SACn,2)-18+150,:)+y(Codraai(i+200*SACn,3)-18+150,:)+y(Codraai(i+200*SACn,4)-18+150,:)+...
y(Codraai(i+200*SACn,5)-18+150,:)+y(Codraai(i+200*SACn,6)-18+150,:)+y(Codraai(i+200*SACn,7)-18+150,:)+y(Codraai(i+200*SACn,8)-18+150,:));
gt(SACn+1,:,:)=gtemp;
i=1:3:18;
dytemp3(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i+200*SACn)'.*(Egluta-ytemp3(i,:))+gtemp(i,:).*(EGABAp-ytemp3(i,:))+gK.*(EK-ytemp3(i,:))+delta.*((ytemp3(i+1,:)-ytemp3(i,:))));
i=2:3:18;
dytemp3(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i+200*SACn)'.*(Egluta-ytemp3(i,:))+gtemp(i,:).*(EGABAp-ytemp3(i,:))+gK.*(EK-ytemp3(i,:))+delta.*((ytemp3(i-1,:)-ytemp3(i,:))+(ytemp3(i+1,:)-ytemp3(i,:))));
i=3:3:18;
dytemp3(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i+200*SACn)'.*(Egluta-ytemp3(i,:))+gtemp(i,:).*(EGABAp-ytemp3(i,:))+gK.*(EK-ytemp3(i,:))+delta.*((ytemp3(i-1,:)-ytemp3(i,:))+(ytemp3(i+17,:)-ytemp3(i,:))));
i=20:3:36;
dytemp3(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i+200*SACn)'.*(Egluta-ytemp3(i,:))+gtemp(i,:).*(EGABAd-ytemp3(i,:))+gK.*(EK-ytemp3(i,:))+delta.*((ytemp3(i-17,:)-ytemp3(i,:))+(ytemp3(i+19+(i-(+20))/3*2,:)-ytemp3(i,:))));
i=39:5:66;
dytemp3(i,:) = 1./tau.*(ggluta_bound.*Sg(:,i+200*SACn)'.*(Egluta-ytemp3(i,:))+gtemp(i,:).*(EGABAd-ytemp3(i,:))+gK.*(EK-ytemp3(i,:))+delta.*((ytemp3(i-19-(i-(+39))/5*2,:)-ytemp3(i,:))));
dyt3(SACn+1,:,:)=dytemp3;
yt3(SACn+1,:,:)=ytemp3;
end
%unpicking all the temps to fill dy
for SACn=0:SACnum-1
dy(SACn*200+101:SACn*200+148,:)=dsgt1(SACn+1,:,:);
dy(SACn*200+151:SACn*200+198,:)=dsgt2(SACn+1,:,:);
y(SACn*200+101:SACn*200+148,:)=sgt1(SACn+1,:,:);
y(SACn*200+151:SACn*200+198,:)=sgt2(SACn+1,:,:);
gGABA(SACn*200+1:SACn*200+66,:)=gt(SACn+1,:,:);
dy(SACn*200+1:SACn*200+66,:)=dyt3(SACn+1,:,:);
y(SACn*200+1:SACn*200+66,:)=yt3(SACn+1,:,:);
end
end
  1 Comment
Sanne
Sanne on 16 Mar 2015
Ok apparantly this is very difficult to do since the parallelization must concentrate on one time step.. does anybody know if this is possible?

Sign in to comment.

Sign in to answer this question.

Answers (0)

Categories

Find more on Loops and Conditional Statements in Help Center and File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!