KJ7LNW / xnec2c

Xnec2c is a high-performance multi-threaded electromagnetic simulation package to model antenna near- and far-field radiation patterns for Linux and UNIX operating systems.
https://www.xnec2c.org/
GNU General Public License v3.0
77 stars 16 forks source link

Trouble simulating a biquad for WiFi #35

Closed rklasen closed 10 months ago

rklasen commented 11 months ago

Hi, thank you for this fascinating software! I'm new to nec in general, so please excuse my inexperience.

I'm trying to simulate a biquad antenna for wifi, so 2.4 - 2.5 GHz. With a lot of trial and error, I produced this nec file:

CM --- NEC2 Input File created or edited by xnec2c 4.4.12 ---
CE --- End Comments ---
GW     1     1  -1.50000E-03  0.00000E+00  1.90000E-02  1.50000E-03  0.00000E+00  1.90000E-02  2.50000E-04
GW     2    16  -1.50000E-03  0.00000E+00  1.90000E-02 -1.50000E-03  1.50000E-03  1.90000E-02  2.50000E-04
GW     3    16  -1.50000E-03  1.50000E-03  1.90000E-02 -2.22739E-02  2.22739E-02  1.90000E-02  2.50000E-04
GW     4    16  -2.22739E-02  2.22739E-02  1.90000E-02  0.00000E+00  4.45477E-02  1.90000E-02  2.50000E-04
GW     5    16   0.00000E+00  4.45477E-02  1.90000E-02  2.22739E-02  2.22739E-02  1.90000E-02  2.50000E-04
GW     6    16   2.22739E-02  2.22739E-02  1.90000E-02  1.50000E-03  1.50000E-03  1.90000E-02  2.50000E-04
GW     7    16   1.50000E-03  1.50000E-03  1.90000E-02  1.50000E-03  0.00000E+00  1.90000E-02  2.50000E-04
GW     8    16   1.50000E-03  0.00000E+00  1.90000E-02  1.50000E-03 -1.50000E-03  1.90000E-02  2.50000E-04
GW     9    16   1.50000E-03 -1.50000E-03  1.90000E-02  2.22739E-02 -2.22739E-02  1.90000E-02  2.50000E-04
GW    10    16   2.22739E-02 -2.22739E-02  1.90000E-02  0.00000E+00 -4.45477E-02  1.90000E-02  2.50000E-04
GW    11    16   0.00000E+00 -4.45477E-02  1.90000E-02 -2.22739E-02 -2.22739E-02  1.90000E-02  2.50000E-04
GW    12    16  -2.22739E-02 -2.22739E-02  1.90000E-02 -1.50000E-03 -1.50000E-03  1.90000E-02  2.50000E-04
GW    13    16  -1.50000E-03 -1.50000E-03  1.90000E-02 -1.50000E-03  0.00000E+00  1.90000E-02  2.50000E-04
GW    14    16  -1.50000E-03  0.00000E+00  1.90000E-02 -1.50000E-03  1.50000E-03  1.90000E-02  2.50000E-04
GE     0     0   0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00
EX     0     0     1      0  1.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00
FR     0    11     0      0  2.40000E+03  1.00000E+01  2.50000E+03  0.00000E+00  0.00000E+00  0.00000E+00
NH     0     0     0      0  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00
NE     0    10     1     10 -1.35000E+00  0.00000E+00 -1.35000E+00  3.00000E-01  0.00000E+00  3.00000E-01
RP     0    19    37   1000  0.00000E+00  0.00000E+00  1.00000E+01  1.00000E+01  0.00000E+00  0.00000E+00
EN     0     0     0      0  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00

Which results in this geometry:

image

However, when I run the simulation, the gain graph is a flat line:

image

Right clicking on the gain graph to see a specific frequency causes the program to freeze.


I'm compiling from source on Arch Linux, all the example antennas that I tested work just fine.

So I actually have multiple questions:

Thanks in advance!

rklasen commented 11 months ago

Oh, so I've just tested a few other example files, namely:

Which are both in the 2-3 GHz range. For both, the software freezes when I try to right click the VSWR graph.

So there seems to be some issue with high frequencies?

KJ7LNW commented 11 months ago

However, when I run the simulation, the gain graph is a flat line

Uncheck View->Clamp VSWR to 10. Also, if you have the monitor resolution, make the window bigger so you can see the vertical graduations.

can xnec2c even simulate 2,45 GHz? All examples were a couple 100 MHz.

Should be fine.

For the antenna construction, I assume all units are meters, is that correct?

Correct.

The simulation is very fast, to the point I'm skeptical if it is even correct (For my use case). Can I set the calculation precision somewhere?

Calculation uses a 64-bit double. That should be fine for most use cases. nec2++ has a long double implementation if you need it. About limitations, see this doc, too: http://www.antentop.org/w4rnl.001/amod3.html

What could cause the software to freeze?

I'm not sure, its strange that it hangs after succeeding a calculation. I'll see what can be done about the lockup, hopefully over the next few days.

FYI, this is what it looks like on my desk:

image

rklasen commented 10 months ago

