From 21ccdc7d0a990626b43fae2bb31366f400460861 Mon Sep 17 00:00:00 2001 From: Thorsten Liebig Date: Wed, 18 Jul 2012 16:59:42 +0200 Subject: [PATCH] TFSF: allow 1D/2D plane wave excitations Signed-off-by: Thorsten Liebig --- FDTD/extensions/engine_ext_tfsf.cpp | 149 +++++++++++------ FDTD/extensions/operator_ext_tfsf.cpp | 226 ++++++++++++++++---------- FDTD/extensions/operator_ext_tfsf.h | 1 + 3 files changed, 236 insertions(+), 140 deletions(-) diff --git a/FDTD/extensions/engine_ext_tfsf.cpp b/FDTD/extensions/engine_ext_tfsf.cpp index db142ba..418ffee 100644 --- a/FDTD/extensions/engine_ext_tfsf.cpp +++ b/FDTD/extensions/engine_ext_tfsf.cpp @@ -58,39 +58,63 @@ void Engine_Ext_TFSF::DoPostVoltageUpdates() { nP = (n+1)%3; nPP = (n+2)%3; + + // lower plane pos[nP] = m_Op_TFSF->m_Start[nP]; - ui_pos = 0; - for (unsigned int i=0;im_numLines[nP];++i) + if (m_Op_TFSF->m_ActiveDir[n][0]) { - pos[nPP] = m_Op_TFSF->m_Start[nPP]; - for (unsigned int j=0;jm_numLines[nPP];++j) + + for (unsigned int i=0;im_numLines[nP];++i) { - // current updates - pos[n] = m_Op_TFSF->m_Start[n]; + pos[nPP] = m_Op_TFSF->m_Start[nPP]; + for (unsigned int j=0;jm_numLines[nPP];++j) + { + // current updates + pos[n] = m_Op_TFSF->m_Start[n]; - m_Eng->SetVolt(nP,pos, m_Eng->GetVolt(nP,pos) - + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][0][0][ui_pos])*m_Op_TFSF->m_VoltAmp[n][0][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][0][0][ui_pos]]] - + m_Op_TFSF->m_VoltDelayDelta[n][0][0][ui_pos] *m_Op_TFSF->m_VoltAmp[n][0][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][0][0][ui_pos]]] ); + m_Eng->SetVolt(nP,pos, m_Eng->GetVolt(nP,pos) + + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][0][0][ui_pos])*m_Op_TFSF->m_VoltAmp[n][0][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][0][0][ui_pos]]] + + m_Op_TFSF->m_VoltDelayDelta[n][0][0][ui_pos] *m_Op_TFSF->m_VoltAmp[n][0][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][0][0][ui_pos]]] ); - m_Eng->SetVolt(nPP,pos, m_Eng->GetVolt(nPP,pos) - + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][0][1][ui_pos])*m_Op_TFSF->m_VoltAmp[n][0][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][0][1][ui_pos]]] - + m_Op_TFSF->m_VoltDelayDelta[n][0][1][ui_pos] *m_Op_TFSF->m_VoltAmp[n][0][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][0][1][ui_pos]]] ); + m_Eng->SetVolt(nPP,pos, m_Eng->GetVolt(nPP,pos) + + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][0][1][ui_pos])*m_Op_TFSF->m_VoltAmp[n][0][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][0][1][ui_pos]]] + + m_Op_TFSF->m_VoltDelayDelta[n][0][1][ui_pos] *m_Op_TFSF->m_VoltAmp[n][0][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][0][1][ui_pos]]] ); - pos[n] = m_Op_TFSF->m_Stop[n]; - - m_Eng->SetVolt(nP,pos, m_Eng->GetVolt(nP,pos) - + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][1][0][ui_pos])*m_Op_TFSF->m_VoltAmp[n][1][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][1][0][ui_pos]]] - + m_Op_TFSF->m_VoltDelayDelta[n][1][0][ui_pos] *m_Op_TFSF->m_VoltAmp[n][1][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][1][0][ui_pos]]] ); - - m_Eng->SetVolt(nPP,pos, m_Eng->GetVolt(nPP,pos) - + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][1][1][ui_pos])*m_Op_TFSF->m_VoltAmp[n][1][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][1][1][ui_pos]]] - + m_Op_TFSF->m_VoltDelayDelta[n][1][1][ui_pos] *m_Op_TFSF->m_VoltAmp[n][1][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][1][1][ui_pos]]] ); - - ++pos[nPP]; - ++ui_pos; + ++pos[nPP]; + ++ui_pos; + } + ++pos[nP]; + } + } + + // upper plane + pos[nP] = m_Op_TFSF->m_Start[nP]; + ui_pos = 0; + if (m_Op_TFSF->m_ActiveDir[n][1]) + { + + for (unsigned int i=0;im_numLines[nP];++i) + { + pos[nPP] = m_Op_TFSF->m_Start[nPP]; + for (unsigned int j=0;jm_numLines[nPP];++j) + { + // current updates + pos[n] = m_Op_TFSF->m_Stop[n]; + + m_Eng->SetVolt(nP,pos, m_Eng->GetVolt(nP,pos) + + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][1][0][ui_pos])*m_Op_TFSF->m_VoltAmp[n][1][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][1][0][ui_pos]]] + + m_Op_TFSF->m_VoltDelayDelta[n][1][0][ui_pos] *m_Op_TFSF->m_VoltAmp[n][1][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][1][0][ui_pos]]] ); + + m_Eng->SetVolt(nPP,pos, m_Eng->GetVolt(nPP,pos) + + (1.0-m_Op_TFSF->m_VoltDelayDelta[n][1][1][ui_pos])*m_Op_TFSF->m_VoltAmp[n][1][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_VoltDelay[n][1][1][ui_pos]]] + + m_Op_TFSF->m_VoltDelayDelta[n][1][1][ui_pos] *m_Op_TFSF->m_VoltAmp[n][1][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_VoltDelay[n][1][1][ui_pos]]] ); + + ++pos[nPP]; + ++ui_pos; + } + ++pos[nP]; } - ++pos[nP]; } } } @@ -118,41 +142,66 @@ void Engine_Ext_TFSF::DoPostCurrentUpdates() unsigned int pos[3]; for (int n=0;n<3;++n) { + if (!m_Op_TFSF->m_ActiveDir[n][0] && !m_Op_TFSF->m_ActiveDir[n][1]) + continue; + nP = (n+1)%3; nPP = (n+2)%3; + + // lower plane pos[nP] = m_Op_TFSF->m_Start[nP]; - ui_pos = 0; - for (unsigned int i=0;im_numLines[nP];++i) + if (m_Op_TFSF->m_ActiveDir[n][0]) { - pos[nPP] = m_Op_TFSF->m_Start[nPP]; - for (unsigned int j=0;jm_numLines[nPP];++j) + for (unsigned int i=0;im_numLines[nP];++i) { - // current updates - pos[n] = m_Op_TFSF->m_Start[n]-1; + pos[nPP] = m_Op_TFSF->m_Start[nPP]; + for (unsigned int j=0;jm_numLines[nPP];++j) + { + // current updates + pos[n] = m_Op_TFSF->m_Start[n]-1; - m_Eng->SetCurr(nP,pos, m_Eng->GetCurr(nP,pos) - + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][0][0][ui_pos])*m_Op_TFSF->m_CurrAmp[n][0][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][0][0][ui_pos]]] - + m_Op_TFSF->m_CurrDelayDelta[n][0][0][ui_pos] *m_Op_TFSF->m_CurrAmp[n][0][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][0][0][ui_pos]]] ); + m_Eng->SetCurr(nP,pos, m_Eng->GetCurr(nP,pos) + + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][0][0][ui_pos])*m_Op_TFSF->m_CurrAmp[n][0][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][0][0][ui_pos]]] + + m_Op_TFSF->m_CurrDelayDelta[n][0][0][ui_pos] *m_Op_TFSF->m_CurrAmp[n][0][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][0][0][ui_pos]]] ); - m_Eng->SetCurr(nPP,pos, m_Eng->GetCurr(nPP,pos) - + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][0][1][ui_pos])*m_Op_TFSF->m_CurrAmp[n][0][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][0][1][ui_pos]]] - + m_Op_TFSF->m_CurrDelayDelta[n][0][1][ui_pos] *m_Op_TFSF->m_CurrAmp[n][0][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][0][1][ui_pos]]] ); + m_Eng->SetCurr(nPP,pos, m_Eng->GetCurr(nPP,pos) + + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][0][1][ui_pos])*m_Op_TFSF->m_CurrAmp[n][0][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][0][1][ui_pos]]] + + m_Op_TFSF->m_CurrDelayDelta[n][0][1][ui_pos] *m_Op_TFSF->m_CurrAmp[n][0][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][0][1][ui_pos]]] ); - pos[n] = m_Op_TFSF->m_Stop[n]; - - m_Eng->SetCurr(nP,pos, m_Eng->GetCurr(nP,pos) - + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][1][0][ui_pos])*m_Op_TFSF->m_CurrAmp[n][1][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][1][0][ui_pos]]] - + m_Op_TFSF->m_CurrDelayDelta[n][1][0][ui_pos] *m_Op_TFSF->m_CurrAmp[n][1][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][1][0][ui_pos]]] ); - - m_Eng->SetCurr(nPP,pos, m_Eng->GetCurr(nPP,pos) - + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][1][1][ui_pos])*m_Op_TFSF->m_CurrAmp[n][1][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][1][1][ui_pos]]] - + m_Op_TFSF->m_CurrDelayDelta[n][1][1][ui_pos] *m_Op_TFSF->m_CurrAmp[n][1][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][1][1][ui_pos]]] ); - - ++pos[nPP]; - ++ui_pos; + ++pos[nPP]; + ++ui_pos; + } + ++pos[nP]; + } + } + + // upper plane + pos[nP] = m_Op_TFSF->m_Start[nP]; + ui_pos = 0; + if (m_Op_TFSF->m_ActiveDir[n][1]) + { + for (unsigned int i=0;im_numLines[nP];++i) + { + pos[nPP] = m_Op_TFSF->m_Start[nPP]; + for (unsigned int j=0;jm_numLines[nPP];++j) + { + // current updates + pos[n] = m_Op_TFSF->m_Stop[n]; + + m_Eng->SetCurr(nP,pos, m_Eng->GetCurr(nP,pos) + + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][1][0][ui_pos])*m_Op_TFSF->m_CurrAmp[n][1][0][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][1][0][ui_pos]]] + + m_Op_TFSF->m_CurrDelayDelta[n][1][0][ui_pos] *m_Op_TFSF->m_CurrAmp[n][1][0][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][1][0][ui_pos]]] ); + + m_Eng->SetCurr(nPP,pos, m_Eng->GetCurr(nPP,pos) + + (1.0-m_Op_TFSF->m_CurrDelayDelta[n][1][1][ui_pos])*m_Op_TFSF->m_CurrAmp[n][1][1][ui_pos]*signal[m_DelayLookup[ m_Op_TFSF->m_CurrDelay[n][1][1][ui_pos]]] + + m_Op_TFSF->m_CurrDelayDelta[n][1][1][ui_pos] *m_Op_TFSF->m_CurrAmp[n][1][1][ui_pos]*signal[m_DelayLookup[1+m_Op_TFSF->m_CurrDelay[n][1][1][ui_pos]]] ); + + ++pos[nPP]; + ++ui_pos; + } + ++pos[nP]; } - ++pos[nP]; } } } diff --git a/FDTD/extensions/operator_ext_tfsf.cpp b/FDTD/extensions/operator_ext_tfsf.cpp index ea3aa0f..5b1965d 100644 --- a/FDTD/extensions/operator_ext_tfsf.cpp +++ b/FDTD/extensions/operator_ext_tfsf.cpp @@ -152,6 +152,15 @@ bool Operator_Ext_TFST::BuildExtension() m_numLines[n] = m_Stop[n]-m_Start[n]+1; m_IncLow[n] = m_PropDir[n]>=0; + if (m_Start[n]==0) + m_ActiveDir[n][0]=false; + else + m_ActiveDir[n][0]=true; + if (m_Stop[n]==m_Op->GetOriginalNumLines(n)-1) + m_ActiveDir[n][1]=false; + else + m_ActiveDir[n][1]=true; + if (m_IncLow[n]) { ui_origin[n] = m_Start[n]; @@ -205,16 +214,22 @@ bool Operator_Ext_TFST::BuildExtension() numP = m_numLines[nP]*m_numLines[nPP]; + if (!m_ActiveDir[n][0] && !m_ActiveDir[n][1]) + continue; + for (int l=0;l<2;++l) for (int c=0;c<2;++c) { - m_VoltDelay[n][l][c]=new unsigned int[numP]; - m_VoltDelayDelta[n][l][c]=new FDTD_FLOAT[numP]; - m_VoltAmp[n][l][c]=new FDTD_FLOAT[numP]; + if (m_ActiveDir[n][l]) + { + m_VoltDelay[n][l][c]=new unsigned int[numP]; + m_VoltDelayDelta[n][l][c]=new FDTD_FLOAT[numP]; + m_VoltAmp[n][l][c]=new FDTD_FLOAT[numP]; - m_CurrDelay[n][l][c]=new unsigned int[numP]; - m_CurrDelayDelta[n][l][c]=new FDTD_FLOAT[numP]; - m_CurrAmp[n][l][c]=new FDTD_FLOAT[numP]; + m_CurrDelay[n][l][c]=new unsigned int[numP]; + m_CurrDelayDelta[n][l][c]=new FDTD_FLOAT[numP]; + m_CurrAmp[n][l][c]=new FDTD_FLOAT[numP]; + } } ui_pos = 0; @@ -225,130 +240,159 @@ bool Operator_Ext_TFST::BuildExtension() { // current updates pos[n] = m_Start[n]; - m_Op->GetYeeCoords(nP,pos,coord,false); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_CurrDelay[n][0][1][ui_pos] = floor(delay); - m_CurrDelayDelta[n][0][1][ui_pos] = delay - floor(delay); - m_CurrAmp[n][0][1][ui_pos] = m_E_Amp[nP]*m_Op->GetEdgeLength(nP,pos); - m_Op->GetYeeCoords(nPP,pos,coord,false); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_CurrDelay[n][0][0][ui_pos] = floor(delay); - m_CurrDelayDelta[n][0][0][ui_pos] = delay - floor(delay); - m_CurrAmp[n][0][0][ui_pos] = m_E_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos); + if (m_ActiveDir[n][0]) + { + m_Op->GetYeeCoords(nP,pos,coord,false); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_CurrDelay[n][0][1][ui_pos] = floor(delay); + m_CurrDelayDelta[n][0][1][ui_pos] = delay - floor(delay); + m_CurrAmp[n][0][1][ui_pos] = m_E_Amp[nP]*m_Op->GetEdgeLength(nP,pos); - --pos[n]; - m_CurrAmp[n][0][0][ui_pos]*=m_Op->GetIV(nP,pos); - m_CurrAmp[n][0][1][ui_pos]*=m_Op->GetIV(nPP,pos); + m_Op->GetYeeCoords(nPP,pos,coord,false); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_CurrDelay[n][0][0][ui_pos] = floor(delay); + m_CurrDelayDelta[n][0][0][ui_pos] = delay - floor(delay); + m_CurrAmp[n][0][0][ui_pos] = m_E_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos); - pos[n] = m_Stop[n]; - m_Op->GetYeeCoords(nP,pos,coord,false); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_CurrDelay[n][1][1][ui_pos] = floor(delay); - m_CurrDelayDelta[n][1][1][ui_pos] = delay - floor(delay); - m_CurrAmp[n][1][1][ui_pos] = m_E_Amp[nP]*m_Op->GetEdgeLength(nP,pos); + --pos[n]; + m_CurrAmp[n][0][0][ui_pos]*=m_Op->GetIV(nP,pos); + m_CurrAmp[n][0][1][ui_pos]*=m_Op->GetIV(nPP,pos); + } - m_Op->GetYeeCoords(nPP,pos,coord,false); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_CurrDelay[n][1][0][ui_pos] = floor(delay); - m_CurrDelayDelta[n][1][0][ui_pos] = delay - floor(delay); - m_CurrAmp[n][1][0][ui_pos] = m_E_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos); + if (m_ActiveDir[n][1]) + { + pos[n] = m_Stop[n]; + m_Op->GetYeeCoords(nP,pos,coord,false); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_CurrDelay[n][1][1][ui_pos] = floor(delay); + m_CurrDelayDelta[n][1][1][ui_pos] = delay - floor(delay); + m_CurrAmp[n][1][1][ui_pos] = m_E_Amp[nP]*m_Op->GetEdgeLength(nP,pos); - m_CurrAmp[n][1][0][ui_pos]*=m_Op->GetIV(nP,pos); - m_CurrAmp[n][1][1][ui_pos]*=m_Op->GetIV(nPP,pos); + m_Op->GetYeeCoords(nPP,pos,coord,false); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_CurrDelay[n][1][0][ui_pos] = floor(delay); + m_CurrDelayDelta[n][1][0][ui_pos] = delay - floor(delay); + m_CurrAmp[n][1][0][ui_pos] = m_E_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos); + + m_CurrAmp[n][1][0][ui_pos]*=m_Op->GetIV(nP,pos); + m_CurrAmp[n][1][1][ui_pos]*=m_Op->GetIV(nPP,pos); + } if (m_IncLow[n]) { - m_CurrAmp[n][0][0][ui_pos]*=-1; - m_CurrAmp[n][1][1][ui_pos]*=-1; + if (m_ActiveDir[n][0]) + m_CurrAmp[n][0][0][ui_pos]*=-1; + if (m_ActiveDir[n][1]) + m_CurrAmp[n][1][1][ui_pos]*=-1; } else { - m_CurrAmp[n][0][1][ui_pos]*=-1; - m_CurrAmp[n][1][0][ui_pos]*=-1; + if (m_ActiveDir[n][0]) + m_CurrAmp[n][0][1][ui_pos]*=-1; + if (m_ActiveDir[n][1]) + m_CurrAmp[n][1][0][ui_pos]*=-1; } if (pos[nP]==m_Stop[nP]) { - m_CurrAmp[n][0][1][ui_pos]=0; - m_CurrAmp[n][1][1][ui_pos]=0; + if (m_ActiveDir[n][0]) + m_CurrAmp[n][0][1][ui_pos]=0; + if (m_ActiveDir[n][1]) + m_CurrAmp[n][1][1][ui_pos]=0; } if (pos[nPP]==m_Stop[nPP]) { - m_CurrAmp[n][0][0][ui_pos]=0; - m_CurrAmp[n][1][0][ui_pos]=0; + if (m_ActiveDir[n][0]) + m_CurrAmp[n][0][0][ui_pos]=0; + if (m_ActiveDir[n][1]) + m_CurrAmp[n][1][0][ui_pos]=0; } // voltage updates pos[n] = m_Start[n]-1; - m_Op->GetYeeCoords(nP,pos,coord,true); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_VoltDelay[n][0][1][ui_pos] = floor(delay); - m_VoltDelayDelta[n][0][1][ui_pos] = delay - floor(delay); - m_VoltAmp[n][0][1][ui_pos] = m_H_Amp[nP]*m_Op->GetEdgeLength(nP,pos,true); + if (m_ActiveDir[n][0]) + { + m_Op->GetYeeCoords(nP,pos,coord,true); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_VoltDelay[n][0][1][ui_pos] = floor(delay); + m_VoltDelayDelta[n][0][1][ui_pos] = delay - floor(delay); + m_VoltAmp[n][0][1][ui_pos] = m_H_Amp[nP]*m_Op->GetEdgeLength(nP,pos,true); - m_Op->GetYeeCoords(nPP,pos,coord,true); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_VoltDelay[n][0][0][ui_pos] = floor(delay); - m_VoltDelayDelta[n][0][0][ui_pos] = delay - floor(delay); - m_VoltAmp[n][0][0][ui_pos] = m_H_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos,true); + m_Op->GetYeeCoords(nPP,pos,coord,true); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_VoltDelay[n][0][0][ui_pos] = floor(delay); + m_VoltDelayDelta[n][0][0][ui_pos] = delay - floor(delay); + m_VoltAmp[n][0][0][ui_pos] = m_H_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos,true); - ++pos[n]; - m_VoltAmp[n][0][0][ui_pos]*=m_Op->GetVI(nP,pos); - m_VoltAmp[n][0][1][ui_pos]*=m_Op->GetVI(nPP,pos); + ++pos[n]; + m_VoltAmp[n][0][0][ui_pos]*=m_Op->GetVI(nP,pos); + m_VoltAmp[n][0][1][ui_pos]*=m_Op->GetVI(nPP,pos); + } pos[n] = m_Stop[n]; - m_Op->GetYeeCoords(nP,pos,coord,true); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_VoltDelay[n][1][1][ui_pos] = floor(delay); - m_VoltDelayDelta[n][1][1][ui_pos] = delay - floor(delay); - m_VoltAmp[n][1][1][ui_pos] = m_H_Amp[nP]*m_Op->GetEdgeLength(nP,pos,true); + if (m_ActiveDir[n][1]) + { + m_Op->GetYeeCoords(nP,pos,coord,true); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_VoltDelay[n][1][1][ui_pos] = floor(delay); + m_VoltDelayDelta[n][1][1][ui_pos] = delay - floor(delay); + m_VoltAmp[n][1][1][ui_pos] = m_H_Amp[nP]*m_Op->GetEdgeLength(nP,pos,true); - m_Op->GetYeeCoords(nPP,pos,coord,true); - dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; - delay = dist*unit/__C0__/dT; - m_maxDelay = max((unsigned int)delay,m_maxDelay); - m_VoltDelay[n][1][0][ui_pos] = floor(delay); - m_VoltDelayDelta[n][1][0][ui_pos] = delay - floor(delay); - m_VoltAmp[n][1][0][ui_pos] = m_H_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos,true); + m_Op->GetYeeCoords(nPP,pos,coord,true); + dist = fabs(coord[0]-origin[0])*m_PropDir[0]+fabs(coord[1]-origin[1])*m_PropDir[1]+fabs(coord[2]-origin[2])*m_PropDir[2]; + delay = dist*unit/__C0__/dT; + m_maxDelay = max((unsigned int)delay,m_maxDelay); + m_VoltDelay[n][1][0][ui_pos] = floor(delay); + m_VoltDelayDelta[n][1][0][ui_pos] = delay - floor(delay); + m_VoltAmp[n][1][0][ui_pos] = m_H_Amp[nPP]*m_Op->GetEdgeLength(nPP,pos,true); - m_VoltAmp[n][1][0][ui_pos]*=m_Op->GetVI(nP,pos); - m_VoltAmp[n][1][1][ui_pos]*=m_Op->GetVI(nPP,pos); + m_VoltAmp[n][1][0][ui_pos]*=m_Op->GetVI(nP,pos); + m_VoltAmp[n][1][1][ui_pos]*=m_Op->GetVI(nPP,pos); + } if (m_IncLow[n]) { - m_VoltAmp[n][1][0][ui_pos]*=-1; - m_VoltAmp[n][0][1][ui_pos]*=-1; + if (m_ActiveDir[n][1]) + m_VoltAmp[n][1][0][ui_pos]*=-1; + if (m_ActiveDir[n][0]) + m_VoltAmp[n][0][1][ui_pos]*=-1; } else { - m_VoltAmp[n][0][0][ui_pos]*=-1; - m_VoltAmp[n][1][1][ui_pos]*=-1; + if (m_ActiveDir[n][0]) + m_VoltAmp[n][0][0][ui_pos]*=-1; + if (m_ActiveDir[n][1]) + m_VoltAmp[n][1][1][ui_pos]*=-1; } if (pos[nP]==m_Stop[nP]) { - m_VoltAmp[n][0][0][ui_pos]=0; - m_VoltAmp[n][1][0][ui_pos]=0; + if (m_ActiveDir[n][0]) + m_VoltAmp[n][0][0][ui_pos]=0; + if (m_ActiveDir[n][1]) + m_VoltAmp[n][1][0][ui_pos]=0; } if (pos[nPP]==m_Stop[nPP]) { - m_VoltAmp[n][0][1][ui_pos]=0; - m_VoltAmp[n][1][1][ui_pos]=0; + if (m_ActiveDir[n][0]) + m_VoltAmp[n][0][1][ui_pos]=0; + if (m_ActiveDir[n][1]) + m_VoltAmp[n][1][1][ui_pos]=0; } ++pos[nPP]; @@ -372,10 +416,12 @@ Engine_Extension* Operator_Ext_TFST::CreateEngineExtention() void Operator_Ext_TFST::ShowStat(ostream &ostr) const { Operator_Extension::ShowStat(ostr); + cout << "Active directions\t: " << "(" << m_ActiveDir[0][0] << "/" << m_ActiveDir[0][1] << ", " << m_ActiveDir[1][0] << "/" << m_ActiveDir[1][1] << ", " << m_ActiveDir[2][0] << "/" << m_ActiveDir[2][1] << ")" << endl; cout << "Propagation direction\t: " << "(" << m_PropDir[0] << ", " << m_PropDir[1] << ", " << m_PropDir[2] << ")" << endl; cout << "E-field amplitude\t: " << "(" << m_E_Amp[0] << ", " << m_E_Amp[1] << ", " << m_E_Amp[2] << ")" << endl; cout << "m_H_Amp amplitude\t: " << "(" << m_H_Amp[0]*__Z0__ << ", " << m_H_Amp[1]*__Z0__ << ", " << m_H_Amp[2]*__Z0__ << ")/Z0" << endl; cout << "Box Dimensions\t\t: " << m_numLines[0] << " x " << m_numLines[1] << " x " << m_numLines[2] << endl; cout << "Max. Delay (TS)\t\t: " << m_maxDelay << endl; - cout << "Memory usage (est.)\t: ~" << m_numLines[0] * m_numLines[1] * m_numLines[2] * 6 * 4 * 4 / 1024 << " kiB" << endl; + int dirs = m_ActiveDir[0][0] + m_ActiveDir[0][1] + m_ActiveDir[1][0] + m_ActiveDir[1][1] + m_ActiveDir[2][0] + m_ActiveDir[2][1] ; + cout << "Memory usage (est.)\t: ~" << m_numLines[0] * m_numLines[1] * m_numLines[2] * dirs * 4 * 4 / 1024 << " kiB" << endl; } diff --git a/FDTD/extensions/operator_ext_tfsf.h b/FDTD/extensions/operator_ext_tfsf.h index 4f42b64..c95f33b 100644 --- a/FDTD/extensions/operator_ext_tfsf.h +++ b/FDTD/extensions/operator_ext_tfsf.h @@ -51,6 +51,7 @@ protected: Excitation* m_Exc; bool m_IncLow[3]; + bool m_ActiveDir[3][2]; // m_ActiveDir[direction][low/high] unsigned int m_Start[3]; unsigned int m_Stop[3]; unsigned int m_numLines[3];