<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#ffffff" text="#000000">
Mark, Jed,<br>
<br>
Thanks a lot for your response.<br>
<br>
<blockquote
cite="mid:648EC7C0-2485-4322-823B-9AB8488063E3@columbia.edu"
type="cite">
<div>It looks like you have a largish \omega\mu\sigma because, as
Jed mentioned, one would not expect GAMG to work well for the
curl-curl part.</div>
</blockquote>
<br>
Well in my case it is the other way around. I solve this equation
for low \omega and possibly very low \sigma. The curl-curl term
dominates always.<br>
<br>
<blockquote
cite="mid:648EC7C0-2485-4322-823B-9AB8488063E3@columbia.edu"
type="cite">
<div><br>
</div>
<div>A few things to try:</div>
<div><br>
</div>
<div>1) '-pc_gamg_type sa' will provide what should be a better
solver for SPD problems.</div>
</blockquote>
<br>
Unfortunately, the convergence with this option is worse for several
operators I've tried. <br>
Would it help to specify coordinates vertices?<br>
<br>
<blockquote
cite="mid:648EC7C0-2485-4322-823B-9AB8488063E3@columbia.edu"
type="cite">
<div><br>
</div>
<div>2) Try:</div>
<div><br>
</div>
<div>
<div>-pc_type hypre ! instead of -pc_type gamg</div>
<div>-pc_hypre_type boomeramg </div>
<div>-pc_hypre_boomeramg_strong_threshold 0.6</div>
</div>
</blockquote>
<br>
When configuring with hypre I get "Cannot use hypre with complex
numbers it is not coded for this capability".<br>
I can reformulate my problem in terms of real matrix:<br>
<br>
If C = A + iB <br>
Cx=b == [A -B; B A] [xr; xi] = [br; bi]<br>
<br>
But I'm not sure that making spectra two times larger would not
influence my condition number?<br>
<br>
<blockquote
cite="mid:648EC7C0-2485-4322-823B-9AB8488063E3@columbia.edu"
type="cite">
<div>
<div><br>
</div>
<div>3) You say 'staggard' but I just see E here. Do you have E
on faces? I forget how staggering works here. If E is cell
centered then you have a system of 3x3 blocks (with the right
ordering) and GAMG might benefit from setting the block size
to tell it this:</div>
<div><br>
</div>
<div>MatSetBlockSize(mat,3); <br>
</div>
</div>
</blockquote>
<br>
I have fields on the edges, but I can formulate another staggering
scheme where fields are on the faces. Would it help?<br>
My matrix is 3x3 block. For instance, if I have 100^3 grid, then
matrix is 3*10^6 and each block is of 10^6 size. <br>
Which block size should I pass to MatSetBlockSize? 10^6 or 3?
Because from its description it is not obvious.<br>
<br>
<blockquote
cite="mid:648EC7C0-2485-4322-823B-9AB8488063E3@columbia.edu"
type="cite">
<div>
<div>And Jed's answer addresses your 2nd question about
null-space. These solvers will degrade as \omega\mu\sigma
gets smaller.</div>
</div>
</blockquote>
<br>
I do observe this for low \omega and low \sigma (e.g. in the air). I
was thinking how could one project out this null-space. <br>
Hitpmair, as Jed pointed, gives some clues. However it requires some
efforts to implement and wanted first to try petsc built in stuff.<br>
<br>
<blockquote
cite="mid:648EC7C0-2485-4322-823B-9AB8488063E3@columbia.edu"
type="cite">
<div>
<div><br>
</div>
<div>Mark</div>
<div><br>
<div>
<div>On Jan 19, 2012, at 5:37 PM, Jed Brown wrote:</div>
<br class="Apple-interchange-newline">
<blockquote type="cite">
<div class="gmail_quote">On Wed, Jan 18, 2012 at 08:22,
Alexander Grayver <span dir="ltr"><<a
moz-do-not-send="true"
href="mailto:agrayver@gfz-potsdam.de">agrayver@gfz-potsdam.de</a>></span>
wrote:<br>
<blockquote class="gmail_quote" style="margin: 0px 0px
0px 0.8ex; border-left: 1px solid rgb(204, 204, 204);
padding-left: 1ex; position: static; z-index: auto;">
Hello petsc team,<br>
<br>
I solve 3D vector Helmholtz equation like following:<br>
<br>
\nabla \times \nabla \times E + i\omega\mu\sigma E =
-J<br>
</blockquote>
<div><br>
</div>
<div>Multigrid methods for curl-curl problems are pretty
specialized. ML and Hypre have support for specific
discretizations, I don't know if they support an
imaginary shift. The PETSc interface to these packages
does not currently support their special Maxwell
interfaces.</div>
<div><br>
</div>
<div>The methods are generally based on edge relaxation
or auxiliary space preconditioning, see Hiptmair or
Arnold, Falk, and Winther for the mathematical
background.</div>
<div> </div>
<blockquote class="gmail_quote" style="margin: 0pt 0pt
0pt 0.8ex; border-left: 1px solid rgb(204, 204, 204);
padding-left: 1ex;">
<br>
I use structured staggered grid and FD. The solution
is a vector that consists of three parts E = {Ex Ey
Ez}. The operator is symmetric matrix with complex
numbers on diagonal.<br>
I'm interested in solving this system with iterative
techniques. I applied newly presented gamg and it
gives promising results, but all I did is just:<br>
-ksp_type tfqmr -pc_type gamg<br>
<br>
I played with different ksp_type and gamg options
which are listed on PCGAMG doc page, but nothing
improved convergence.<br>
Could you please guide me a bit through usage of this
technique?<br>
The precise questions are:<br>
1. Do I have to do something to say petsc that my
equation is a vector equation? Is it important for
gamg pc?<br>
2. Should I take into account null-space using
KSPSetNullSpace? Since it is well known that as \omega
or \sigma get small, null-space of geometric term
(curl curl operator) starts to dominate and system
gets more ill-conditioned.<br>
3. Which options for gamg may improve convergence in
my case?<br>
<br>
<br>
Thanks a lot in advance.<br>
<font color="#888888">
<br>
-- <br>
Regards,<br>
Alexander<br>
<br>
</font></blockquote>
</div>
<br>
</blockquote>
</div>
<br>
</div>
</div>
</blockquote>
<br>
<br>
<pre class="moz-signature" cols="72">--
Regards,
Alexander</pre>
</body>
</html>