Hi, thanks for the detailed answer.

I've increased the Window size and removed the VSWR clamp, and my windows now look exactly like yours:

image

(I was really hoping the antenna would perform better, but that's hardly a software issue :sweat_smile: )

So now, only the freezing remains. Does that happen to you as well?


For the sake of trying, I've recreated the same antenna ten times larger, and decreased the stimulation frequency tenfold. My hope was that this would give comparable results:

image

The radiation pattern looks similar at least and the VSWR graph is also similar, but the gain plots looks slightly different (I assume dielectric constants matter here).

But most importantly: The software doesn't freeze in this case, I can click in the VSWR graph and inspect the radiation pattern for each simulated frequency.


Thanks again for the quick reply!


EDIT: Now that I have your and my plots side-by-side, I see that your gain plot is actually different from mine, even at the same frequency.

KJ7LNW commented 10 months ago

I was really hoping the antenna would perform better, but that's hardly a software issue 😅 [...] For the sake of trying, I've recreated the same antenna ten times larger, and decreased the stimulation frequency tenfold. My hope was that this would give comparable results [...]

I seem to remember that you want at least 10 "segments" per wavelength. So if you are at 13cm (2.4ghz) then you need about 1 segment/cm of GW wire length. For more accuracy you can increase the segment count to ~100/wavelength. You have GW N 16 ... in most places which is fine if each wire is <16cm. (I didn't actually check the lengths, though)

So now, only the freezing remains. Does that happen to you as well?

It does. Since it can be reproduced, maybe we can find a solution.

... recreated the same antenna ten times larger [and] The software doesn't freeze in this case, I can click in the VSWR graph and inspect the radiation pattern for each simulated frequency.

Interesting. I also had wondered if it was a geometry issue, and not a frequency issue. It looks like it spins 100% CPU inside Set_Interaction_Matrix() so possibly there is a bug there. I didn't write that part of the code but I'll look into it and see if we can figure out the cause.

KJ7LNW commented 10 months ago

I bisected and found this:

first bad commit: [1cc146c4e269bc6964017a53c8ea0fd3dba4d2bf] Disable EM calculation if the selected frequency step has already been calculated.

Now I need to figure out what the bug is...

rklasen commented 10 months ago

I had an error in my geometry indeed. I've tested my same nec file with nec2++ like you suggested, which hinted me towards overlapping wires. I've mistakenly created a wire segment twice, which explains the abysmal VSWR graph for a biquad which should have correct dimensions.

I've read (parts of) the NEC-2 Manual Part III, which also states that any segment should be smaller than lamdba/10 but preferably not smaller than lambda/1000 (Page 7, "Extremely short segments, less than about 10^-3 Lambda, should also be avoided since the similarity of the constant and cosine components of the current expansion leads to numerical inaccuracy.")

So anyway, I corrected my geometry but kept the geometry scaling of 10 times and the frequency division of 10 times to get this:

image

I think I made an error earlier, am now confident that this scaling is a suitable workaround for the frequency issue.


Nonetheless, I am also convinced that the issue is somehow related to frequency and not geometry. In the examples, there is file called 13cm_Yagi.nec, for completeness here:

CM --- NEC2 Input File created or edited by xnec2c 3.5 ---
CM Yagi antenna for 2.4 GHz
CE --- End Comments ---
GW     1    23  0.00000E+00  0.00000E+00 -2.62500E-02  0.00000E+00  0.00000E+00  2.62500E-02  1.50000E-03
GW     2    23 -1.30000E-02  0.00000E+00 -2.87500E-02 -1.30000E-02  0.00000E+00  2.87500E-02  1.50000E-03
GW     3    21  1.20000E-02  0.00000E+00 -2.48000E-02  1.20000E-02  0.00000E+00  2.48000E-02  1.50000E-03
GW     4    20  3.18000E-02  0.00000E+00 -2.44000E-02  3.18000E-02  0.00000E+00  2.44000E-02  1.50000E-03
GW     5    20  5.93000E-02  0.00000E+00 -2.39000E-02  5.93000E-02  0.00000E+00  2.39000E-02  1.50000E-03
GW     6    20  8.93000E-02  0.00000E+00 -2.35000E-02  8.93000E-02  0.00000E+00  2.35000E-02  1.50000E-03
GW     7    20  1.24000E-01  0.00000E+00 -2.32500E-02  1.24000E-01  0.00000E+00  2.32500E-02  1.50000E-03
GW     8    20  1.61000E-01  0.00000E+00 -2.31000E-02  1.61000E-01  0.00000E+00  2.31000E-02  1.50000E-03
GW     9    20  2.00000E-01  0.00000E+00 -2.28500E-02  2.00000E-01  0.00000E+00  2.28500E-02  1.50000E-03
GW    10    20  2.40000E-01  0.00000E+00 -2.28500E-02  2.40000E-01  0.00000E+00  2.28500E-02  1.50000E-03
GW    11    20  2.80000E-01  0.00000E+00 -2.28500E-02  2.80000E-01  0.00000E+00  2.28500E-02  1.50000E-03
GM     0     0  0.00000E+00  0.00000E+00  0.00000E+00 -1.35000E-01  0.00000E+00  0.00000E+00  0.00000E+00
GE     0     0  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00
EX     0     1    12      0  1.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00
FR     0    41     0      0  2.00000E+03  2.00000E+01  2.80000E+03  0.00000E+00  0.00000E+00  0.00000E+00
RP     0    19    37      0  0.00000E+00  0.00000E+00  1.00000E+01  1.00000E+01  0.00000E+00  0.00000E+00
EN     0     0     0      0  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00

The file can be opened and frequencies below 2147 MHz can be clicked on in the VSWR graph without issue:

image

However, when you go slightly higher, the green bar in the VSWR graph disappears, and the program freezes:

image


The noteworthy thing here is that the frequency is set to a completely wrong value just prior to the freezing:

image

I have some limited experience with linear algebra, and I seem to recall that numeric matrix inversions significantly reduce the calculation precision. For example, a double has a usable precision of about 15 decimal places, but when I was dealing with 4x4 matrices of doubles, numeric inversion of these matrices resulted in only 10 usable decimal places.

Put plainly, if $B$ is the numerically (as opposed to analytically) calculated inverse of the matrix $A$, then:

$$ A B = \mathbb{1} $$

Is only true for the first 10 digits. If the frequency is internally handled in Hz, we're already at 2.45e9 Hz, so probably close enough for numeric issues to be relevant.

I haven't read enough of xnec2c code to see if this is the case here too, it's just the first thing that came to mind when I read the function name Set_Interaction_Matrix().

rklasen commented 10 months ago

So from what I can tell, the integrator function here is the one that goes into an endless loop.

I've added a simple iteration counter to exit even if the integral doesn't converge.

int iter = 0;
  int maxIter = 1e1;

  while (TRUE && iter < maxIter) {
    iter++;

That alleviates the freezing, but obviously it doesn't fix the underlying fault. I've tried replacing all doubles with long doubles (and exchanged all log to logl, cos to cosl etc) and added a test_long function accordingly, but that doesn't suffice yet.


I'm a little confused to see a hand-rolled integrator here, since the main app allows selection of the math lib (i.e. intel MKL or NEC2). If we could substitute this integration for a different, numerically more robust method, it may solve the issue.

KJ7LNW commented 10 months ago

I'm a little confused to see a hand-rolled integrator here, since the main app allows selection of the math lib (i.e. intel MKL or NEC2). If we could substitute this integration for a different, numerically more robust method, it may solve the issue.

True. This was originally written in FORTRAN in the 80's and ported to C by Neoklis 5B4AZ in the early 2000's. I worked with Neoklis for a while on an automated antenna optimizer and have started maintaining the xnec2c software for bugs and minor features. Neoklis is the one who hand-ported this to C from Fortran so long ago.

The most CPU-expensive function in NEC2 is factoring and solving with Gaussian Elimination. That xnec2c code now (optionally) supports OpenBLAS, ATLAS and MKL using zgetrf and zgetrs which can be several times faster.

If you understand the math well enough to implement this in OpenBLAS, ATLAS, or Intel MKL for other parts of the code then that would be awesome. See src/mathlib.c and the mathfuncs array:

// To add a new mathfunc:
//   * Update the enum in mathlib.h if the calling convention differs
//   * Add function pointer typedefs in mathlib.h
//   * Add a row here
//   * Implement the new func below to use current_mathlib.
//   * Add a prototype for the new func in mathlib.h
static char *mathfuncs[] = {
    [MATHLIB_ZGETRF] = "zgetrf",
    [MATHLIB_ZGETRS] = "zgetrs"
};
KJ7LNW commented 10 months ago

Also good point about the frequency selection being completely wrong. That might be the cause of the issue.

KJ7LNW commented 10 months ago

Thanks for your help troubleshooting this. If you rebuild from master it should be working for you. Please re-open the issue if you have trouble.

Also check out https://github.com/KJ7LNW/xnec2c-optimize/ if you'd like to optimize your antenna geometry with the Simplex optimization algorithm.

(If you are interested in replacing intx() with a BLAS implementation, or any other improvement, then patches are welcome.)

Thanks again.

-Eric

KJ7LNW commented 10 months ago

(Fixed commit description and force-pushed as a129ee1f7f9dae7242fa2de1043500443ad6b883)

rklasen commented 10 months ago

I can confirm that your fix is working. That means the integrator intx() is not at fault, and was instead simply given a negative frequency, which it could not integrate over and hence went into an endless loop.

I can add a convergence check via PR like I outlined in an earlier comment, if we want to handle edge cases like these, since I think all for/while loops should have a certain exit point.

Other than that, the integrator seems stable enough, so I would prefer not to touch it. From what I can tell this integrator is "only" used for display and interpolation, and not the actual NEC calculations, so it wouldn't affect the radiation patterns etc anyway.

So, again: thanks for your swift fix, and for this software in general. I will also check out the antenna optimizer you liked, seems like I could really benefit from that.

KJ7LNW commented 10 months ago

You're welcome, I'm glad its working for you!