Skip to content

Commit

Permalink
Merge pull request #7 from Amazingkivas/update_branch
Browse files Browse the repository at this point in the history 
Updating PML
  • Loading branch information
@Amazingkivas
Amazingkivas authored Apr 17, 2024
2 parents 3a4db10 + 6b55fd2 commit 210ea1c
Show file tree
Hide file tree
Showing 10 changed files with 289 additions and 263 deletions.
Binary file removed PlotScript/Animations/animation_Bx.gif
View file
Binary file not shown.
Binary file removed PlotScript/Animations/animation_By.gif
View file
Binary file not shown.
Binary file removed PlotScript/Animations/animation_Bz.gif
View file
Binary file not shown.
Binary file removed PlotScript/Animations/animation_Ey.gif
View file
Binary file not shown.
Binary file removed PlotScript/Animations/animation_Ez.gif
View file
Binary file not shown.
Binary file modified PlotScript/src/Release/sample.exe
View file
Binary file not shown.
38 changes: 19 additions & 19 deletions include/FDTD.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,6 @@ class Field
private:
int Ni, Nj, Nk;
std::vector<double> field;
double nul = 0.0;

public:
Field(const int _Ni, const int _Nj, const int _Nk);
Expand All @@ -33,36 +32,34 @@ class FDTD
private:
Field Ex, Ey, Ez, Bx, By, Bz;

Field Exy, Exz, Eyx, Eyz, Ezx, Ezy;
Field Bxy, Bxz, Byx, Byz, Bzx, Bzy;
Field EsigmaX, EsigmaY, EsigmaZ, BsigmaX, BsigmaY, BsigmaZ;

std::vector<Field> Jx;
std::vector<Field> Jy;
std::vector<Field> Jz;

/*std::vector<double> EsigmaX;
std::vector<double> EsigmaY;
std::vector<double> EsigmaZ;
std::vector<double> BsigmaX;
std::vector<double> BsigmaY;
std::vector<double> BsigmaZ;*/
Field Exy, Exz, Eyx, Eyz, Ezx, Ezy;
Field Bxy, Bxz, Byx, Byz, Bzx, Bzy;
Field EsigmaX, EsigmaY, EsigmaZ,
BsigmaX, BsigmaY, BsigmaZ;

Parameters parameters;
PMLparameters pml;
double pml_percent;
double pml_percent, dt;
int pml_size_i, pml_size_j, pml_size_k;
double dt;

void update_E(int bounds_i[2], int bounds_j[2], int bounds_k[2], int t);
void update_B(int bounds_i[2], int bounds_j[2], int bounds_k[2], int t);
void update_B(int bounds_i[2], int bounds_j[2], int bounds_k[2]);

double PMLcoef(double sigma, double constant);
double PMLcoef(double sigma);
void update_E_PML(int bounds_i[2], int bounds_j[2], int bounds_k[2]);
void update_B_PML(int bounds_i[2], int bounds_j[2], int bounds_k[2]);
void set_sigma_x(int bounds_i[2], int bounds_j[2], int bounds_k[2], double SGm[2], std::function<int(int, int, int)> dist);
void set_sigma_y(int bounds_i[2], int bounds_j[2], int bounds_k[2], double SGm[2], std::function<int(int, int, int)> dist);
void set_sigma_z(int bounds_i[2], int bounds_j[2], int bounds_k[2], double SGm[2], std::function<int(int, int, int)> dist);

void set_sigma_x(int bounds_i[2], int bounds_j[2], int bounds_k[2],
double SGm, std::function<int(int, int, int)> dist);
void set_sigma_y(int bounds_i[2], int bounds_j[2], int bounds_k[2],
double SGm, std::function<int(int, int, int)> dist);
void set_sigma_z(int bounds_i[2], int bounds_j[2], int bounds_k[2],
double SGm, std::function<int(int, int, int)> dist);

double calc_max_value(Field& field);

public:
FDTD(Parameters _parameters, double _dt, double _pml_percent);
Expand All @@ -71,4 +68,7 @@ class FDTD
std::vector<Field>& get_current(Component);

std::vector<std::vector<Field>> update_fields(const int time);

double calc_reflection(Field E_start[3], Field B_start[3],
Field E_final[3], Field B_final[3]);
};
19 changes: 3 additions & 16 deletions include/Structures.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,17 +5,15 @@
namespace FDTDconst
{
const double C = 3e10;
const double R = 0.01;
const double EPS0 = 1.60217656535e-12;
const double PI = 3.1415926535;
const double MU0 = 1e6 * EPS0;
const double R = 1e-12;
const double EPS0 = 1.0;
const double MU0 = EPS0;
const double N = 4.0;
}

