High Performance Computing

in Power System

function U=Parareal(F,G,T,u0,N,K);

dT=T/N; TT=0:dT:T;

U{1}(1,:)=u0;

for n=1:N

Go(n+1,:)=G(TT(n),TT(n+1),U{1}(n,:));

U{1}(n+1,:)=Go(n+1,:);

end;

for k=1:K

for n=1:N

Fn(n+1,:)=F(TT(n),TT(n+1),U{k}(n,:));

end;

U{k+1}(1,:)=u0;

for n=1:N

Gn(n+1,:)=G(TT(n),TT(n+1),U{k+1}(n,:));

U{k+1}(n+1,:)=Fn(n+1,:)+Gn(n+1,:)-Go(n+1,:);

end;

Go=Gn;

end;

function[t,u]=ForwardEuler(f,tspan,u0,N); dt=(tspan(2)-tspan(1))/N;

t=(tspan(1):dt:tspan(2))’;

u(1,:)=u0(:);

for n=1:N,

u(n+1,:)=u(n,:)+dt*f(t(n),u(n,:));

end;

function u=SForwardEuler(f,t0,t1,u0,n); [t,u]=ForwardEuler(f,[t0 t1],u0,n); u=u(end,:);

sigma=10;r=28;b=8/3;

f=@(t,x) [sigma*(x(2)-x(1)) r*x(1)-x(2)-x(1)*x(3) x(1)*x(2)-b*x(3)];

MF=10; MG=1;

F=@(t0,t1,u0) SForwardEuler(f,t0,t1,u0,MF);

G=@(t0,t1,u0) SForwardEuler(f,t0,t1,u0,MG);

K=20; u0=[20;5;-5]; T=5; N=500;

U=Parareal(F,G,T,u0,N,K);

[t,u]=ForwardEuler(f,[0 T],u0,MF*N);

TT=0:T/N:T;

for k=1:K

plot3(u(:,1),u(:,2),u(:,3),’-b’...

,U{k}(:,1),U{k}(:,2),U{k}(:,3),’.’);

axis([-20 30 -30 40 -10 60]); view([-13,8]); xlabel(’x’); ylabel(’y’); zlabel(’z’);

grid on

pause

end

Page updated

Google Sites

Report abuse