Average gain computation is wrong #72
Description
The following example:
#!/usr/bin/python3
import PyNEC
nec = PyNEC.nec_context()
geo = nec.get_geometry()
geo.wire(1, 9, 0, 0, 2, 0, 0, 7, .003, 1, 1)
nec.geometry_complete(1)
nec.set_extended_thin_wire_kernel(True)
nec.fr_card(0, 1, 30, 0)
nec.ex_card(0, 1, 5, 0, 1, 0, 0, 0, 0, 0)
nec.gn_card(1, 0, 0, 0, 0, 0, 0, 0)
#nec.rp_card(0, 10, 2, 0, 0, 0, 1301, 0, 0, 10, 90, 0, 0)
nec.rp_card(0, 91, 361, 0, 0, 0, 1301, 0, 0, 1, 1, 0, 0)
p = nec.get_radiation_pattern(0)
print ("avg gain: %.5f solid angle: %.5f"
% (p.get_average_power_gain(), p.get_average_power_solid_angle())
)
Which is a vertical dipole from the nec user documentation e.g. from https://www.nec2.org/other/nec2prt3.pdf example 3 on page 92 but change of wire radius from 30cm (!) to 3mm. I'm also using 1° stepped phi and theta angles instead of 10° stepping for theta and just 0 and 90° for phi (see commented line in code), just 2 values for phi is possible because the vertical dipole has a completely symmetric pattern.
Since this uses ideal ground it should produce an average gain fairly close to 2.0
In my example (which uses python bindings but I'm quite sure the problem is in the C++ code) this prints:
avg gain: 2.05305 solid angle: 2.00000
which is off by more than 2%.
The nec2c code prints for this example with the following input:
CE
GW 1 9 0. 0. 2. 0. 0. 7. .003
GE 1
EK
FR 0 1 0 0 30.
EX 0 1 5 0 1.
GN 1
RP 0 91 361 1301 0. 0. 1. 1.
EN
AVERAGE POWER GAIN: 1.9901E+00 - SOLID ANGLE USED IN AVERAGING: (+2.0000)*PI STERADIANS
when I replace the RP card with:
RP 0 10 2 1301 0. 0. 10. 90.
(which is equivalent to the commented line in the code) I get
AVERAGE POWER GAIN: 1.9876E+00 - SOLID ANGLE USED IN AVERAGING: (+0.5000)*PI STERADIANS
which is still close enough to 2.0
If I enable the commented line in the code (and comment the current rp_card line instead) I'm getting:
avg gain: 3.88009 solid angle: 0.50000
which is way off.
My observation is that the average gain is overestimated systematically. I note this in another project I'm working on.
Hope you have some time looking into this
Ralf