zcy
/
openEMS
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

dispersive material extension: prepared for higher order dispersion

pull/1/head
Thorsten Liebig 2012-04-24 11:05:19 +02:00
parent fb5b330d25
commit 796fd83f7b
7 changed files with 288 additions and 192 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

110
FDTD/extensions/engine_ext_dispersive.cpp
View File

@ -22,117 +22,141 @@
Engine_Ext_Dispersive::Engine_Ext_Dispersive(Operator_Ext_Dispersive* op_ext_disp) : Engine_Extension(op_ext_disp)
{
m_Op_Ext_Disp = op_ext_disp;
int order = m_Op_Ext_Disp->m_Order;
curr_ADE = new FDTD_FLOAT**[order];
volt_ADE = new FDTD_FLOAT**[order];
for (int o=0;o<order;++o)
{
curr_ADE[o] = new FDTD_FLOAT*[3];
volt_ADE[o] = new FDTD_FLOAT*[3];
for (int n=0; n<3; ++n)
{
if (m_Op_Ext_Disp->m_curr_ADE_On==true)
if (m_Op_Ext_Disp->m_curr_ADE_On[o]==true)
{
curr_ADE[n] = new FDTD_FLOAT[m_Op_Ext_Disp->m_LM_Count];
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
curr_ADE[n][i]=0.0;
curr_ADE[o][n] = new FDTD_FLOAT[m_Op_Ext_Disp->m_LM_Count[o]];
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count[o]; ++i)
curr_ADE[o][n][i]=0.0;
}
else
curr_ADE[n] = NULL;
if (m_Op_Ext_Disp->m_volt_ADE_On==true)
curr_ADE[o][n] = NULL;
if (m_Op_Ext_Disp->m_volt_ADE_On[o]==true)
{
volt_ADE[n] = new FDTD_FLOAT[m_Op_Ext_Disp->m_LM_Count];
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
volt_ADE[n][i]=0.0;
volt_ADE[o][n] = new FDTD_FLOAT[m_Op_Ext_Disp->m_LM_Count[o]];
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count[o]; ++i)
volt_ADE[o][n][i]=0.0;
}
else
volt_ADE[n] = NULL;
volt_ADE[o][n] = NULL;
}
}
}
Engine_Ext_Dispersive::~Engine_Ext_Dispersive()
{
if (curr_ADE==NULL && volt_ADE==NULL)
return;
for (int o=0;o<m_Op_Ext_Disp->m_Order;++o)
{
for (int n=0; n<3; ++n)
{
delete[] curr_ADE[n];
curr_ADE[n] = NULL;
delete[] volt_ADE[n];
volt_ADE[n] = NULL;
delete[] curr_ADE[o][n];
delete[] volt_ADE[o][n];
}
delete[] curr_ADE[o];
delete[] volt_ADE[o];
}
delete[] curr_ADE;
curr_ADE=NULL;
delete[] volt_ADE;
volt_ADE=NULL;
}
void Engine_Ext_Dispersive::Apply2Voltages()
{
if (m_Op_Ext_Disp->m_volt_ADE_On==false) return;
for (int o=0;o<m_Op_Ext_Disp->m_Order;++o)
{
if (m_Op_Ext_Disp->m_volt_ADE_On[o]==false) continue;
unsigned int **pos = m_Op_Ext_Disp->m_LM_pos;
unsigned int **pos = m_Op_Ext_Disp->m_LM_pos[o];
//switch for different engine types to access faster inline engine functions
switch (m_Eng->GetType())
{
case Engine::BASIC:
{
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count.at(o); ++i)
{
m_Eng->Engine::SetVolt(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[0][i]);
m_Eng->Engine::SetVolt(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[1][i]);
m_Eng->Engine::SetVolt(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[2][i]);
m_Eng->Engine::SetVolt(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][0][i]);
m_Eng->Engine::SetVolt(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][1][i]);
m_Eng->Engine::SetVolt(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][2][i]);
}
break;
}
case Engine::SSE:
{
Engine_sse* eng_sse = (Engine_sse*)m_Eng;
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count.at(o); ++i)
{
eng_sse->Engine_sse::SetVolt(0,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[0][i]);
eng_sse->Engine_sse::SetVolt(1,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[1][i]);
eng_sse->Engine_sse::SetVolt(2,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[2][i]);
eng_sse->Engine_sse::SetVolt(0,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][0][i]);
eng_sse->Engine_sse::SetVolt(1,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][1][i]);
eng_sse->Engine_sse::SetVolt(2,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][2][i]);
}
break;
}
default:
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count.at(o); ++i)
{
m_Eng->SetVolt(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetVolt(0,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[0][i]);
m_Eng->SetVolt(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetVolt(1,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[1][i]);
m_Eng->SetVolt(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetVolt(2,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[2][i]);
m_Eng->SetVolt(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetVolt(0,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][0][i]);
m_Eng->SetVolt(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetVolt(1,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][1][i]);
m_Eng->SetVolt(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetVolt(2,pos[0][i],pos[1][i],pos[2][i]) - volt_ADE[o][2][i]);
}
break;
}
}
}
void Engine_Ext_Dispersive::Apply2Current()
{
if (m_Op_Ext_Disp->m_curr_ADE_On==false) return;
for (int o=0;o<m_Op_Ext_Disp->m_Order;++o)
{
if (m_Op_Ext_Disp->m_curr_ADE_On[o]==false) continue;
unsigned int **pos = m_Op_Ext_Disp->m_LM_pos;
unsigned int **pos = m_Op_Ext_Disp->m_LM_pos[o];
//switch for different engine types to access faster inline engine functions
switch (m_Eng->GetType())
{
case Engine::BASIC:
{
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count.at(o); ++i)
{
m_Eng->Engine::SetCurr(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[0][i]);
m_Eng->Engine::SetCurr(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[1][i]);
m_Eng->Engine::SetCurr(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[2][i]);
m_Eng->Engine::SetCurr(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][0][i]);
m_Eng->Engine::SetCurr(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][1][i]);
m_Eng->Engine::SetCurr(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->Engine::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][2][i]);
}
break;
}
case Engine::SSE:
{
Engine_sse* eng_sse = (Engine_sse*)m_Eng;
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count.at(o); ++i)
{
eng_sse->Engine_sse::SetCurr(0,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[0][i]);
eng_sse->Engine_sse::SetCurr(1,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[1][i]);
eng_sse->Engine_sse::SetCurr(2,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[2][i]);
eng_sse->Engine_sse::SetCurr(0,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][0][i]);
eng_sse->Engine_sse::SetCurr(1,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][1][i]);
eng_sse->Engine_sse::SetCurr(2,pos[0][i],pos[1][i],pos[2][i], eng_sse->Engine_sse::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][2][i]);
}
break;
}
default:
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Disp->m_LM_Count.at(o); ++i)
{
m_Eng->SetCurr(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetCurr(0,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[0][i]);
m_Eng->SetCurr(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetCurr(1,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[1][i]);
m_Eng->SetCurr(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetCurr(2,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[2][i]);
m_Eng->SetCurr(0,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetCurr(0,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][0][i]);
m_Eng->SetCurr(1,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetCurr(1,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][1][i]);
m_Eng->SetCurr(2,pos[0][i],pos[1][i],pos[2][i], m_Eng->GetCurr(2,pos[0][i],pos[1][i],pos[2][i]) - curr_ADE[o][2][i]);
}
break;
}
}
}

