On Mon, Feb 20, 2012 at 4:52 PM, Patrick Alken <span dir="ltr"><<a href="mailto:patrick.alken@colorado.edu">patrick.alken@colorado.edu</a>></span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000">
Hello all,<br>
<br>
I am having great difficulty solving a 3D finite difference
equation in spherical coordinates. I am solving the equation in a
spherical shell region S(a,b), with the boundary conditions being
that the function is 0 on both boundaries (r = a and r = b). I
haven't imposed any boundary conditions on theta or phi which may be
a reason its not converging. The phi boundary condition would be
that the function is periodic in phi, but I don't know if this needs
to be put into the matrix somehow?<br></div></blockquote><div><br></div><div>1) The periodicity appears in the definition of the FD derivative in phi. Since this is Cartesian, you can use a DA in 3D, and make one</div>
<div>direction periodic.</div><div><br></div><div>2) Don't you have a coordinate singularity at the pole? This is why every code I know of uses something like a Ying-Yang grid.</div><div><br></div><div> Matt</div><div>
</div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div bgcolor="#FFFFFF" text="#000000">
I nondimensionalized the equation before solving which helped a
little bit. I've also scaled the matrix and RHS vectors by their
maximum element to make all entries <= 1.<br>
<br>
I've tried both direct and iterative solvers. The direct solvers
give a fairly accurate solution for small grids but seem unstable
for larger grids. The PETSc iterative solvers converge for very
small grids but for medium to large grids don't converge at all.<br>
<br>
When running with the command (for a small grid):<br>
<br>
<b>> ./main -ksp_converged_reason -ksp_monitor_true_residual
-pc_type svd -pc_svd_monitor</b><br>
<br>
I get the output:<br>
<br>
SVD: condition number 5.929088512946e+03, 0 of 1440 singular
values are (nearly) zero<br>
SVD: smallest singular values: 2.742809162118e-04
2.807446554985e-04 1.548488288425e-03 1.852332719983e-03
2.782708934678e-03<br>
SVD: largest singular values : 1.590835571953e+00
1.593368145758e+00 1.595771695877e+00 1.623691828398e+00
1.626235829632e+00<br>
0 KSP preconditioned resid norm 2.154365616645e+03 true resid norm
8.365589263063e+00 ||r(i)||/||b|| 1.000000000000e+00<br>
1 KSP preconditioned resid norm 4.832753933427e-10 true resid norm
4.587845792963e-12 ||r(i)||/||b|| 5.484187244549e-13<br>
Linear solve converged due to CONVERGED_RTOL iterations 1<br>
<br>
When plotting the output of this SVD solution, it looks pretty good,
but svd isn't practical for larger grids.<br>
<br>
Using the command (on the same grid):<br>
<br>
<b>> ./main -ksp_converged_reason -ksp_monitor_true_residual
-ksp_compute_eigenvalues -ksp_gmres_restart 1000 -pc_type none</b><br>
<br>
The output is attached. There do not appear to be any 0 eigenvalues.
The solution here is much less accurate than the SVD case since it
didn't converge.<br>
<br>
I've also tried the -ksp_diagonal_scale -ksp_diagonal_scale_fix
options which don't help very much.<br>
<br>
Any advice on how to trouble shoot this would be greatly
appreciated.<br>
<br>
Some things I've checked already:<br>
<br>
1) there aren't any 0 rows in the matrix<br>
2) using direct solvers on very small grids seems to give decent
solutions<br>
3) there don't appear to be any 0 singular values or eigenvalues<br>
<br>
Perhaps the matrix has a null space, but I don't know how I would
find out what the null space is? Is there a tutorial on how to do
this?<br>
<br>
Thanks in advance!<br>
<br>
</div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>