<div class="gmail_quote">On Sun, Aug 14, 2011 at 11:33, Paul Anton Letnes <span dir="ltr">&lt;<a href="mailto:paul.anton.letnes@gmail.com">paul.anton.letnes@gmail.com</a>&gt;</span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div id=":72">Ingve is my supervisor - so it&#39;s a good guess! But, no. BiCGSTAB works well for his code, by the way. He gets a massive speedup compared to direct solvers.<br></div></blockquote><div><br></div><div>Heh, okay. So it works when the problem is of the kind I&#39;ve seen before and the theory I&#39;m familiar with says it should work. Comforting perhaps, but doesn&#39;t solve your problem.</div>
<div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;"><div id=":72">
<br>
I am working on a different equation called the Reduced Rayleigh Equation (RRE) for 2 dimensions, which has not been solved directly before (AFAIK).  Have a look here:</div></blockquote><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div id=":72">
<a href="http://www.tandfonline.com/doi/abs/10.1088/0959-7174/6/3/006" target="_blank">http://www.tandfonline.com/doi/abs/10.1088/0959-7174/6/3/006</a><br>
The form of the solution is given in Eq. 5, and the equation itself is given in Eq. 8 with some definitions in Eq. 9-11. As far as I can tell there is only 2 terms, so it&#39;s a first kind equation. But by all means, feel free to correct me on that one. The cited paper uses a completely different approach (perturbation theory where 0. order is the flat surface) to solving the RRE, by the way. Our approach is more &quot;brute force&quot; - discretize the integral on a 2D grid in \vec{q}_\parallel space and solve the linear equation system.<br>
</div></blockquote><div><br></div><div>Hmm, I would have to understand the spectral properties of this kernel better. Is involves products instead of differences (like the more conventional Green&#39;s function you showed earlier) and is clearly non-compact, so my earlier comment about trying to invert a compact operator does not apply.</div>
<div><br></div><div>But you also likely won&#39;t get any nice decay in the spectrum with which to use unpreconditioned Krylov methods.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div id=":72">
<br>
The idea is that the RRE is less memory intensive than the rigorous approach in the paper you cited. However, if we can&#39;t avoid the LU factorization, the RRE approach will be much more CPU intensive, and of lesser interest than we would have hoped...<br>

<br>
We have submitted a first paper to ArXiv; it will appear on monday, I think.<br>
<div class="im"><br>
&gt; Note that the system has the form<br>
&gt;<br>
&gt; J_H(x_\parallel | \omega) = J_H(x_\parallel | \omega)_{inc} + \int (...) G(x | x&#39;) J_H(x_\parallel&#39; | \omega) + \int (...) G(x | x&#39;) J_E(x_\parallel&#39; | \omega)<br>
&gt;<br>
&gt; which is the form of a second order integral equation. I assume the incident field J_H(...)_{inc} is known in this equation. If you dropped the term on the left hand side in this equation, you would have a Fredholm integral equation of the first kind to &quot;solve&quot;, which is problematic at a mathematical level due to ill-posedness.<br>

<br>
</div>But since it&#39;s of the second kind, it&#39;s better posed, right?<br></div></blockquote><div><br></div><div>Right, because the eigenvalues decay rapidly to 1 instead of decaying rapidly to 0.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div id=":72">
<div class="im"><br>
&gt; I have downloaded and attempted to use a different BiCGSTAB code. It converges, but only after several hundred (about 400) for a very small (not physically interesting) problem. It would appear that if we are to get good performance, some form of preconditioning is necessary.<br>

&gt;<br>
&gt; Do the eigenvalues decay quickly? Can you plot some eigenvalues? They should decay rapidly to a positive value like (with appropriate scaling) 1.<br>
<br>
</div>I&#39;ll have a look at the eigenvalues. Is your suggestion to plot the magnitude of the complex eigenvectors as a function of their index, after sorting them? I guess that must be what you mean by &quot;plotting&quot;.</div>
</blockquote></div><br><div>Let&#39;s start with a scatter plot of the eigenvalues. Can you do a problem that is representative of the physics in less than, say, 1000 degrees of freedom? If so, I would just use Matlab (or Octave, etc).  You want to be able to plot the eigenvector associated with a chosen eigenvalue in some way that is meaningful to you. We want to see if the wavelength (in terms of the variables you are discretizing over) of the modes has some useful correlation with the size of the associated eigenvalues. If so, we may be able to build some sort of multigrid preconditioner.</div>