9
FDTD/extensions/engine_ext_dispersive.h
View File

@ -36,11 +36,16 @@ public:
protected:
Operator_Ext_Dispersive* m_Op_Ext_Disp;
//! Dispersive order
int m_Order;
//! ADE currents
FDTD_FLOAT *curr_ADE[3];
// Array setup: curr_ADE[N_order][direction][mesh_pos]
FDTD_FLOAT ***curr_ADE;
//! ADE voltages
FDTD_FLOAT *volt_ADE[3];
// Array setup: volt_ADE[N_order][direction][mesh_pos]
FDTD_FLOAT ***volt_ADE;
};
#endif // ENGINE_EXT_DISPERSIVE_H

99
FDTD/extensions/engine_ext_lorentzmaterial.cpp
View File

@ -22,6 +22,7 @@
Engine_Ext_LorentzMaterial::Engine_Ext_LorentzMaterial(Operator_Ext_LorentzMaterial* op_ext_lorentz) : Engine_Ext_Dispersive(op_ext_lorentz)
{
m_Op_Ext_Lor = op_ext_lorentz;
m_Order = m_Op_Ext_Lor->GetDispersionOrder();
}
Engine_Ext_LorentzMaterial::~Engine_Ext_LorentzMaterial()
@ -31,113 +32,119 @@ Engine_Ext_LorentzMaterial::~Engine_Ext_LorentzMaterial()
void Engine_Ext_LorentzMaterial::DoPreVoltageUpdates()
{
if (m_Op_Ext_Lor->m_volt_ADE_On==false) return;
for (int o=0;o<m_Order;++o)
{
if (m_Op_Ext_Lor->m_volt_ADE_On==false) continue;
unsigned int **pos = m_Op_Ext_Lor->m_LM_pos;
unsigned int **pos = m_Op_Ext_Lor->m_LM_pos[o];
//switch for different engine types to access faster inline engine functions
switch (m_Eng->GetType())
{
case Engine::BASIC:
{
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count.at(o); ++i)
{
volt_ADE[0][i] *= m_Op_Ext_Lor->v_int_ADE[0][i];
volt_ADE[0][i] += m_Op_Ext_Lor->v_ext_ADE[0][i] * m_Eng->Engine::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][0][i] *= m_Op_Ext_Lor->v_int_ADE[o][0][i];
volt_ADE[o][0][i] += m_Op_Ext_Lor->v_ext_ADE[o][0][i] * m_Eng->Engine::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[1][i] *= m_Op_Ext_Lor->v_int_ADE[1][i];
volt_ADE[1][i] += m_Op_Ext_Lor->v_ext_ADE[1][i] * m_Eng->Engine::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][1][i] *= m_Op_Ext_Lor->v_int_ADE[o][1][i];
volt_ADE[o][1][i] += m_Op_Ext_Lor->v_ext_ADE[o][1][i] * m_Eng->Engine::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[2][i] *= m_Op_Ext_Lor->v_int_ADE[2][i];
volt_ADE[2][i] += m_Op_Ext_Lor->v_ext_ADE[2][i] * m_Eng->Engine::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][2][i] *= m_Op_Ext_Lor->v_int_ADE[o][2][i];
volt_ADE[o][2][i] += m_Op_Ext_Lor->v_ext_ADE[o][2][i] * m_Eng->Engine::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]);
}
break;
}
case Engine::SSE:
{
Engine_sse* eng_sse = (Engine_sse*)m_Eng;
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count.at(o); ++i)
{
volt_ADE[0][i] *= m_Op_Ext_Lor->v_int_ADE[0][i];
volt_ADE[0][i] += m_Op_Ext_Lor->v_ext_ADE[0][i] * eng_sse->Engine_sse::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][0][i] *= m_Op_Ext_Lor->v_int_ADE[o][0][i];
volt_ADE[o][0][i] += m_Op_Ext_Lor->v_ext_ADE[o][0][i] * eng_sse->Engine_sse::GetVolt(0,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[1][i] *= m_Op_Ext_Lor->v_int_ADE[1][i];
volt_ADE[1][i] += m_Op_Ext_Lor->v_ext_ADE[1][i] * eng_sse->Engine_sse::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][1][i] *= m_Op_Ext_Lor->v_int_ADE[o][1][i];
volt_ADE[o][1][i] += m_Op_Ext_Lor->v_ext_ADE[o][1][i] * eng_sse->Engine_sse::GetVolt(1,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[2][i] *= m_Op_Ext_Lor->v_int_ADE[2][i];
volt_ADE[2][i] += m_Op_Ext_Lor->v_ext_ADE[2][i] * eng_sse->Engine_sse::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][2][i] *= m_Op_Ext_Lor->v_int_ADE[o][2][i];
volt_ADE[o][2][i] += m_Op_Ext_Lor->v_ext_ADE[o][2][i] * eng_sse->Engine_sse::GetVolt(2,pos[0][i],pos[1][i],pos[2][i]);
}
break;
}
default:
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count.at(o); ++i)
{
volt_ADE[0][i] *= m_Op_Ext_Lor->v_int_ADE[0][i];
volt_ADE[0][i] += m_Op_Ext_Lor->v_ext_ADE[0][i] * m_Eng->GetVolt(0,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][0][i] *= m_Op_Ext_Lor->v_int_ADE[o][0][i];
volt_ADE[o][0][i] += m_Op_Ext_Lor->v_ext_ADE[o][0][i] * m_Eng->GetVolt(0,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[1][i] *= m_Op_Ext_Lor->v_int_ADE[1][i];
volt_ADE[1][i] += m_Op_Ext_Lor->v_ext_ADE[1][i] * m_Eng->GetVolt(1,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][1][i] *= m_Op_Ext_Lor->v_int_ADE[o][1][i];
volt_ADE[o][1][i] += m_Op_Ext_Lor->v_ext_ADE[o][1][i] * m_Eng->GetVolt(1,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[2][i] *= m_Op_Ext_Lor->v_int_ADE[2][i];
volt_ADE[2][i] += m_Op_Ext_Lor->v_ext_ADE[2][i] * m_Eng->GetVolt(2,pos[0][i],pos[1][i],pos[2][i]);
volt_ADE[o][2][i] *= m_Op_Ext_Lor->v_int_ADE[o][2][i];
volt_ADE[o][2][i] += m_Op_Ext_Lor->v_ext_ADE[o][2][i] * m_Eng->GetVolt(2,pos[0][i],pos[1][i],pos[2][i]);
}
break;
}
}
}
void Engine_Ext_LorentzMaterial::DoPreCurrentUpdates()
{
if (m_Op_Ext_Lor->m_curr_ADE_On==false) return;
for (int o=0;o<m_Order;++o)
{
if (m_Op_Ext_Lor->m_curr_ADE_On==false) continue;
unsigned int **pos = m_Op_Ext_Lor->m_LM_pos;
unsigned int **pos = m_Op_Ext_Lor->m_LM_pos[o];
//switch for different engine types to access faster inline engine functions
switch (m_Eng->GetType())
{
case Engine::BASIC:
{
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count.at(o); ++i)
{
curr_ADE[0][i] *= m_Op_Ext_Lor->i_int_ADE[0][i];
curr_ADE[0][i] += m_Op_Ext_Lor->i_ext_ADE[0][i] * m_Eng->Engine::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][0][i] *= m_Op_Ext_Lor->i_int_ADE[o][0][i];
curr_ADE[o][0][i] += m_Op_Ext_Lor->i_ext_ADE[o][0][i] * m_Eng->Engine::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[1][i] *= m_Op_Ext_Lor->i_int_ADE[1][i];
curr_ADE[1][i] += m_Op_Ext_Lor->i_ext_ADE[1][i] * m_Eng->Engine::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][1][i] *= m_Op_Ext_Lor->i_int_ADE[o][1][i];
curr_ADE[o][1][i] += m_Op_Ext_Lor->i_ext_ADE[o][1][i] * m_Eng->Engine::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[2][i] *= m_Op_Ext_Lor->i_int_ADE[2][i];
curr_ADE[2][i] += m_Op_Ext_Lor->i_ext_ADE[2][i] * m_Eng->Engine::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][2][i] *= m_Op_Ext_Lor->i_int_ADE[o][2][i];
curr_ADE[o][2][i] += m_Op_Ext_Lor->i_ext_ADE[o][2][i] * m_Eng->Engine::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]);
}
break;
}
case Engine::SSE:
{
Engine_sse* eng_sse = (Engine_sse*)m_Eng;
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count.at(o); ++i)
{
curr_ADE[0][i] *= m_Op_Ext_Lor->i_int_ADE[0][i];
curr_ADE[0][i] += m_Op_Ext_Lor->i_ext_ADE[0][i] * eng_sse->Engine_sse::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][0][i] *= m_Op_Ext_Lor->i_int_ADE[o][0][i];
curr_ADE[o][0][i] += m_Op_Ext_Lor->i_ext_ADE[o][0][i] * eng_sse->Engine_sse::GetCurr(0,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[1][i] *= m_Op_Ext_Lor->i_int_ADE[1][i];
curr_ADE[1][i] += m_Op_Ext_Lor->i_ext_ADE[1][i] * eng_sse->Engine_sse::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][1][i] *= m_Op_Ext_Lor->i_int_ADE[o][1][i];
curr_ADE[o][1][i] += m_Op_Ext_Lor->i_ext_ADE[o][1][i] * eng_sse->Engine_sse::GetCurr(1,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[2][i] *= m_Op_Ext_Lor->i_int_ADE[2][i];
curr_ADE[2][i] += m_Op_Ext_Lor->i_ext_ADE[2][i] * eng_sse->Engine_sse::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][2][i] *= m_Op_Ext_Lor->i_int_ADE[o][2][i];
curr_ADE[o][2][i] += m_Op_Ext_Lor->i_ext_ADE[o][2][i] * eng_sse->Engine_sse::GetCurr(2,pos[0][i],pos[1][i],pos[2][i]);
}
break;
}
default:
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count; ++i)
for (unsigned int i=0; i<m_Op_Ext_Lor->m_LM_Count.at(o); ++i)
{
curr_ADE[0][i] *= m_Op_Ext_Lor->i_int_ADE[0][i];
curr_ADE[0][i] += m_Op_Ext_Lor->i_ext_ADE[0][i] * m_Eng->GetCurr(0,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][0][i] *= m_Op_Ext_Lor->i_int_ADE[o][0][i];
curr_ADE[o][0][i] += m_Op_Ext_Lor->i_ext_ADE[o][0][i] * m_Eng->GetCurr(0,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[1][i] *= m_Op_Ext_Lor->i_int_ADE[1][i];
curr_ADE[1][i] += m_Op_Ext_Lor->i_ext_ADE[1][i] * m_Eng->GetCurr(1,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][1][i] *= m_Op_Ext_Lor->i_int_ADE[o][1][i];
curr_ADE[o][1][i] += m_Op_Ext_Lor->i_ext_ADE[o][1][i] * m_Eng->GetCurr(1,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[2][i] *= m_Op_Ext_Lor->i_int_ADE[2][i];
curr_ADE[2][i] += m_Op_Ext_Lor->i_ext_ADE[2][i] * m_Eng->GetCurr(2,pos[0][i],pos[1][i],pos[2][i]);
curr_ADE[o][2][i] *= m_Op_Ext_Lor->i_int_ADE[o][2][i];
curr_ADE[o][2][i] += m_Op_Ext_Lor->i_ext_ADE[o][2][i] * m_Eng->GetCurr(2,pos[0][i],pos[1][i],pos[2][i]);
}
break;
}
}
}

