<div dir="ltr">Dear Petsc users,<div><br></div><div>We are working on a matrix-free solver for NS equations. We first developed the solver in the SNES framework. But because we want to take control of the Newton iterations and reuse some of our existing code, we stop using SNES and begin developing the matrix-free solver in KSP framework.</div>
<div><br>But we have an issue when using MFFD and VecSetValuesBlocked together.</div><div>My code looks like</div><div><br></div><div>call MatMFFDSetFunction(pet_mat_mf, flowsol_ng_mf, ctx, ierpetsc)<br>call MatMFFDSetBase(pet_mat_mf, pet_solu_snes, PETSC_NULL_OBJECT, ierpetsc)<br>
<br></div><div class="gmail_extra">The vec pet_solu_snes are set to blocksize 5.<br><br></div><div class="gmail_extra">flowsol_ng_mf looks like:<br></div><div class="gmail_extra">subroutine flowsol_ng_mf ( ctxx, pet_solu_snesa, pet_rhs_snesa, ierrpet )<br>
...............<br> call VecGetBlockSize(pet_rhs_snesa,i,ierpetsc)<br> write(*,*) "pet_rhs_snesa",i<br>.............<br> call VecSetValuesBlocked ( pet_rhs_snesa, ndperl, pirow2, rhsloc, ADD_VALUES, ierpetsc )<br>
..........<br></div><div class="gmail_extra">end<br><br></div><div class="gmail_extra">I'm worried about the block size of pet_rhs_snesa is not properly initialized. So I checked the blocksize in the function evaluation.<br>
<br></div><div class="gmail_extra">The output with single CPU looks like:<br><br> 0 KSP preconditioned resid norm 7.818306841541e+00 true resid norm 9.619278462343e-03 ||r(i)||/||b|| 1.000000000000e+00<br> pet_rhs_snesa 1<br>
pet_rhs_snesa 1<br> pet_rhs_snesa 5<br> 1 KSP preconditioned resid norm 1.739193723080e+00 true resid norm 2.564460053330e-03 ||r(i)||/||b|| 2.665958848545e-01<br> pet_rhs_snesa 1<br> pet_rhs_snesa 5<br>
2 KSP preconditioned resid norm 3.816590845074e-01 true resid norm 8.222646506071e-04 ||r(i)||/||b|| 8.548090730777e-02<br> pet_rhs_snesa 1<br> pet_rhs_snesa 5<br><br></div><div class="gmail_extra">I think the block size is not initialized correctly for the first call, which leads to the misuse of VecSetValuesBlocked. <br>
<br></div><div class="gmail_extra">One way to circumvent is to create a temporary vec with correct blocksize, assemble it with VecSetValuesBlocked and then copy back. But Is there other solutions? Thank you very much.<br>
<br><br><br></div><div class="gmail_extra">Song<br></div><div class="gmail_extra"><br><br></div><div class="gmail_extra"><br><div class="gmail_quote">
On Sat, Mar 8, 2014 at 8:51 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">
<div>On Sat, Mar 8, 2014 at 7:45 PM, Song Gao <span dir="ltr"><<a href="mailto:song.gao2@mail.mcgill.ca" target="_blank">song.gao2@mail.mcgill.ca</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"><div dir="ltr">Thank you all. <br><br>Most of the difficulties come from the current structure of the program. Its size and "age" make any radical modification a challenge. The large uses of global variables and deprecated "goto" statements call for a revision, which however is unlikely to occur, against our better judgement...<br>
<br>That being said, the current tools provided byPETSc are sufficient, however not necessarily convenient for our purposes. As a wishful thinking we can say that the implementation of PETSc SNES features into existing codes would be easier if the outer Newton-like loop could be managed by the original code, out of the SNES context. On the other hand, we realize that this might be in contrast with the requirements of PETSc itself. <br>
</div></blockquote><div><br></div></div><div>I guess I wanted to know specifically why you cannot just dump all of the other things in the legacy code into SNESSetUpdate()</div><div>or SNESSetConvergenceTest(). Was there something about the way we structured the hooks that made this impossible?</div>
<div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr">The application of the Eistenstat-Walker method together with a Matrix-Free approach are of key importance. Thus, as kindly suggested by Jed Brown, implementing our own version of EW scheme might turn out to be the only way. We would be happy to provide more details if you think that they might be helpful for the future developments of PETSc.<div class="gmail_extra">
<br><br><div class="gmail_quote">On Fri, Mar 7, 2014 at 2:39 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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">
<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><div>On Fri, Mar 7, 2014 at 1:21 PM, Song Gao <span dir="ltr"><<a href="mailto:song.gao2@mail.mcgill.ca" target="_blank">song.gao2@mail.mcgill.ca</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"><div dir="ltr"><div><div>Hello,<br><br></div>We are working on a legacy codes which solves the NS equations. The codes take control of newton iterations itself and use KSP to solve the linear system. </div>
<div><br></div><div>
We modified the codes, used SNES to control the newton iteration and changed the code to matrix free fashion. But the legacy codes did a lot of other things between two newton iterations (such as output solution, update variables....). I know we could use linesearchpostcheck but it is difficulty to do that correctly. Therefore, we decide to go back to the KSP framework but still use matrix free.<br>
</div></div></blockquote><div><br></div></div><div>What makes it difficult?</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr"><div></div>When using SNES, we always use the runtime option -snes_ksp_ew, we observe that for some test cases, the residual stalls without -snes_ksp_ew, but converges with-snes_ksp_ew. So I'm thinking if it is possible to use -snes_ksp_ew in KSP? Thanks in advance.<br>
<div><br></div><div><br></div></div>
</blockquote></div></div><span><font color="#888888"><br><br clear="all"><span><font color="#888888"><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
</font></span></font></span></div></div>
</blockquote></div><br></div></div>
</blockquote></div></div><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
</div></div></div>
</blockquote></div><br></div></div>