namespace FDTDstruct
{
enum class Component { EX, EY, EZ, BX, BY, BZ, JX, JY, JZ };

enum class Axis { X, Y, Z };

struct SelectedFields {
Expand All @@ -33,26 +31,15 @@ namespace FDTDstruct
double period_z = static_cast<double>(m) * FDTDconst::C;
};

struct PMLparameters {
double sigma_low_xy;
double sigma_low_yz;
double sigma_low_zx;
double sigma_up_xy;
double sigma_up_yz;
double sigma_up_zx;
};

struct Parameters {
int Ni;
int Nj;
int Nk;

double ax;
double bx;

double ay;
double by;

double az;
double bz;

Expand Down
43 changes: 21 additions & 22 deletions samples/sample.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,8 @@ Axis get_axis(Component field_E, Component field_B)
return selected_axis;
}

void plane_wave(std::vector<char*> arguments, std::vector<double> numbers, char* outfile_path)
void plane_wave(std::vector<char*> arguments,
std::vector<double> numbers, char* outfile_path)
{
// Saving data from the console
Component selected_field_1 = static_cast<Component>(std::atoi(arguments[1]));
Expand All @@ -53,17 +54,13 @@ void plane_wave(std::vector<char*> arguments, std::vector<double> numbers, char*
double sizes_z[2] = { numbers[7], numbers[8] }; // { az, bz }
int iterations_num = static_cast<int>(numbers[9]);
double time = numbers[10];
double time_step = 0.4; // time / static_cast<double>(iterations_num); // dt
double time_step = time / static_cast<double>(iterations_num); // dt

// Initialization of the initializing function and the true solution function
std::function<double(double, double[2])> initial_func = [](double x, double size[2])
{
return sin(2.0 * M_PI * (x - size[0]) / (size[1] - size[0]));
};
std::function<double(double, double, double[2])> true_func = [](double x, double t, double size[2])
{
return sin(2.0 * M_PI * (x - size[0] - FDTDconst::C * t) / (size[1] - size[0]));
};

// Determination of the wave propagation axis
Axis axis = get_axis(flds_selected[0], flds_selected[1]);
Expand All @@ -80,17 +77,17 @@ void plane_wave(std::vector<char*> arguments, std::vector<double> numbers, char*
grid_sizes[0],
grid_sizes[1],
grid_sizes[2],
0.0,
sizes_x[0],
sizes_x[1],
0.0,
sizes_y[0],
sizes_y[1],
0.0,
sizes_z[0],
sizes_z[1],
FDTDconst::C, //(sizes_x[1] - sizes_x[0]) / static_cast<double>(grid_sizes[0]),
FDTDconst::C, //(sizes_y[1] - sizes_y[0]) / static_cast<double>(grid_sizes[1]),
FDTDconst::C, //(sizes_z[1] - sizes_z[0]) / static_cast<double>(grid_sizes[2])
(sizes_x[1] - sizes_x[0]) / static_cast<double>(grid_sizes[0]),
(sizes_y[1] - sizes_y[0]) / static_cast<double>(grid_sizes[1]),
(sizes_z[1] - sizes_z[0]) / static_cast<double>(grid_sizes[2])
};
FDTD method(params, time_step, 0.1);
FDTD method(params, time_step, 0.0);

// Meaningful calculations
Test_FDTD test(params);
Expand All @@ -104,15 +101,13 @@ void plane_wave(std::vector<char*> arguments, std::vector<double> numbers, char*
}

std::vector<std::vector<Field>> data = method.update_fields(iterations_num);

//std::cout << test.get_max_abs_error(method.get_field(fld_tested), fld_tested, true_func, time) << std::endl;

if (write_flag == true)
{
// Writing the results to a file
try
{
write_spherical(data, Axis::X, iterations_num, outfile_path);
write_plane(method, axis, outfile_path);
}
catch (const std::exception& except)
{
Expand All @@ -125,18 +120,22 @@ void spherical_wave(int n, int it, char* base_path = "")
{
CurrentParameters cur_param
{
2,
8,
4,
0.4,
0.2,
};
double T = cur_param.period;
double Tx = cur_param.period_x;
double Ty = cur_param.period_y;
double Tz = cur_param.period_z;
cur_param.iterations = static_cast<int>(static_cast<double>(cur_param.period) / cur_param.dt);
std::function<double(double, double, double, double)> cur_func = [T, Tx, Ty, Tz](double x, double y, double z, double t)
std::function<double(double, double, double, double)> cur_func
= [T, Tx, Ty, Tz](double x, double y, double z, double t)
{
return sin(2.0 * M_PI * t / T) * pow(cos(2.0 * M_PI * x / Tx), 2.0) * pow(cos(2.0 * M_PI * y / Ty), 2.0) * pow(cos(2.0 * M_PI * z / Tz), 2.0);
return sin(2.0 * M_PI * t / T)
* pow(cos(2.0 * M_PI * x / Tx), 2.0)
* pow(cos(2.0 * M_PI * y / Ty), 2.0)
* pow(cos(2.0 * M_PI * z / Tz), 2.0);
};

// Initialization of the structures and method
Expand Down Expand Up @@ -182,8 +181,8 @@ int main(int argc, char* argv[])
if (argc == 1)
{
#ifdef __USE_SPHERICAL_WAVE__
int N = 50; // 100;
int Iterations = 250; // 200;
int N = 70;
int Iterations = 300;
spherical_wave(N, Iterations, "../../PlotScript/");
#endif

Expand Down
Loading

0 comments on commit 210ea1c

Please sign in to comment.