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

added Curve for excitation and dump excitation in operator-dump

pull/1/head
Thorsten Liebig 2010-03-27 23:05:45 +01:00
parent 98f72a855e
commit ea2823377d
3 changed files with 109 additions and 21 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

93
FDTD/operator.cpp
View File

@ -106,7 +106,7 @@ bool Operator::SnapToMesh(double* dcoord, unsigned int* uicoord, bool lower)
else if (dcoord[n]>discLines[n][numLines[n]-1]) {ok=false;uicoord[n]=numLines[n]-1; if (lower) uicoord[n]=numLines[n]-2;}
else if (dcoord[n]==discLines[n][numLines[n]-1]) {uicoord[n]=numLines[n]-1; if (lower) uicoord[n]=numLines[n]-2;}
else
for (unsigned int i=1;i<numLines[n]-1;++i)
for (unsigned int i=1;i<numLines[n];++i)
{
if (dcoord[n]<discLines[n][i])
{
@ -224,9 +224,9 @@ void Operator::ShowSize()
cout << "-----------------------------" << endl;
}
void Operator::CalcGaussianPulsExcitation(double f0, double fc)
bool Operator::CalcGaussianPulsExcitation(double f0, double fc)
{
if (dT==0) return;
if (dT==0) return false;
ExciteLength = (unsigned int)(2.0 * 9.0/(2.0*PI*fc) / dT);
cerr << "Operator::CalcGaussianPulsExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl;
@ -238,14 +238,16 @@ void Operator::CalcGaussianPulsExcitation(double f0, double fc)
ExciteSignal[n] = cos(2.0*PI*f0*(n*dT-9.0/(2.0*PI*fc)))*exp(-1*pow(2.0*PI*fc*n*dT/3.0-3,2));
// cerr << ExciteSignal[n] << endl;
}
return true;
}
void Operator::CalcSinusExcitation(double f0, int nTS)
bool Operator::CalcSinusExcitation(double f0, int nTS)
{
if (dT==0) return;
if (nTS<=0) return;
if (dT==0) return false;
if (nTS<=0) return false;
ExciteLength = (unsigned int)(nTS);
cerr << "Operator::CalcSinusExcitation: Length of the excite signal: " << ExciteLength << " timesteps" << endl;
delete[] ExciteSignal;
ExciteSignal = new FDTD_FLOAT[ExciteLength+1];
ExciteSignal[0]=0.0;
@ -254,6 +256,7 @@ void Operator::CalcSinusExcitation(double f0, int nTS)
ExciteSignal[n] = sin(2.0*PI*f0*n*dT);
// cerr << ExciteSignal[n] << endl;
}
return true;
}
void Operator::DumpOperator2File(string filename)
@ -266,10 +269,18 @@ void Operator::DumpOperator2File(string filename)
return;
}
string names[] = {"vv", "vi", "iv" , "ii"};
FDTD_FLOAT**** array[] = {vv,vi,iv,ii};
FDTD_FLOAT**** exc = Create_N_3DArray(numLines);
for (unsigned int n=0;n<E_Exc_Count;++n)
{
exc[E_Exc_dir[n]][E_Exc_index[0][n]][E_Exc_index[1][n]][E_Exc_index[2][n]] = E_Exc_amp[n];
}
ProcessFields::DumpMultiVectorArray2VTK(file, names , array , 4, discLines, numLines);
string names[] = {"vv", "vi", "iv" , "ii", "exc"};
FDTD_FLOAT**** array[] = {vv,vi,iv,ii,exc};
ProcessFields::DumpMultiVectorArray2VTK(file, names , array , 5, discLines, numLines);
Delete_N_3DArray(exc,numLines);
file.close();
}
@ -774,6 +785,70 @@ bool Operator::CalcEFieldExcitation()
}
}
//special treatment for primitives of type curve (treated as wires) see also Calc_PEC
double p1[3];
double p2[3];
double deltaN=0.0;
int n;
struct Grid_Path path;
CSPropElectrode* elec=NULL;
CSProperties* prop=NULL;
vector<CSProperties*> vec_prop = CSX->GetPropertyByType(CSProperties::ELECTRODE);
for (size_t p=0;p<vec_prop.size();++p)
{
prop = vec_prop.at(p);
elec = prop->ToElectrode();
for (size_t n=0;n<prop->GetQtyPrimitives();++n)
{
CSPrimitives* prim = prop->GetPrimitive(n);
CSPrimCurve* curv = prim->ToCurve();
if (curv)
{
for (size_t i=1;i<curv->GetNumberOfPoints();++i)
{
curv->GetPoint(i-1,p1);
curv->GetPoint(i,p2);
path = FindPath(p1,p2);
for (size_t t=0;t<path.dir.size();++t)
{
n = path.dir.at(t);
pos[0] = path.posPath[0].at(t);
pos[1] = path.posPath[1].at(t);
pos[2] = path.posPath[2].at(t);
MainOp->SetPos(pos[0],pos[1],pos[2]);
deltaN=fabs(MainOp->GetIndexDelta(n,pos[n]));
coord[0] = discLines[0][pos[0]];
coord[1] = discLines[1][pos[1]];
coord[2] = discLines[2][pos[2]];
coord[n] += 0.5*deltaN;
// cerr << n << " " << coord[0] << " " << coord[1] << " " << coord[2] << endl;
if (elec!=NULL)
{
if ((elec->GetActiveDir(n)) && (pos[n]<(int)numLines[n]-1))
{
amp = elec->GetWeightedExcitation(n,coord)*deltaN*gridDelta;
if (amp!=0)
{
vExcit.push_back(amp);
vDelay.push_back((unsigned int)(elec->GetDelay()/dT));
vDir.push_back(n);
vIndex[0].push_back(pos[0]);
vIndex[1].push_back(pos[1]);
vIndex[2].push_back(pos[2]);
}
if (elec->GetExcitType()==1) //hard excite
{
vv[n][pos[0]][pos[1]][pos[2]] = 0;
vi[n][pos[0]][pos[1]][pos[2]] = 0;
}
}
}
}
}
}
}
}
E_Exc_Count = vExcit.size();
cerr << "Operator::CalcEFieldExcitation: Found number of excitations points: " << E_Exc_Count << endl;
if (E_Exc_Count==0)

