D:/demona/2.0/manual/doc/2-31set/matlab/cont.m

Go to the documentation of this file.

clear all
clc
n=65;

global lambda
global unk
global re

lambda=0.1;
re=1000;
unk=2;

N=n*n*unk;




maxNewton=10;
tol=1e-5;

maxit=1000;


% initial guess for the solution vector
x=zeros(N,1);

f0=myF(x,n);

initial_norm=norm(f0)

jac=sparse(N,N);

for i=1:maxNewton
    
    disp('Newton iter')
    i
    
    %initialize the update vector
    delta=zeros(N,1);
        
    Jac=myJac(x,n);
        tol=1e-2; %*norm(f0)
%    delta = gmres(@(v) mfmult(v,x,n),-f0,min(100,N),tol,min(1000,N));

%delta = gmres(Jac,-f0,min(100,N),tol,min(1000,N));
     delta = -Jac\f0;
    
    x=x+delta;
    
    f0=myF(x,n);
      
    if (norm(f0)<=tol)
        
        disp('Newton converged')
        break
    end
    
end

s=zeros(N/unk,1);
v=zeros(N/unk,1);

for i=1:N/unk
    s(i)=x(2*i-1);
    v(i)=x(2*i);
end


S=reshape(s,n,n);
subplot(1,2,1)
contourf(S');
V=reshape(v,n,n);
subplot(1,2,2)
contourf(V');
%disp('maximum temperature:')
% max(x)

%pause

disp('re=100 is solved')


Fp=myFp(x,n);
Jac=myJac(x,n);

delta=-Jac\Fp;

%delta(133:134)

reOld=re;
re=1000;

x=x+delta*(re-reOld);
x(133:134)



f0=myF(x,n);

initial_norm=norm(f0)

for i=1:maxNewton
    
    disp('Newton iter')
    i
    
    %initialize the update vector
    delta=zeros(N,1);
        
    Jac=myJac(x,n);
%    delta = gmres(@(v) mfmult(v,x,n),-f0,min(100,N),tol,min(1000,N));
    tol=1e-2*norm(f0)
delta = gmres(Jac,-f0,min(100,N),tol,min(1000,N));
%     delta = -Jac\f0;
    
    x=x+delta;
    
    f0=myF(x,n);
      
    if (norm(f0)<=tol)
        
        disp('Newton converged')
        break
    end
    
end
figure
s=zeros(N/unk,1);
v=zeros(N/unk,1);

for i=1:N/unk
    s(i)=x(2*i-1);
    v(i)=x(2*i);
end


S=reshape(s,n,n);
subplot(1,2,1)
contourf(S');
V=reshape(v,n,n);
subplot(1,2,2)
contourf(V');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% clear all
% clc
% n=65;
% 
% N=n*n;
% 
% global lambda
% lambda=0.1;
% 
% maxNewton=10;
% tol=1e-5;
% 
% maxit=1000;
% 
% 
% % initial guess for the solution vector
% x=zeros(N,1);
% 
% f0=myF(x,n);
% 
% initial_norm=norm(f0)
% 
% jac=sparse(N,N);
% 
% for i=1:maxNewton
%     
%     disp('Newton iter')
%     i
%     
%     %initialize the update vector
%     delta=zeros(N,1);
%         
%     Jac=myJac(x,n);
% %    delta = gmres(@(v) mfmult(v,x,n),-f0,min(100,N),tol,min(1000,N));
%      delta = -Jac\f0;
%     
%     x=x+delta;
%     
%     f0=myF(x,n);
%       
%     if (norm(f0)<=tol)
%         
%         disp('Newton converged')
%         break
%     end
%     
% end
% 
% X=reshape(x,n,n);
% surf(X');
% disp('maximum temperature:')
% max(x)
% 
% 
% disp('0.1 is solved')
% 
% 
% Fp=myFp(x,n);
% Jac=myJac(x,n);
% 
% delta=-Jac\Fp;
% 
% delta(67)
% 
% lambdaOld=lambda;
% lambda=6.8;
% 
% x=x+delta*(lambda-lambdaOld);
% x(67)
% 
% 
% 
% f0=myF(x,n);
% 
% initial_norm=norm(f0)
% 
% for i=1:maxNewton
%     
%     disp('Newton iter')
%     i
%     
%     %initialize the update vector
%     delta=zeros(N,1);
%         
%     Jac=myJac(x,n);
% %    delta = gmres(@(v) mfmult(v,x,n),-f0,min(100,N),tol,min(1000,N));
%      delta = -Jac\f0;
%     
%     x=x+delta;
%     
%     f0=myF(x,n);
%       
%     if (norm(f0)<=tol)
%         
%         disp('Newton converged')
%         break
%     end
%     
% end
% 
% 
% 
% 

---