42
FDTD/extensions/operator_ext_dispersive.cpp
View File

@ -21,30 +21,44 @@
Operator_Ext_Dispersive::Operator_Ext_Dispersive(Operator* op) : Operator_Extension(op)
{
m_curr_ADE_On = false;
m_volt_ADE_On = false;
m_curr_ADE_On = NULL;
m_volt_ADE_On = NULL;
m_LM_pos[0]=NULL;
m_LM_pos[1]=NULL;
m_LM_pos[2]=NULL;
m_LM_pos=NULL;
m_curr_ADE_On=NULL;
m_volt_ADE_On=NULL;
m_Order = 0;
}
Operator_Ext_Dispersive::~Operator_Ext_Dispersive()
{
delete[] m_LM_pos[0];
delete[] m_LM_pos[1];
delete[] m_LM_pos[2];
delete[] m_curr_ADE_On;
delete[] m_volt_ADE_On;
m_curr_ADE_On=NULL;
m_volt_ADE_On=NULL;
m_LM_pos[0]=NULL;
m_LM_pos[1]=NULL;
m_LM_pos[2]=NULL;
for (int n=0;n<m_Order;++n)
{
delete[] m_LM_pos[n][0];
delete[] m_LM_pos[n][1];
delete[] m_LM_pos[n][2];
}
delete[] m_LM_pos;
m_LM_pos=NULL;
m_Order=0;
m_LM_Count.clear();
}
void Operator_Ext_Dispersive::ShowStat(ostream &ostr) const
{
Operator_Extension::ShowStat(ostr);
string On_Off[2] = {"Off", "On"};
ostr << " Active cells\t\t: " << m_LM_Count << endl;
ostr << " Voltage ADE is \t: " << On_Off[m_volt_ADE_On] << endl;
ostr << " Current ADE is \t: " << On_Off[m_curr_ADE_On] << endl;
ostr << " Max. Dispersion Order N = " << m_Order << endl;
for (int i=0;i<m_Order;++i)
{
ostr << " N=" << i << ":\t Active cells\t\t: " << m_LM_Count[i] << endl;
ostr << " N=" << i << ":\t Voltage ADE is \t: " << On_Off[m_volt_ADE_On[i]] << endl;
ostr << " N=" << i << ":\t Current ADE is \t: " << On_Off[m_curr_ADE_On[i]] << endl;
}
}

