/*** Translated to the C language by N. Kyriazis 20 Aug 2003 ***
Program NEC(input,tape5=input,output,tape11,tape12,tape13,tape14,
tape15,tape16,tape20,tape21)
Numerical Electromagnetics Code (NEC2) developed at Lawrence
Livermore lab., Livermore, CA. (contact G. Burke at 415-422-8414
for problems with the NEC code. For problems with the vax implem-
entation, contact J. Breakall at 415-422-8196 or E. Domning at 415
422-5936)
file created 4/11/80.
***********Notice**********
This computer code material was prepared as an account of work
sponsored by the United States government. Neither the United
States nor the United States Department Of Energy, nor any of
their employees, nor any of their contractors, subcontractors,
or their employees, makes any warranty, express or implied, or
assumes any legal liability or responsibility for the accuracy,
completeness or usefulness of any information, apparatus, product
or process disclosed, or represents that its use would not infringe
privately-owned rights.
*******************************************************************/
#include "nec2c.h"
/* common /netcx/ */
extern netcx_t netcx;
/* common /vsorc/ */
extern vsorc_t vsorc;
/* common /data/ */
extern data_t data;
/* common /crnt/ */
extern crnt_t crnt;
/* pointers to input/output files */
extern FILE *input_fp, *output_fp, *plot_fp;
/*-------------------------------------------------------------------*/
/* subroutine netwk solves for structure currents for a given */
/* excitation including the effect of non-radiating networks if */
/* present. */
void netwk (complex double *cm, int *ip, complex double *einc)
{
int *ipnt = NULL, *nteqa = NULL, *ntsca = NULL;
int jump1, jump2, nteq = 0, ntsc = 0, nseg2, irow2 = 0, j, ndimn;
int neqz2, neqt, irow1 = 0, i, nseg1, isc1 = 0, isc2 = 0;
double asmx, asa, pwr, y11r, y11i, y12r, y12i, y22r, y22i;
complex double *vsrc = NULL, *rhs = NULL, *cmn = NULL;
complex double *rhnt = NULL, *rhnx = NULL, ymit, vlt, cux;
neqz2 = netcx.neq2;
if (neqz2 == 0)
neqz2 = 1;
netcx.pin = 0.;
netcx.pnls = 0.;
neqt = netcx.neq + netcx.neq2;
ndimn = j = (2 * netcx.nonet + vsorc.nsant);
/* Allocate network buffers */
if (netcx.nonet != 0)
{
mem_alloc ((void *) &rhs, data.np3m * sizeof (complex double));
i = j * sizeof (complex double);
mem_alloc ((void *) &rhnt, i);
mem_alloc ((void *) &rhnx, i);
mem_alloc ((void *) &cmn, i * j);
i = j * sizeof (int);
mem_alloc ((void *) &ntsca, i);
mem_alloc ((void *) &nteqa, i);
mem_alloc ((void *) &ipnt, i);
mem_alloc ((void *) &vsrc, vsorc.nsant * sizeof (complex double));
}
else if (netcx.masym != 0)
{
i = j * sizeof (int);
mem_alloc ((void *) &ipnt, i);
}
if (netcx.ntsol == 0)
{
/* compute relative matrix asymmetry */
if (netcx.masym != 0)
{
irow1 = 0;
if (netcx.nonet != 0)
{
for (i = 0; i < netcx.nonet; i++)
{
nseg1 = netcx.iseg1[i];
for (isc1 = 0; isc1 < 2; isc1++)
{
if (irow1 == 0)
{
ipnt[irow1] = nseg1;
nseg1 = netcx.iseg2[i];
irow1++;
continue;
}
for (j = 0; j < irow1; j++)
if (nseg1 == ipnt[j])
break;
if (j == irow1)
{
ipnt[irow1] = nseg1;
irow1++;
}
nseg1 = netcx.iseg2[i];
} /* for( isc1 = 0; isc1 < 2; isc1++ ) */
} /* for( i = 0; i < netcx.nonet; i++ ) */
} /* if( netcx.nonet != 0) */
if (vsorc.nsant != 0)
{
for (i = 0; i < vsorc.nsant; i++)
{
nseg1 = vsorc.isant[i];
if (irow1 == 0)
{
ipnt[irow1] = nseg1;
irow1++;
continue;
}
for (j = 0; j < irow1; j++)
if (nseg1 == ipnt[j])
break;
if (j == irow1)
{
ipnt[irow1] = nseg1;
irow1++;
}
} /* for( i = 0; i < vsorc.nsant; i++ ) */
} /* if( vsorc.nsant != 0) */
if (irow1 >= 2)
{
for (i = 0; i < irow1; i++)
{
isc1 = ipnt[i] - 1;
asmx = data.si[isc1];
for (j = 0; j < neqt; j++)
rhs[j] = CPLX_00;
rhs[isc1] = CPLX_10;
solves (
cm,
ip,
rhs,
netcx.neq,
1,
data.np,
data.n,
data.mp,
data.m);
cabc (rhs);
for (j = 0; j < irow1; j++)
{
isc1 = ipnt[j] - 1;
cmn[j + i * ndimn] = rhs[isc1] / asmx;
}
} /* for( i = 0; i < irow1; i++ ) */
asmx = 0.;
asa = 0.;
for (i = 1; i < irow1; i++)
{
isc1 = i;
for (j = 0; j < isc1; j++)
{
cux = cmn[i + j * ndimn];
pwr = cabsl ((cux - cmn[j + i * ndimn]) / cux);
asa += pwr * pwr;
if (pwr < asmx)
continue;
asmx = pwr;
nteq = ipnt[i];
ntsc = ipnt[j];
} /* for( j = 0; j < isc1; j++ ) */
} /* for( i = 1; i < irow1; i++ ) */
asa = sqrtl (asa * 2. / (double) (irow1 * (irow1 - 1)));
output_fp,
"\n\n"
" MAXIMUM RELATIVE ASYMMETRY OF THE DRIVING POINT "
"ADMITTANCE\n"
" MATRIX IS %10.3E FOR SEGMENTS %d AND %d\n"
" RMS RELATIVE ASYMMETRY IS %10.3E",
asmx,
nteq,
ntsc,
asa);
} /* if( irow1 >= 2) */
} /* if( netcx.masym != 0) */
/* solution of network equations */
if (netcx.nonet != 0)
{
for (i = 0; i < ndimn; i++)
{
rhnx[i] = CPLX_00;
for (j = 0; j < ndimn; j++)
cmn[j + i * ndimn] = CPLX_00;
}
/* sort network and source data and */
/* assign equation numbers to segments */
nteq = 0;
ntsc = 0;
for (j = 0; j < netcx.nonet; j++)
{
nseg1 = netcx.iseg1[j];
nseg2 = netcx.iseg2[j];
if (netcx.ntyp[j] <= 1)
{
y11r = netcx.x11r[j];
y11i = netcx.x11i[j];
y12r = netcx.x12r[j];
y12i = netcx.x12i[j];
y22r = netcx.x22r[j];
y22i = netcx.x22i[j];
}
else
{
y22r = TP * netcx.x11i[j] / data.wlam;
y12r = 0.;
y12i = 1. / (netcx.x11r[j] * sinl (y22r));
y11r = netcx.x12r[j];
y11i = -y12i * cosl (y22r);
y22r = netcx.x22r[j];
y22i = y11i + netcx.x22i[j];
y11i = y11i + netcx.x12i[j];
if (netcx.ntyp[j] != 2)
{
y12r = -y12r;
y12i = -y12i;
}
} /* if( netcx.ntyp[j] <= 1) */
jump1 = FALSE;
if (vsorc.nsant != 0)
{
for (i = 0; i < vsorc.nsant; i++)
if (nseg1 == vsorc.isant[i])
{
isc1 = i;
jump1 = TRUE;
break;
}
} /* if( vsorc.nsant != 0) */
jump2 = FALSE;
if (!jump1)
{
isc1 = -1;
if (nteq != 0)
{
for (i = 0; i < nteq; i++)
if (nseg1 == nteqa[i])
{
irow1 = i;
jump2 = TRUE;
break;
}
} /* if( nteq != 0) */
if (!jump2)
{
irow1 = nteq;
nteqa[nteq] = nseg1;
nteq++;
}
} /* if( ! jump1 ) */
else
{
if (ntsc != 0)
{
for (i = 0; i < ntsc; i++)
{
if (nseg1 == ntsca[i])
{
irow1 = ndimn - (i + 1);
jump2 = TRUE;
break;
}
}
} /* if( ntsc != 0) */
if (!jump2)
{
irow1 = ndimn - (ntsc + 1);
ntsca[ntsc] = nseg1;
vsrc[ntsc] = vsorc.vsant[isc1];
ntsc++;
}
} /* if( ! jump1 ) */
jump1 = FALSE;
if (vsorc.nsant != 0)
{
for (i = 0; i < vsorc.nsant; i++)
{
if (nseg2 == vsorc.isant[i])
{
isc2 = i;
jump1 = TRUE;
break;
}
}
} /* if( vsorc.nsant != 0) */
jump2 = FALSE;
if (!jump1)
{
isc2 = -1;
if (nteq != 0)
{
for (i = 0; i < nteq; i++)
if (nseg2 == nteqa[i])
{
irow2 = i;
jump2 = TRUE;
break;
}
} /* if( nteq != 0) */
if (!jump2)
{
irow2 = nteq;
nteqa[nteq] = nseg2;
nteq++;
}
} /* if( ! jump1 ) */
else
{
if (ntsc != 0)
{
for (i = 0; i < ntsc; i++)
if (nseg2 == ntsca[i])
{
irow2 = ndimn - (i + 1);
jump2 = TRUE;
break;
}
} /* if( ntsc != 0) */
if (!jump2)
{
irow2 = ndimn - (ntsc + 1);
ntsca[ntsc] = nseg2;
vsrc[ntsc] = vsorc.vsant[isc2];
ntsc++;
}
} /* if( ! jump1 ) */
/* fill network equation matrix and right hand side vector with
*/
/* network short-circuit admittance matrix coefficients. */
if (isc1 == -1)
{
cmn[irow1 + irow1 * ndimn] -=
cmplx (y11r, y11i) * data.si[nseg1 - 1];
cmn[irow1 + irow2 * ndimn] -=
cmplx (y12r, y12i) * data.si[nseg1 - 1];
}
else
{
rhnx[irow1] +=
cmplx (y11r, y11i) * vsorc.vsant[isc1] / data.wlam;
rhnx[irow2] +=
cmplx (y12r, y12i) * vsorc.vsant[isc1] / data.wlam;
}
if (isc2 == -1)
{
cmn[irow2 + irow2 * ndimn] -=
cmplx (y22r, y22i) * data.si[nseg2 - 1];
cmn[irow2 + irow1 * ndimn] -=
cmplx (y12r, y12i) * data.si[nseg2 - 1];
}
else
{
rhnx[irow1] +=
cmplx (y12r, y12i) * vsorc.vsant[isc2] / data.wlam;
rhnx[irow2] +=
cmplx (y22r, y22i) * vsorc.vsant[isc2] / data.wlam;
}
} /* for( j = 0; j < netcx.nonet; j++ ) */
/* add interaction matrix admittance */
/* elements to network equation matrix */
for (i = 0; i < nteq; i++)
{
for (j = 0; j < neqt; j++)
rhs[j] = CPLX_00;
irow1 = nteqa[i] - 1;
rhs[irow1] = CPLX_10;
solves (
cm,
ip,
rhs,
netcx.neq,
1,
data.np,
data.n,
data.mp,
data.m);
cabc (rhs);
for (j = 0; j < nteq; j++)
{
irow1 = nteqa[j] - 1;
cmn[i + j * ndimn] += rhs[irow1];
}
} /* for( i = 0; i < nteq; i++ ) */
/* factor network equation matrix */
factr (nteq, cmn, ipnt, ndimn);
} /* if( netcx.nonet != 0) */
} /* if( netcx.ntsol != 0) */
if (netcx.nonet != 0)
{
/* add to network equation right hand side */
/* the terms due to element interactions */
for (i = 0; i < neqt; i++)
rhs[i] = einc[i];
solves (cm, ip, rhs, netcx.neq, 1, data.np, data.n, data.mp, data.m);
cabc (rhs);
for (i = 0; i < nteq; i++)
{
irow1 = nteqa[i] - 1;
rhnt[i] = rhnx[i] + rhs[irow1];
}
/* solve network equations */
solve (nteq, cmn, ipnt, rhnt, ndimn);
/* add fields due to network voltages to electric fields */
/* applied to structure and solve for induced current */
for (i = 0; i < nteq; i++)
{
irow1 = nteqa[i] - 1;
einc[irow1] -= rhnt[i];
}
solves (cm, ip, einc, netcx.neq, 1, data.np, data.n, data.mp, data.m);
cabc (einc);
if (netcx.nprint == 0)
{
output_fp,
"\n\n\n"
" "
"--------- STRUCTURE EXCITATION DATA AT NETWORK CONNECTION POINTS "
"--------");
output_fp,
"\n"
" TAG SEG. VOLTAGE (VOLTS) CURRENT (AMPS) "
" IMPEDANCE (OHMS) ADMITTANCE (MHOS) POWER\n"
" NO. NO. REAL IMAGINARY REAL IMAGINARY "
" REAL IMAGINARY REAL IMAGINARY (WATTS)");
}
for (i = 0; i < nteq; i++)
{
irow1 = nteqa[i] - 1;
vlt = rhnt[i] * data.si[irow1] * data.wlam;
cux = einc[irow1] * data.wlam;
ymit = cux / vlt;
netcx.zped = vlt / cux;
irow2 = data.itag[irow1];
pwr = .5 * creall (vlt * conjl (cux));
netcx.pnls = netcx.pnls - pwr;
if (netcx.nprint == 0)
output_fp,
"\n"
" %4d %5d %11.4E %11.4E %11.4E %11.4E"
" %11.4E %11.4E %11.4E %11.4E %11.4E",
irow2,
irow1 + 1,
print_creall (vlt),
print_cimagl (vlt),
print_creall (cux),
print_cimagl (cux),
print_creall (netcx.zped),
print_cimagl (netcx.zped),
print_creall (ymit),
print_cimagl (ymit),
pwr);
}
if (ntsc != 0)
{
for (i = 0; i < ntsc; i++)
{
irow1 = ntsca[i] - 1;
vlt = vsrc[i];
cux = einc[irow1] * data.wlam;
ymit = cux / vlt;
netcx.zped = vlt / cux;
irow2 = data.itag[irow1];
pwr = .5 * creall (vlt * conjl (cux));
netcx.pnls = netcx.pnls - pwr;
if (netcx.nprint == 0)
output_fp,
"\n"
" %4d %5d %11.4E %11.4E %11.4E %11.4E"
" %11.4E %11.4E %11.4E %11.4E %11.4E",
irow2,
irow1 + 1,
print_creall (vlt),
print_cimagl (vlt),
print_creall (cux),
print_cimagl (cux),
print_creall (netcx.zped),
print_cimagl (netcx.zped),
print_creall (ymit),
print_cimagl (ymit),
pwr);
} /* for( i = 0; i < ntsc; i++ ) */
} /* if( ntsc != 0) */
} /* if( netcx.nonet != 0) */
else
{
/* solve for currents when no networks are present */
solves (cm, ip, einc, netcx.neq, 1, data.np, data.n, data.mp, data.m);
cabc (einc);
ntsc = 0;
}
if ((vsorc.nsant + vsorc.nvqd) == 0)
return;
output_fp,
"\n\n\n"
" "
"- - - ANTENNA INPUT PARAMETERS - - -\n");
output_fp,
"\n"
" TAG SEG. VOLTAGE (VOLTS) CURRENT (AMPS) "
" IMPEDANCE (OHMS) ADMITTANCE (MHOS) POWER\n"
" NO. NO. REAL IMAG. REAL IMAG. "
" REAL IMAG. REAL IMAG. (WATTS)");
if (vsorc.nsant != 0)
{
for (i = 0; i < vsorc.nsant; i++)
{
isc1 = vsorc.isant[i] - 1;
vlt = vsorc.vsant[i];
if (ntsc == 0)
{
cux = einc[isc1] * data.wlam;
irow1 = 0;
}
else
{
for (j = 0; j < ntsc; j++)
if (ntsca[j] == isc1 + 1)
break;
irow1 = ndimn - (j + 1);
cux = rhnx[irow1];
for (j = 0; j < nteq; j++)
cux -= cmn[j + irow1 * ndimn] * rhnt[j];
cux = (einc[isc1] + cux) * data.wlam;
irow1++;
} /* if( ntsc == 0) */
ymit = cux / vlt;
netcx.zped = vlt / cux;
pwr = .5 * creall (vlt * conjl (cux));
netcx.pin = netcx.pin + pwr;
if (irow1 != 0)
netcx.pnls = netcx.pnls + pwr;
irow2 = data.itag[isc1];
output_fp,
"\n"
"%6d%6d%12.5E%12.5E%12.5E%12.5E"
"%12.5E%12.5E%12.5E%12.5E%12.5E",
irow2,
isc1 + 1,
print_creall (vlt),
print_cimagl (vlt),
print_creall (cux),
print_cimagl (cux),
print_creall (netcx.zped),
print_cimagl (netcx.zped),
print_creall (ymit),
print_cimagl (ymit),
pwr);
} /* for( i = 0; i < vsorc.nsant; i++ ) */
} /* if( vsorc.nsant != 0) */
if (vsorc.nvqd == 0)
return;
for (i = 0; i < vsorc.nvqd; i++)
{
isc1 = vsorc.ivqd[i] - 1;
vlt = vsorc.vqd[i];
cux = cmplx (crnt.air[isc1], crnt.aii[isc1]);
ymit = cmplx (crnt.bir[isc1], crnt.bii[isc1]);
netcx.zped = cmplx (crnt.cir[isc1], crnt.cii[isc1]);
pwr = data.si[isc1] * TP * .5;
cux = (cux - ymit * sinl (pwr) + netcx.zped * cosl (pwr)) * data.wlam;
ymit = cux / vlt;
netcx.zped = vlt / cux;
pwr = .5 * creall (vlt * conjl (cux));
netcx.pin = netcx.pin + pwr;
irow2 = data.itag[isc1];
output_fp,
"\n"
"%6d%6d%12.5E%12.5E%12.5E%12.5E"
"%12.5E%12.5E%12.5E%12.5E%12.5E",
irow2,
isc1 + 1,
print_creall (vlt),
print_cimagl (vlt),
print_creall (cux),
print_cimagl (cux),
print_creall (netcx.zped),
print_cimagl (netcx.zped),
print_creall (ymit),
print_cimagl (ymit),
pwr);
} /* for( i = 0; i < vsorc.nvqd; i++ ) */
/* Free network buffers */
free_ptr ((void *) &ipnt);
free_ptr ((void *) &nteqa);
free_ptr ((void *) &ntsca);
free_ptr ((void *) &vsrc);
free_ptr ((void *) &rhs);
free_ptr ((void *) &cmn);
free_ptr ((void *) &rhnt);
free_ptr ((void *) &rhnx);
return;
}
/*-----------------------------------------------------------------------*/