using Matlab in science problem

I use Finite Difference Method (FDM) in two dimension in my thesis. I want to write my program for FDM method in two dimension in MATLAB, but I couldn't find suitable function for my program, for example I want to solve bellow equation by MATLAB:
alfaf * F(ix ,j)
+ deltaf* F(ix-1,j)
+ miuf * F(ix,j-1)
- landaf* F(ix+1,j)
- Af * F(ix,j+1)
= 0

4 commentaires

Rafael Freire
Rafael Freire le 22 Juil 2011
witch one of the values is a constant and witch one is a variable?
Doug Hull
Doug Hull le 22 Juil 2011
You need to give more information here. I reformated the question to be readable. What have you tried, which are the unknowns. Do you know how to write a function in MATLAB?
Farnaz Tabatabaei
Farnaz Tabatabaei le 23 Juil 2011
Dear Friend;
Please find my complete program that has many errors as bellow answer.
Regards
Farnaz
Farnaz Tabatabaei
Farnaz Tabatabaei le 23 Juil 2011
With many thanks for your kind reply to me, please find the program that I wrote as attached. Please be informed that my equation is a kind of partial differential equation (like laplace equation). F functions are Variable & alfaf, deltaf, miuf, Af are my constant, but they are depended on degrees.
my programm:
% edit program
% example 13900328
clc
clear all
N=input('N=');
% created Boundry condition
kk=0;
for ix=1:N-1;
for j=1:N-1;
kk=kk+1;
%kissi=input('Please input kissi: ');
%Etta=input('Please input Etta: ');
o=0.0000001;
%L=input('Please input L: ');
L=1;
%Omega=input('Please input Omega: ');
Omega=1;
%miu0=input('Please input miu0: ');
h=-pi+0.00000001;
miu0=0.74;
%epsilon0=input('Please input epsilon0: ');
epsilon0=1;
%Betta=input('Please input Betta: ');
Betta=1;
C=3*10^8;
DeltaEtta = L / (N-1);
DeltaKissi = 2 * pi /(N-1);
Etta(ix)=o+(ix-1)*DeltaEtta;
Kissi(j)=h+(j-1)*DeltaKissi;
Etta=Etta(ix);
Kissi=Kissi(j);
k(ix,j) = (2*L.^2).*(((Omega.^2).*epsilon0.*miu0.*(L^2).*(Etta^2).*(cosh(Etta)-sinh(Etta))-(Betta^2).*sinh(Etta).*(cosh(Etta)-cos(Kissi)).*(1-cosh(Etta).*cos(Kissi)))./((Omega^2).*epsilon0.*miu0.*(L^2)*sinh(Etta.^2).*(((Omega.^2)).*epsilon0.*miu0.*(L^2).*sinh(Etta.^2)-2*(Betta.^2).*(cosh(Etta)-cos(Kissi)).^2)+(Betta.^4).*(cosh(Etta)-cos(Kissi)).^4));
q(ix,j) = ((2*L.^2.*(Betta).^2.*sinh(Etta).^2.*sin(Kissi).*(cosh(Etta)-cos(Kissi)))./((Omega).^2.*epsilon0.*miu0.*L.^2.*sinh(Etta).^2.*((Omega).^2.*epsilon0.*miu0.*L.^2.*sinh(Etta).^2.-2.*(Betta).^2.*(cosh(Etta)-cos(Kissi)).^2)+(Betta).^4.*(cosh(Etta)-cos(Kissi)).^4));
GammaT2(ix,j) = ((Omega).^2.*epsilon0.*miu0-(Betta).^2.*((cosh(Etta)-cos(Kissi)).^.2./L.^2.*sinh(Etta)^2));
alfaprim(ix,j) = ((((Betta*sin(Kissi)*(cosh(Etta)-cos(Kissi))^2)/L^3*sinh(Etta))* ...
(-k(ix,j)-(q(ix,j)*(1-cosh(Etta)*cos(Kissi))/(sinh(Etta)*sin(Kissi))-(2/(GammaT2(ix,j))* ...
(DeltaEtta))-(q(ix,j)*(cosh(Etta)-cos(Kissi))/(sin(Kissi)*(DeltaEtta)))+(k(ix,j) * ...
(cosh(Etta)-cos(Kissi))/(sin(Kissi)*(DeltaKissi)))-2*(1-cosh(Etta)* ...
cos(Kissi))/(GammaT2(ix,j))* sin(Kissi)*(DeltaKissi)*sinh(Etta)*(DeltaKissi)))));
deltaprim(ix,j) = (((Betta*sin(Kissi)*(cosh(Etta)-cos(Kissi))^2)/L^3*sinh(Etta))* ...
(2/((GammaT2(ix,j))*(DeltaEtta))+q(ix,j)*(cosh(Etta)-cos(Kissi))/(sin(Kissi)*(DeltaEtta))));
miuprim(ix,j) = ((-(Betta*(cosh(Etta)-cos(Kissi))^2)/L^3*sinh(Etta))* ...
(k(ix,j)*(cosh(Etta)-cos(Kissi))/(DeltaKissi))-2*(1-cosh(Etta)*cos(Kissi)) ...
/((GammaT2(ix,j))*sinh(Etta)*(DeltaKissi)));
alfa(ix,j) = ((((Omega)*(cosh(Etta)-cos(Kissi))^2)/L^2)*((-2/((GammaT2(ix,j))*(DeltaEtta)^2)) ...
-(2/((GammaT2(ix,j))*(DeltaKissi)^2))+(k(ix,j)/(DeltaEtta))+(1-cosh(Etta)*cos(Kissi))/ ...
((GammaT2(ix,j))*sinh(Etta)*(cosh(Etta)-cos(Kissi))*(DeltaEtta))+q(ix,j)/(DeltaKissi)- ...
(sin(Kissi)/(GammaT2(ix,j)*DeltaKissi*cosh(Etta)-cos(Kissi)))));
Delta(ix,j) = ((((Omega)*(cosh(Etta)-cos(Kissi))^2)/L^2)*(1/(DeltaEtta)^2-k(ix,j)/(DeltaEtta)- ...
(1-cosh(Etta)*cos(Kissi))/((GammaT2(ix,j))*(DeltaEtta)*sinh(Etta)*(cosh(Etta) - cos(Kissi)))));
landa(ix,j) = (((Omega)*(cosh(Etta)-cos(Kissi))^2)/((L^2)*(GammaT2(ix,j))*(DeltaEtta)^2));
A(ix,j) = (((Omega)*(cosh(Etta)-cos(Kissi))^2)/((L^2)*(GammaT2(ix,j))*(DeltaKissi)^2));
miu(ix,j)=(((Omega)*(cosh(Etta)-cos(Kissi))^2)/((L^2)*(GammaT2(ix,j))*((1/(DeltaKissi)^2)- ...
(q(ix,j) /(DeltaKissi))+(sin(Kissi)/(DeltaKissi)*(GammaT2(ix,j))*(cosh(Etta)-cos(Kissi))))))
alfaf=alfaprim-i.*alfa;
Deltaf=deltaprim-i.*Delta;
miuf=miuprim-i.*miu;
landaf=i.*landa;
Af=i.*A;
%B=rand(N-1,1);
% AA=rand(N-1,N-1);
%F=rand(N-1,1);
F(ix,j)=(-1./alfaf).*(Deltaf.*F(ix-1,j)+miuf.*F(ix+1,j)-Af.*F(ix,j+1));
B(kk,1)=-subs(F(ix,j),{F(ix,j),F(ix+1,j)},{0,0});
AA(ix,j)=subs(F(ix,j),{F(ix+1,j)},{0,0})-B(kk,1);
% Boundry Condition
if (ix==1)
F(ix,j)=2*i;
if (ix==N)
F(ix,j)=3*i;
if (j==1 || j==N)
F(ix,1)=1+i*1.5;
F(ix,N)=1+i*1.5;
alfaf*F(ix,j)+Deltaf*F(ix-1,j)+miuf*F(ix,j-1)-landaf*F(ix+1,j)-Af*F(ix,j+1)==0;
end
end
end
end
end
Thank you & Best Regards
Farnaz

Connectez-vous pour commenter.

Réponses (0)

Catégories

En savoir plus sur Sciences dans Centre d'aide et File Exchange

Produits

Tags

Aucun tag saisi pour le moment.

Community Treasure Hunt

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

Start Hunting!

Translated by