19
FDTD/extensions/operator_ext_dispersive.h
View File

@ -20,6 +20,7 @@
//#include "operator.h"
#include "operator_extension.h"
#include "vector"
//! Abstract base class for all dispersive material models, based on an ADE (additional differential equation)
class Operator_Ext_Dispersive : public Operator_Extension
@ -28,6 +29,8 @@ class Operator_Ext_Dispersive : public Operator_Extension
public:
virtual ~Operator_Ext_Dispersive();
virtual int GetDispersionOrder() {return m_Order;}
virtual string GetExtensionName() const {return string("Dispersive Material Abstract Base class");}
virtual void ShowStat(ostream &ostr) const;
@ -35,13 +38,17 @@ public:
protected:
Operator_Ext_Dispersive(Operator* op);
//! Lorentz material count
unsigned int m_LM_Count;
//! Index with dispersive material
unsigned int *m_LM_pos[3];
//! Dispersive order
int m_Order;
bool m_curr_ADE_On;
bool m_volt_ADE_On;
//! Dispersive material count
vector<unsigned int> m_LM_Count;
//! Index with dispersive material
// Array setup: m_LM_pos[N_order][direction][mesh_pos]
unsigned int ***m_LM_pos;
bool *m_curr_ADE_On;
bool *m_volt_ADE_On;
};
#endif // OPERATOR_EXT_DISPERSIVE_H

