-
Notifications
You must be signed in to change notification settings - Fork 2
/
Band_Diagram_main_1.m
114 lines (109 loc) · 4.29 KB
/
Band_Diagram_main_1.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
clear; clc; %close all;
number_of_time_steps=2^14; no_of_k=30;
courant_factor = 0.99;
max_source_freq = 8e10;
n_grid_pt_per_lattice_const=100;
% Constants
eps_0 = 8.854187817e-12; mu_0 = 4*pi*1e-7; c = 1/sqrt(mu_0*eps_0);
% number_of_frequencies = 2^10; % for post_pocess without fft
syms triangular square TE TM PPC MPC PMPC Magnetized_PMPC % PPC type 1 == plasma rods in vaccum background
%triangular=0; square=1; TE=0; TM=1; PPC=0; MPC=1; PMPC=2; Magnetized_PMPC=3; % PPC type 1 == plasma rods in vaccum background
unit_cell = triangular;
mode = TM;
Case = PMPC; ppcType = 1; %type 1 == plasma rods in vaccum background
if Case == PPC
a = 3e-3;
r = sqrt(.5*a^2/pi);
% Plasma parameters
omega_p = 2*pi*28.7e9;
collision_freq = 20e9;
%
MaxDetactablePeak = 1e-30;
elseif Case == PMPC
a=1;
r=0.1; %sqrt(.5*a^2/pi);%0.2;0.2e-2; %
omega_p_Array = (0.1:0.1:2)*pi*c/a;
collision_freq = 0;
end
% q=1.6e-19; m=9.1e-31; %kg
% B0=0.5; %Tesla
% omega_cyclotron=q*B0/m;
% complemetary material
eps_r = 1;
%% ---------------------------------
r_by_a = [0:0.05:0.15 0.2:0.025:0.5];
BandSq_X(length(r_by_a),length(omega_p_Array),30)=0;
BandSq_Gamma=BandSq_X*0;
BandSq_J=BandSq_X*0;
% figure(2);
for plasmaFreqInd = 1:length(omega_p_Array)
omega_p = omega_p_Array(plasmaFreqInd);
omega_str='000000000';
omega_str(1:length(num2str(omega_p*1e-9)))=num2str(omega_p*1e-9);
omega_folder_string = [char(unit_cell) '\' char(mode) '\plasmaFreq_' omega_str '_GHz']; omega_folder_string(omega_folder_string=='.') = '_';
mkdir(omega_folder_string);
end
compute_k_space_parameters;
plasmaFreqInd = 1; rInd=7;
if exist('status1.mat','file'), load status1; end
plasmaFreqInd0=plasmaFreqInd;
rInd0 = rInd;
for rInd=rInd0:12
r = r_by_a(rInd);
r_str='00000';
r_str(1:length(num2str(r)))=num2str(r); r_str(r_str=='.') = '_';
r_folder_string = [char(unit_cell) '\' char(mode) '\r_by_a_' r_str];
mkdir(r_folder_string);
initialize_parameters;
for plasmaFreqInd = plasmaFreqInd0:length(omega_p_Array)
close all
omega_p = omega_p_Array(plasmaFreqInd);
omega_str='000000000';
omega_str(1:length(num2str(omega_p*1e-9)))=num2str(omega_p*1e-9); omega_str(omega_str=='.') = '_';
omega_folder_string = [char(unit_cell) '\' char(mode) '\plasmaFreq_' omega_str '_GHz'];
% mkdir(omega_folder_string);
if mode==TE
disp('mode = TE')
if Case == PPC
time_march_PPC_TE
elseif Case == PMPC
time_march_PMPC_TE
end
else
disp('mode = TM')
time_march_PMPC_TM
end
BandSqOriginal = BandSq;
save([omega_folder_string, '\BandSqOriginal_r_by_a_', r_str],'BandSqOriginal','-MAT');
save([r_folder_string, '\BandSqOriginal_PlasmaFreq_', omega_str, '_GHz'],'BandSqOriginal','-MAT');
%% --- Plotting results
BandSq = BandSqOriginal;
% process_BandSq_2
figure()
if Case == PPC
plot(0:length(kx)-1,BandSq*c/a*1e-9,'x','LineWidth',1,'color',[0 0 0]);
axis([0 length(kx)-1 0 80]); axis1=axis;
else
plot(0:length(kx)-1,BandSq,'x','LineWidth',1,'color',[0 0 0]);
axis([0 length(kx)-1 0 4]); axis1=axis;
end
xlabel('Wave vector','FontSize',20);
ylabel('Normalised Frequency','FontSize',20);
if mode == TE
title('Tranverse Electric (TEz) Photonic Band Structure')
else
title('Trasverse Magnetic (TMz) Photonic Band Structure')
end
%line(0:length(kx)-1,BandSq);
line([vertPos(1)-1 vertPos(1)-1],[0 axis1(4)])
line([vertPos(2)-1 vertPos(2)-1],[0 axis1(4)])
%%
BandSq_Gamma(rInd,plasmaFreqInd,:)=BandSq(end,:); % the first frequency in the lowest band is generally missing
BandSq_J(rInd,plasmaFreqInd,:)=BandSq(vertPos(1),:);
BandSq_X(rInd,plasmaFreqInd,:)=BandSq(vertPos(2),:);
saveas(gcf,[omega_folder_string '\r_by_a_' r_str '.jpg'],'jpg')
saveas(gcf,[r_folder_string '\PlasmaFreq_' omega_str '_GHz.jpg'],'jpg')
save status1
end
plasmaFreqInd0=1;
end