حل معادله پواسون دو بعدی به روش ماتریس

حل معادله پواسون دو بعدی به روش ماتریسی

در این مورد با استفاده از روابط زیر

و جایگذاری آنها در معادله پواسون دوبعدی به معادله ماتریسی مشابه رابطه زیر میرسیم:

مهمترین مشکل تشکیل ماتریس ضرایب در حالت دو بعدیست که به نسبت حالت یک بعدی کار مشکلییست. اسکریپت زیر معادله پواسون را به روش ماتریسی برای مسئله ای با همان شرایط اولیه که در بالا انرا به روش SOR حل کردیم، حل می نماید.

ایراد روش ماتریس بخصوص در مواردی که با ابعاد بیش از یک بعد سرو کار داریم بحث حافظه مورد نیاز است. در چنین مواردی بهتر است از روشهایی با دقت بالاتر استفاده نماییم تا مجبور نیاشیم ابعاد ماتریسها را خیلی بزرگ نماییم. در مورد معادلاتی نظیر معادله پواسون در ابعاد بیش از یک بعد روش المان محدود جایگزین مناسبی برای روش تفاضل محدود است چرا که  از دقت بالاتری برخوردار است. ایراد دیگر روش تفاضل متناهی در مواجه با مواردی است که با نواحی مرزی با شکلهای پیچیده تر روبرو هستیم که در این موارد این روش عملا کارایی ندارد. در چنین مواردی روش المان را باید بکار برد. در بخشهای بعدی روش استفاده از جعبه ابزار (toolbox) المان محدود را برای حل چنین مسائلی شرح خواهم داد.

clc

clear; % clear workspace to start new

close; % close previous figures to start new

%——————

Nx=50;       % number of mesh- x dirextion

Ny=50;       % number of mesh- y dirextion

nx=Nx-2;     % number of mesh- x dirextion without boundry point clc

clear; % clear workspace to start new

close; % close previous figures to start new

%——————

Nx=50;       % number of mesh- x dirextion

Ny=50;       % number of mesh- y dirextion

nx=Nx-2;     % number of mesh- x dirextion without boundry point x=0 & x=Lx

ny=Ny-2;     % number of mesh- y dirextion without boundry point y=0 & x=Ly

x=linspace(0,1,Nx);

y=linspace(0,1,Ny);

dx=x(2)-x(1);

dy=y(2)-y(1);

X=repmat(x(2:end-1),1,ny)’;

Y=repmat(y(2:end-1),nx,1);

Y=Y(:);

rho=X;

ep=1;

a=diag(ones(1,nx));

b=diag(ones(1,nx-1),1);

c=diag(ones(1,nx-1),-1);

A1=(-2*a+b+c)/dx^2-2*a/dy^2;

A=A1;

for i=1:ny-1

A=blkdiag(A,A1);

end

A=A+diag(ones(1,nx*(ny-1)),nx)/dy^2+diag(ones(1,nx*(ny-1)),-nx)/dy^2;

phi0=zeros(1,nx*ny)’;

phi0(1:nx:end)=1/dx^2;

phi0(nx:nx:end)=1/dx^2;

phi0(1:nx)=0;

phi0(end-nx+1:end)=0;

B=-rho/ep-phi0;

phi=A\B;

phi=reshape(phi,nx,ny);

surf(phi)

مقاله متلب,مطلب,متلب,مقاله برق,مقاله قدرت,مقاله مطلب,مقاله سیمولینک,دانلود متلب,دانلود مقاله متلب,مقالهmatlab ,آموزش متلب,مطلب,متلب,آموزش برق,آموزش قدرت,آموزش مطلب,آموزش سیمولینک,دانلود متلب,دانلود آموزش متلب,آموزشmatlab ,پروژه متلب,مطلب,متلب,پروژه برق,پروژه قدرت,پروژه مطلب,پروژه سیمولینک,دانلود متلب,دانلود پروژه متلب,پروژهmatlab ,