104
FDTD/extensions/operator_ext_lorentzmaterial.cpp
View File

@ -20,28 +20,48 @@
Operator_Ext_LorentzMaterial::Operator_Ext_LorentzMaterial(Operator* op) : Operator_Ext_Dispersive(op)
{
for (int n=0; n<3; ++n)
{
v_int_ADE[n] = NULL;
v_ext_ADE[n] = NULL;
i_int_ADE[n] = NULL;
i_ext_ADE[n] = NULL;
}
v_int_ADE = NULL;
v_ext_ADE = NULL;
i_int_ADE = NULL;
i_ext_ADE = NULL;
}
Operator_Ext_LorentzMaterial::~Operator_Ext_LorentzMaterial()
{
for (int i=0;i<m_Order;++i)
{
for (int n=0; n<3; ++n)
{
delete[] v_int_ADE[n];
v_int_ADE[n] = NULL;
delete[] v_ext_ADE[n];
v_ext_ADE[n] = NULL;
delete[] i_int_ADE[n];
i_int_ADE[n] = NULL;
delete[] i_ext_ADE[n];
i_ext_ADE[n] = NULL;
if (m_volt_ADE_On[i])
{
delete[] v_int_ADE[i][n];
delete[] v_ext_ADE[i][n];
}
if (m_curr_ADE_On[i])
{
delete[] i_int_ADE[i][n];
delete[] i_ext_ADE[i][n];
}
}
if (m_volt_ADE_On[i])
{
delete[] v_int_ADE[i];
delete[] v_ext_ADE[i];
}
if (m_curr_ADE_On[i])
{
delete[] i_int_ADE[i];
delete[] i_ext_ADE[i];
}
}
delete[] v_int_ADE;
delete[] v_ext_ADE;
delete[] i_int_ADE;
delete[] i_ext_ADE;
v_int_ADE = NULL;
v_ext_ADE = NULL;
i_int_ADE = NULL;
i_ext_ADE = NULL;
}
bool Operator_Ext_LorentzMaterial::BuildExtension()
@ -61,6 +81,11 @@ bool Operator_Ext_LorentzMaterial::BuildExtension()
vector<double> i_int[3];
vector<double> i_ext[3];
vector<unsigned int> v_pos[3];
m_Order = 1;
m_volt_ADE_On = new bool[1];
m_volt_ADE_On[0]=false;
m_curr_ADE_On = new bool[1];
m_curr_ADE_On[0]=false;
for (pos[0]=0; pos[0]<numLines[0]; ++pos[0])
{
@ -87,11 +112,11 @@ bool Operator_Ext_LorentzMaterial::BuildExtension()
if (w_plasma>0)
{
b_pos_on = true;
m_volt_ADE_On = true;
m_volt_ADE_On[0] = true;
L_D[n] = 1/(w_plasma*w_plasma*m_Op->EC_C[n][index]);
}
t_relax = mat->GetEpsRelaxTimeWeighted(n,coord);
if ((t_relax>0) && m_volt_ADE_On)
if ((t_relax>0) && m_volt_ADE_On[0])
{
R_D[n] = L_D[n]/t_relax;
}
@ -114,11 +139,11 @@ bool Operator_Ext_LorentzMaterial::BuildExtension()
if (w_plasma>0)
{
b_pos_on = true;
m_curr_ADE_On = true;
m_curr_ADE_On[0] = true;
C_D[n] = 1/(w_plasma*w_plasma*m_Op->EC_L[n][index]);
}
t_relax = mat->GetMueRelaxTimeWeighted(n,coord);
if ((t_relax>0) && m_curr_ADE_On)
if ((t_relax>0) && m_curr_ADE_On[0])
{
G_D[n] = C_D[n]/t_relax;
}
@ -158,33 +183,46 @@ bool Operator_Ext_LorentzMaterial::BuildExtension()
}
//copy all vectors into the array's
m_LM_Count = v_pos[0].size();
m_LM_Count.push_back(v_pos[0].size());
m_LM_pos = new unsigned int**[1];
m_LM_pos[0] = new unsigned int*[3];
v_int_ADE = new FDTD_FLOAT**[1];
v_ext_ADE = new FDTD_FLOAT**[1];
i_int_ADE = new FDTD_FLOAT**[1];
i_ext_ADE = new FDTD_FLOAT**[1];
v_int_ADE[0] = new FDTD_FLOAT*[3];
v_ext_ADE[0] = new FDTD_FLOAT*[3];
i_int_ADE[0] = new FDTD_FLOAT*[3];
i_ext_ADE[0] = new FDTD_FLOAT*[3];
for (int n=0; n<3; ++n)
{
m_LM_pos[n] = new unsigned int[m_LM_Count];
for (unsigned int i=0; i<m_LM_Count; ++i)
m_LM_pos[n][i] = v_pos[n].at(i);
m_LM_pos[0][n] = new unsigned int[m_LM_Count.at(0)];
for (unsigned int i=0; i<m_LM_Count.at(0); ++i)
m_LM_pos[0][n][i] = v_pos[n].at(i);
if (m_volt_ADE_On)
{
v_int_ADE[n] = new FDTD_FLOAT[m_LM_Count];
v_ext_ADE[n] = new FDTD_FLOAT[m_LM_Count];
v_int_ADE[0][n] = new FDTD_FLOAT[m_LM_Count.at(0)];
v_ext_ADE[0][n] = new FDTD_FLOAT[m_LM_Count.at(0)];
for (unsigned int i=0; i<m_LM_Count; ++i)
for (unsigned int i=0; i<m_LM_Count.at(0); ++i)
{
v_int_ADE[n][i] = v_int[n].at(i);
v_ext_ADE[n][i] = v_ext[n].at(i);
v_int_ADE[0][n][i] = v_int[n].at(i);
v_ext_ADE[0][n][i] = v_ext[n].at(i);
}
}
if (m_curr_ADE_On)
{
i_int_ADE[n] = new FDTD_FLOAT[m_LM_Count];
i_ext_ADE[n] = new FDTD_FLOAT[m_LM_Count];
i_int_ADE[0][n] = new FDTD_FLOAT[m_LM_Count.at(0)];
i_ext_ADE[0][n] = new FDTD_FLOAT[m_LM_Count.at(0)];
for (unsigned int i=0; i<m_LM_Count; ++i)
for (unsigned int i=0; i<m_LM_Count.at(0); ++i)
{
i_int_ADE[n][i] = i_int[n].at(i);
i_ext_ADE[n][i] = i_ext[n].at(i);
i_int_ADE[0][n][i] = i_int[n].at(i);
i_ext_ADE[0][n][i] = i_ext[n].at(i);
}
}

11
FDTD/extensions/operator_ext_lorentzmaterial.h
View File

@ -34,16 +34,17 @@ public:
virtual bool IsCylinderCoordsSave(bool closedAlpha, bool R0_included) const {return false;}
virtual std::string GetExtensionName() const {return std::string("Lorentz Dispersive Material Extension");}
virtual std::string GetExtensionName() const {return std::string("Drude/Lorentz Dispersive Material Extension");}
virtual void ShowStat(ostream &ostr) const;
protected:
FDTD_FLOAT *v_int_ADE[3];
FDTD_FLOAT *v_ext_ADE[3];
FDTD_FLOAT *i_int_ADE[3];
FDTD_FLOAT *i_ext_ADE[3];
//ADE update coefficients, array setup: coeff[N_order][direction][mesh_pos_index]
FDTD_FLOAT ***v_int_ADE;
FDTD_FLOAT ***v_ext_ADE;
FDTD_FLOAT ***i_int_ADE;
FDTD_FLOAT ***i_ext_ADE;
};
#endif // OPERATOR_EXT_LORENTZMATERIAL_H