# Numerical Methods for Partial Differential Equations

4 次查看（过去 30 天）
Hi, I am following this source for numerically solving an own PDE:
Here, the file polarcPDE.m shows:
close all;
%clear all;
%Setting up the domain
r1 = 0;
r2 = 1;
th1 = 0;
th2 = 2*pi;
%Setting the step size
dr = 0.01;
dth = 2*pi/90;
dr2 = dr*dr;
dth2 = dth*dth;
r = r1+dr:dr:r2;
th = th1:dth:th2;
N = size(r,2) - 1;
M = size(th,2) - 1;
u2 = zeros(N+1,M+1);
usol = zeros(N+1,M+1);
for i = 1:M+1
u2(N+1,i) = cos(2*th(i));
end
%Plotting mesh with BCs
[R,TH] = ndgrid(r,th);
[X,Y] = pol2cart(TH,R);
subplot(1,3,1);
mesh(X,Y,u2);
title('Domain with the Boundary condition cos(2*theta)');
for i=1:N+1
u2(i,1) = r(i)*r(i);
u2(i,M+1) = r(i)*r(i);
end
%To start the iteration loop
err = 1000;
%Number of iterations
k=0;
u = u2;
tol = 1e-6;
while err > tol
for i = 2:N
for j = 2:M
rip = (r(i)+r(i+1))/2;
rim = (r(i)+r(i-1))/2;
term1 = (-1/(r(i)*dr2))*(rim*u(i-1,j) + rip*u(i+1,j));
term2 = (-1/(r(i)*r(i)*dth2))*(u(i,j-1)+u(i,j+1));
term3 = (-1/(r(i)*dr2))*(rip+rim) - 2/(r(i)*r(i)*dth2);
u(i,j) = (term1+term2)/term3;
end
end
err = max(max(abs(u-u2)));
u2 = u;
k = k+1;
end
%Printing the exact solution profile
subplot(1,3,2);
mesh(X,Y,u2);
title('Approximate Solution');
for i = 1:N+1
for j = 1:M+1
usol(i,j) = r(i)*r(i)*cos(2*th(j));
end
end
subplot(1,3,3);
mesh(X,Y,usol);
title('Exact solution');
error = max(max(abs(u2-usol)));
I am trying to use this layout to change the original PDE given in the information on top of the code, to fit this PDE: However, there is a problem. I cannot find any place in the code where the respective polar PDE is actually given?
In any case, is it possible to change this code to fit the given PDE?
Some start- initial conditions are as follows:
u(x,0)=cos(x)
u(0,y)=cos(y)
u'(x,0)=-sin(x)
Thanks
##### 4 个评论显示 3更早的评论隐藏 3更早的评论
Sergio Manzetti 2018-3-9

Thanks! Unfortunately, nothing is given there about how to change the functional form to MATLAB codes for radial problems.

### 类别

Help CenterFile Exchange 中查找有关 Eigenvalue Problems 的更多信息

### Community Treasure Hunt

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

Start Hunting!