Skip to content
Snippets Groups Projects
Commit caec3f24 authored by Matteo Cicuttin's avatar Matteo Cicuttin
Browse files

Added EH boundary condition.

parent 82e4a308
No related branches found
No related tags found
No related merge requests found
...@@ -62,6 +62,9 @@ eval_boundary_sources(const model& mod, const parameter_loader& mpl, ...@@ -62,6 +62,9 @@ eval_boundary_sources(const model& mod, const parameter_loader& mpl,
case boundary_condition::H_FIELD: case boundary_condition::H_FIELD:
break; break;
case boundary_condition::EH_FIELD:
break;
case boundary_condition::IMPEDANCE: case boundary_condition::IMPEDANCE:
break; break;
...@@ -164,6 +167,24 @@ eval_boundary_sources_new(const model& mod, const parameter_loader& mpl, ...@@ -164,6 +167,24 @@ eval_boundary_sources_new(const model& mod, const parameter_loader& mpl,
case boundary_condition::SURFACE_CURRENT: case boundary_condition::SURFACE_CURRENT:
break; break;
case boundary_condition::EH_FIELD: {
for (auto& [local, global] : sd.face_numbers)
{
auto fEH = [&](const point_3d& pt) -> std::pair<vec3d, vec3d> {
return mpl.eval_boundary_source_eh(mtag, pt, state.curr_time);
};
matxd prj = project_EH_to_face(e, local, fEH);
auto num_bf = prj.rows();
srcfield.Ex.segment(global*num_bf, num_bf) = prj.col(0);
srcfield.Ey.segment(global*num_bf, num_bf) = prj.col(1);
srcfield.Ez.segment(global*num_bf, num_bf) = prj.col(2);
srcfield.Hx.segment(global*num_bf, num_bf) = prj.col(3);
srcfield.Hy.segment(global*num_bf, num_bf) = prj.col(4);
srcfield.Hz.segment(global*num_bf, num_bf) = prj.col(5);
}
} break; /* case boundary_condition::EH_FIELD */
case boundary_condition::PLANE_WAVE_E: { case boundary_condition::PLANE_WAVE_E: {
for (auto& [local, global] : sd.face_numbers) for (auto& [local, global] : sd.face_numbers)
{ {
......
...@@ -43,6 +43,7 @@ enum class boundary_condition ...@@ -43,6 +43,7 @@ enum class boundary_condition
PLANE_WAVE_H, PLANE_WAVE_H,
E_FIELD, E_FIELD,
H_FIELD, H_FIELD,
EH_FIELD,
SURFACE_CURRENT SURFACE_CURRENT
}; };
...@@ -80,6 +81,7 @@ public: ...@@ -80,6 +81,7 @@ public:
bool volume_has_source(int) const; bool volume_has_source(int) const;
vec3d eval_volume_source(int, const point_3d&, double) const; vec3d eval_volume_source(int, const point_3d&, double) const;
vec3d eval_boundary_source(int, const point_3d&, double) const; vec3d eval_boundary_source(int, const point_3d&, double) const;
std::pair<vec3d,vec3d> eval_boundary_source_eh(int, const point_3d&, double) const;
vec3d eval_interface_source(int, const point_3d&, double) const; vec3d eval_interface_source(int, const point_3d&, double) const;
bool has_analytical_solution(void) const; bool has_analytical_solution(void) const;
......
...@@ -105,6 +105,7 @@ build_source_compression_tables(const model& mod, solver_state_gpu& state, const ...@@ -105,6 +105,7 @@ build_source_compression_tables(const model& mod, solver_state_gpu& state, const
switch ( mpl.boundary_type(bd.gmsh_tag()) ) switch ( mpl.boundary_type(bd.gmsh_tag()) )
{ {
case boundary_condition::PLANE_WAVE_E: case boundary_condition::PLANE_WAVE_E:
case boundary_condition::EH_FIELD:
for (size_t i = 0; i < num_fluxes; i++) for (size_t i = 0; i < num_fluxes; i++)
state.bndsrcs_decomp_table_cpu.push_back(gfnum*num_fluxes+i); state.bndsrcs_decomp_table_cpu.push_back(gfnum*num_fluxes+i);
break; break;
......
...@@ -489,6 +489,7 @@ material_params_gpu::copyin(const material_params& mp) ...@@ -489,6 +489,7 @@ material_params_gpu::copyin(const material_params& mp)
#define BCOND_TYPE_HPLANEW "plane_wave_H" #define BCOND_TYPE_HPLANEW "plane_wave_H"
#define BCOND_TYPE_EFIELD "E_field" #define BCOND_TYPE_EFIELD "E_field"
#define BCOND_TYPE_HFIELD "H_field" #define BCOND_TYPE_HFIELD "H_field"
#define BCOND_TYPE_EHFIELD "EH_field"
#define BCOND_TYPE_SURFCURR "surface_current" #define BCOND_TYPE_SURFCURR "surface_current"
#define SEC_IFCONDS "ifaceconds" #define SEC_IFCONDS "ifaceconds"
...@@ -581,7 +582,7 @@ parameter_loader::validate_boundary_condition(const model&, int tag) const ...@@ -581,7 +582,7 @@ parameter_loader::validate_boundary_condition(const model&, int tag) const
if (bc_kind == BCOND_TYPE_IMPEDANCE) if (bc_kind == BCOND_TYPE_IMPEDANCE)
return true; return true;
if (bc_kind == BCOND_TYPE_EPLANEW) if (bc_kind == BCOND_TYPE_EPLANEW || bc_kind == BCOND_TYPE_EHFIELD)
{ {
if (bc[BC_SOURCE].valid()) if (bc[BC_SOURCE].valid())
return true; return true;
...@@ -792,6 +793,9 @@ parameter_loader::boundary_type(int tag) const ...@@ -792,6 +793,9 @@ parameter_loader::boundary_type(int tag) const
if (kind_str == BCOND_TYPE_HFIELD) if (kind_str == BCOND_TYPE_HFIELD)
return boundary_condition::H_FIELD; return boundary_condition::H_FIELD;
if (kind_str == BCOND_TYPE_EHFIELD)
return boundary_condition::EH_FIELD;
if (kind_str == BCOND_TYPE_SURFCURR) if (kind_str == BCOND_TYPE_SURFCURR)
return boundary_condition::SURFACE_CURRENT; return boundary_condition::SURFACE_CURRENT;
...@@ -846,6 +850,15 @@ parameter_loader::eval_boundary_source(int tag, const point_3d& pt, double t) co ...@@ -846,6 +850,15 @@ parameter_loader::eval_boundary_source(int tag, const point_3d& pt, double t) co
return ret; return ret;
} }
std::pair<vec3d, vec3d>
parameter_loader::eval_boundary_source_eh(int tag, const point_3d& pt, double t) const
{
vec3d retE, retH;
sol::tie(retE(0), retE(1), retE(2), retH(0), retH(1), retH(2)) =
lua[SEC_BCONDS][tag][BC_SOURCE](tag, pt.x(), pt.y(), pt.z(), t);
return std::make_pair(retE, retH);
}
vec3d vec3d
parameter_loader::eval_interface_source(int tag, const point_3d& pt, double t) const parameter_loader::eval_interface_source(int tag, const point_3d& pt, double t) const
{ {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment