Commit 18e0810b authored by conmccoid's avatar conmccoid
Browse files

Extrap: working on AltS examples

parent 856c6234
......@@ -5,6 +5,11 @@ function [u1,u2] = ALGO_extrap_AltS_v1_20210708(F1,F2,u2,N)
% alternating Schwarz method. The functionals F1 and F2 compute the
% solutions on the respective subdomains by taking the solution on the
% previous subdomain for its boundary data.
%
% nb: This implementation is not designed for efficiency, as each call of
% F1 and F2 will contain extraneous calculations. An efficient
% implenentation would not subdivide the algorithm into base components,
% but this makes it easier to try out different problems.
for i=1:N
u1=F1(u2);
......
function [u1,u2] = ALGO_extrap_AltS_v2_20210713(A,b,f,B1,B2,ind1,ind2)
% ALTS performs alternating Schwarz, in a compact implementation
% [u1,u2] = AltS(?)
m=1; N=10;
I=speye(size(A));
R1=I(ind1,:); P1=I(:,ind1); A1=R1*A*P1; b1=R1*A*b;
R2=I(ind2,:); P2=I(:,ind2); A2=R2*A*P2; b2=R2*A*b;
B1=R1*A*B1*P2; B2=R2*A*B2*P1;
for i=1:N
for j=1:m
u1 = A1 \ (f(u1) - B1*u2 - b1);
u2 = A2 \ (f(u2) - B2*u1 - b2);
end
end
\ No newline at end of file
% EX_extrap_1DLaplace_v1_20210713
% 1DLAPLACE solves the 1D Laplace problem on [0,1] using alternating Schwarz
N=64; x=1:N-1; x=x'/N; M=N-1;
L=ones(M,1); L=spdiags([L,-2*L,L],0:2,M,M+2)*N^2;
B=[-L(:,1),sparse(M,M),-L(:,end)]; L=L(:,2:end-1);
u=ones(N+1,1); f=1;
vexact=-x.^2/2 + x/2 + 1;
v=SUB_extrap_SubdomainSolve_v1_20210713(L,B,f,u,'linear');
figure(1)
plot(x,v,'r*-',x,vexact,'k')
% Subdomains: [0,a] and [b,1]
a=0.6; b=0.4;
ind1=x<a; ind2=x>b; l1=sum(ind1); l2=M-sum(ind2);
I=speye(M); P1=I(:,ind1); R1=I(ind1,:); B1=sparse(l1+1,l1+1,1,M,M)*P2;
P2=I(:,ind2); R2=I(ind2,:); B2=sparse(l2,l2,1,M,M)*P1;
F1=@(u2) SUB_extrap_SubdomainSolve_v1_20210713(R1*L*P1,-R1*L*B1,R1*(f-B*u),u2,'linear');
F2=@(u1) SUB_extrap_SubdomainSolve_v1_20210713(R2*L*P2,-R2*L*B2,R2*(f-B*u),u1,'linear');
gamma=-5:0.01:5; G=zeros(size(gamma));
for k=1:length(gamma)
u2=gamma(k)*ones(sum(ind2),1);
[u1,u2] = ALGO_extrap_AltS_v1_20210708(F1,F2,u2,1);
u2=P2*u2; G(k)=u2(l1+1);
end
plot(gamma,G,'k.--')
\ No newline at end of file
function v = SUB_extrap_SubdomainSolve_v1_20210713(A,B,f,u,type)
% SUBDOMAINSOLVE solves problem on a given subdomain
% v=SubdomainSolve(A,B,f,u,type) gives the solution v to the equation
% Av=Bu-f. If type='linear' then this is solved directly. If
% type='nonlinear' then f depends on v and the equation is solved
% iteratively. (add a choice of iterative solvers?)
if strcmp(type,'linear')
m=1; f=@(~) f;
elseif strcmp(type,'nonlinear')
m=20;
end
Bu=B*u; v=zeros(size(u));
for i=1:m
v = A \ (Bu - f(v));
end
\ No newline at end of file
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment