FDTD C mã nguồn từ cuốn sách

bạn có thể tham khảo liên kết:
http://www.edaboard.com/viewtopic.php?t=74057
ebook này là do điện từ mô phỏng FDTD (do Sullivan)

và liên kết:
http://www.edaboard.com/viewtopic.php?t=149266
ebook này là điện từ của FDTD (do kunz)

và mã MATLAB từ liên kết:
http://www.mathworks.com/matlabcentral/fileexchange/loadCategory.do
bạn phải tìm kiếm về FDTD MATLAB mã sau đó bắt đầu tải về.
1D này là ví dụ về phương pháp FDTD:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%
% Scott Hudson, WSU Tri-Cities
% 1D điện-sự khác biệt hữu hạn thời gian tên miền (FDTD) chương trình.
Ey% giả định và các thành phần lĩnh vực Hz tuyên truyền theo hướng x.
% Fields, permittivity, tính thấm, và dẫn
% là chức năng của x.Hãy thử thay đổi giá trị của "tiểu sử".
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%

đóng tất cả;
rõ ràng tất cả;

L = 5;% miền dài trong mét
N = 505;% # không gian mẫu ở miền
Tiêu thạch = 800;% # của lặp để thực hiện
fs = 300e6;% nguồn tần số trong Hz
ds = L / N;% không gian bước trong mét
dt = ds/300e6;% "ma thuật thời gian bước"
eps0 = 8.854e-12;% permittivity không gian Việt
mu0 = pi * 4e-7;% tính thấm không gian Việt
x = linspace (0, L, N);% x phối hợp của các mẫu không gian
showWKB = 0;% nếu = 1 sau đó hiển thị WKB appromination lúc kết thúc

% quy mô các yếu tố cho E và H
ae = những (N, 1) * dt / (ds * eps0);
am = cái (N, 1) * dt / (ds * mu0);
như là = cái (N, 1);
epsr = cái (N, 1);
Mur = cái (N, 1);
sigma = zeros (N, 1);%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%
% Ở đây chúng tôi chỉ định epsilon, sigma, và hồ sơ mu.Tôi đã
% được xác định trước một số ví dụ thú vị.Hãy thử cấu hình = 1,2,3,4,5,6 theo thứ tự.
% Bạn có thể định nghĩa epsr (i), Mur (i) (tương permittivity và tính thấm)
% và sigma (i) (dẫn) để được bất cứ điều gì bạn muốn.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%
profile = 1;
for i = 1: N
epsr (i) = 1;
Mur (i) = 1;
W1 = 0,5;
W2 = 1,5;
nếu (profile == 1)% điện môi cửa sổ
if (abs (x (i)-L / 2) <0,5) epsr (i) = 4; cuối
cuối
nếu (profile == 2)% điện môi cửa sổ với chuyển đổi suôn sẻ
if (abs (x (i)-L / 2) <1,5) epsr (i) = 1 3 * (1 cos (pi * (abs (x (i)-L / 2)-W1) / (W2 -W1))) / 2; cuối
if (abs (x (i)-L / 2) <0,5) epsr (i) = 4; cuối
cuối
nếu (profile == 3)% dielectric gián đoạn
if (x (i)> L / 2) epsr (i) = 9; cuối
cuối
nếu (profile == 4)% dielectric disontinuity với lớp kết hợp 1/4-wave
if (x (i)> (L/2-0.1443)) epsr (i) = 3; cuối
if (x (i)> L / 2) epsr (i) = 9; cuối
cuối
nếu (profile == 5)% tiến hành một nửa không gian
if (x (i)> L / 2) sigma (i) = 0,005; cuối
cuối
nếu (profile == 6)% sinusoidal điện môi
epsr (i) = 1 sin (2 * pi * x (i) / L) ^ 2;
cuối
cuối
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%ae = ae. / epsr;
am = am. / Mur;
ae = ae. / (1 dt * (sigma. / epsr) / (2 * eps0));
như = (1-dt * (sigma. / epsr) / (2 * eps0 ))./( 1 dt * (sigma. / epsr) / (2 * eps0));

% lô các permittivity, tính thấm, và dẫn các cấu hình
Con số (1)
subplot (3,1,1);
plot (x, epsr);
Lưới điện trên;
trục ([3 * ds L min (epsr) * 0,9 max (epsr) * 1,1]);
Tiêu đề permittivity tương đối ( '');
subplot (3,1,2);
plot (x, Mur);
Lưới điện trên;
trục ([3 * ds L min (Mur) * 0,9 max (Mur) * 1,1]);
Tiêu đề permeabiliity tương đối ( '');
subplot (3,1,3);
plot (x, sigma);
Lưới điện trên;
trục ([3 * ds L min (sigma) * 0.9-0.001 tối đa (sigma) * 1,1 0,001]);
Tiêu đề ( 'dẫn');

% khởi tạo các lĩnh vực về không
Hz = zeros (N, 1);
Ey = zeros (N, 1);
con số (2);
thiết lập (GCF, 'doublebuffer', 'ngày');% đặt trên đôi đệm cho đồ họa mượt mà
lô (Ey);
Lưới điện trên;

cho ITER = 1: tiêu thạch
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%
The 10% hoặc hơn các dòng tiếp theo của mã là nơi chúng tôi thực sự tích hợp của Maxwell
% phương trình.Tất cả các phần còn lại của chương trình là cơ bản sổ sách kế toán và vẽ.
% "mịn bật" nguồn sinusoidal
Ey (3) = Ey (3) 2 * (1-exp (- ((ITER-1) / 50) ^ 2)) * sin (2 * pi * fs * dt * ITER);
Hz (1) = Hz (2);% hấp thụ điều kiện ranh giới cho trái-tuyên truyền sóng
for i = 2: N-1% cập nhật H trường
Hz (i) = Hz (i)-am (i) * (Ey (i 1)-Ey (i));
cuối
Ey (N) = Ey (N-1);% hấp thụ điều kiện ranh giới cho quyền tuyên truyền sóng
for i = 2: N-1% cập nhật E trường
Ey (i) = như (i) * Ey (i)-ae (i) * (Hz (i)-Hz (i-1));
cuối
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%
con số (2)
giữ off
plot (x, Ey, 'b');
trục ([3 * ds L -2 2]);
Lưới điện trên;
Tiêu đề ( 'E (màu xanh) và 377 * H (đỏ)');
giữ lấy
plot (x, 377 * Hz, 'r');
xlabel ( 'x (m)');
tạm dừng (0);
ITER
cuối

% WKB dự đoán cho các lĩnh vực
nếu (showWKB == 1)
giai đoạn = cumsum ((epsr) ^ 0,5.) * ds;
beta0 = 2 * pi * fs / (300e6);
Lý thuyết = sin (2 * pi * fs * (tiêu thạch 4) * dt-beta0 * giai đoạn). / (epsr. ^ 0,25);
đầu vào ( 'bấm Enter để hiển thị WKB lý thuyết');
plot (x, lý thuyết, 'sinh');
Lý thuyết = sin (2 * pi * fs * (tiêu thạch 4) * dt-beta0 * giai đoạn) .* (epsr. ^ 0,25);
plot (x, lý thuyết, 'r.');
Tiêu đề ( 'E (màu xanh), 377 * H (đỏ), WKB lý thuyết (điểm)');
cuối

 

Welcome to EDABoard.com

Sponsor

Back
Top