4
FDTD/operator.h
View File

@ -38,9 +38,9 @@ public:
virtual int CalcECOperator();
//! Calculate an excitation with center of f0 and the half bandwidth fc
virtual void CalcGaussianPulsExcitation(double f0, double fc);
virtual bool CalcGaussianPulsExcitation(double f0, double fc);
//! Calculate a sinusoidal excitation with frequency f0 and a duration of nTS number of timesteps
virtual void CalcSinusExcitation(double f0, int nTS);
virtual bool CalcSinusExcitation(double f0, int nTS);
virtual void ApplyElectricBC(bool* dirs); //applied by default to all boundaries
virtual void ApplyMagneticBC(bool* dirs);

33
openems.cpp
View File

@ -185,26 +185,39 @@ int openEMS::SetupFDTD(const char* file)
FDTD_Op = new Operator();
if (FDTD_Op->SetGeometryCSX(&CSX)==false) return(-1);
if (DebugMat)
{
FDTD_Op->DumpMaterial2File("material_dump.vtk");
}
FDTD_Op->CalcECOperator();
if (DebugOp)
{
FDTD_Op->DumpOperator2File("operator_dump.vtk");
}
if (Excit_Type==0)
FDTD_Op->CalcGaussianPulsExcitation(f0,fc);
{
if (!FDTD_Op->CalcGaussianPulsExcitation(f0,fc))
{
cerr << "openEMS: excitation setup failed!!" << endl;
exit(2);
}
}
else if (Excit_Type==1)
FDTD_Op->CalcSinusExcitation(f0,NrTS);
{
if (!FDTD_Op->CalcSinusExcitation(f0,NrTS))
{
cerr << "openEMS: excitation setup failed!!" << endl;
exit(2);
}
}
else
{
cerr << "openEMS: Excitation type is unknown" << endl;
exit(-1);
}
if (DebugMat)
{
FDTD_Op->DumpMaterial2File("material_dump.vtk");
}
if (DebugOp)
{
FDTD_Op->DumpOperator2File("operator_dump.vtk");
}
time_t OpDoneTime=time(NULL);
FDTD_Op->ShowSize();