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

changed node-interpolation in H-field dumps 

average the H-fields not the currents
pull/1/head
Sebastian Held 2010-07-19 09:04:13 +02:00
parent 1a93650fa0
commit e3904c0f18
1 changed files with 30 additions and 20 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

50
FDTD/processfields_td.cpp
View File

@ -43,8 +43,6 @@ void ProcessFieldsTD::DumpNodeInterpol(string filename)
FDTD_FLOAT**** H_T = Create_N_3DArray<FDTD_FLOAT>(numLines); FDTD_FLOAT**** H_T = Create_N_3DArray<FDTD_FLOAT>(numLines);
unsigned int pos[3] = {start[0],start[1],start[2]}; unsigned int pos[3] = {start[0],start[1],start[2]};
unsigned int OpPos[3]; unsigned int OpPos[3];
double delta;
// cerr << "processing e-fields... " << endl;
for (pos[0]=0;pos[0]<numLines[0];++pos[0]) for (pos[0]=0;pos[0]<numLines[0];++pos[0])
{ {
OpPos[0]=start[0]+pos[0]*subSample[0]; OpPos[0]=start[0]+pos[0]*subSample[0];
@ -54,27 +52,39 @@ void ProcessFieldsTD::DumpNodeInterpol(string filename)
for (pos[2]=0;pos[2]<numLines[2];++pos[2]) for (pos[2]=0;pos[2]<numLines[2];++pos[2])
{ {
OpPos[2]=start[2]+pos[2]*subSample[2]; OpPos[2]=start[2]+pos[2]*subSample[2];
//in x //in x
delta = Op->GetMeshDelta(0,OpPos,true); H_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(0,OpPos) / Op->GetMeshDelta(0,OpPos,true);
if (delta) OpPos[1]++;
{ H_T[0][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(0,OpPos) / Op->GetMeshDelta(0,OpPos,true);
H_T[0][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(0,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(0,OpPos[0],OpPos[1]+1,OpPos[2]) + Eng->GetCurr(0,OpPos[0],OpPos[1],OpPos[2]+1) + Eng->GetCurr(0,OpPos[0],OpPos[1]+1,OpPos[2]+1); OpPos[2]++;
H_T[0][pos[0]][pos[1]][pos[2]] /= (4*delta); H_T[0][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(0,OpPos) / Op->GetMeshDelta(0,OpPos,true);
} OpPos[1]--;
H_T[0][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(0,OpPos) / Op->GetMeshDelta(0,OpPos,true);
OpPos[2]--;
H_T[0][pos[0]][pos[1]][pos[2]] /= 4.0;
//in y //in y
delta = Op->GetMeshDelta(1,OpPos,true); H_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(1,OpPos) / Op->GetMeshDelta(1,OpPos,true);
if (delta) OpPos[0]++;
{ H_T[1][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(1,OpPos) / Op->GetMeshDelta(1,OpPos,true);
H_T[1][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(1,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(1,OpPos[0]+1,OpPos[1],OpPos[2]) + Eng->GetCurr(1,OpPos[0],OpPos[1],OpPos[2]+1) + Eng->GetCurr(1,OpPos[0]+1,OpPos[1],OpPos[2]+1); OpPos[2]++;
H_T[1][pos[0]][pos[1]][pos[2]] /= (4*delta);//*Op->gridDelta); H_T[1][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(1,OpPos) / Op->GetMeshDelta(1,OpPos,true);
} OpPos[0]--;
H_T[1][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(1,OpPos) / Op->GetMeshDelta(1,OpPos,true);
OpPos[2]--;
H_T[1][pos[0]][pos[1]][pos[2]] /= 4.0;
//in z //in z
delta = Op->GetMeshDelta(2,OpPos,true); H_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(2,OpPos) / Op->GetMeshDelta(2,OpPos,true);
if (delta) OpPos[1]++;
{ H_T[2][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(2,OpPos) / Op->GetMeshDelta(2,OpPos,true);
H_T[2][pos[0]][pos[1]][pos[2]] = Eng->GetCurr(2,OpPos[0],OpPos[1],OpPos[2]) + Eng->GetCurr(2,OpPos[0],OpPos[1]+1,OpPos[2]) + Eng->GetCurr(2,OpPos[0]+1,OpPos[1],OpPos[2]) + Eng->GetCurr(2,OpPos[0]+1,OpPos[1]+1,OpPos[2]); OpPos[0]++;
H_T[2][pos[0]][pos[1]][pos[2]] /= (4*delta); H_T[2][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(2,OpPos) / Op->GetMeshDelta(2,OpPos,true);
} OpPos[1]--;
H_T[2][pos[0]][pos[1]][pos[2]] += Eng->GetCurr(2,OpPos) / Op->GetMeshDelta(2,OpPos,true);
OpPos[0]--;
H_T[2][pos[0]][pos[1]][pos[2]] /= 4.0;
} }
} }
} }