June 2016 - petsc-dev - mailman3.cels.anl.gov

Re: [petsc-dev] Problem with assembly of Jacobian in ex12 & ex77
by Sander Arens 28 Jun '16

28 Jun '16

What are you trying to do with it? I don't think there's more documentation on how to use PetscFE than what's already in the online documentation and in the examples. What do you mean with "grasp the concept"? Thanks, Sander On 27 June 2016 at 17:41, <rickcha(a)googlemail.com> wrote: > That did do the trick with the examples, thank you. I am still struggling > with my own simple example however. Is there any documentation available on > how to use petscFE? I could not find anything besides the three examples. > And I yet have to grasp the concept. > > > Thanks in advance, > M. Hartig > > On 27 Jun 2016, at 15:49, Sander Arens <Sander.Arens(a)ugent.be> wrote: > > I think you forgot to set the order of the finite element spaces, which > default to 0. Does it work if you try for example ex77 with the following > options found here? > <https://bitbucket.org/petsc/petsc/src/8bc5349a152c81f25b59a9592e1e5936c0186…> > > Thanks, > Sander > > On 27 June 2016 at 15:28, <rickcha(a)googlemail.com> wrote: > >> Dear Petsc developer team, >> >> while learning the ways of the petscFE I came across a curiosity I’d like >> to point out to you. The Jacobians that are created from the petscFE object >> seem to have only 0-entries (apart from the obviously empty regions) in >> examples 12 and 77. That might be due to a limited understanding from my >> side on how exactly petscFE and the mentioned examples work, but it struck >> me as unusual. >> >> Thanks in advance, >> >> M. Hartig > > > >

1 0

Re: [petsc-dev] TSComputeRHSJacobianConstant issue
by Barry Smith 28 Jun '16

28 Jun '16

Fixed in barry/fix-tscomputerhsjacobianconstant/maint merged to next for nightly testing Had I any brains at all I would have had Dominic use git bisect to determine the cause (which was Barry who added an optimization that broke Jed's cool trick for constant Jacobians) and this would have been resolved long ago. Curse on Jed for not adding his test cases, that he did write at least, to the automatic testing list in the makefile. Barry > On Mar 23, 2016, at 5:54 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > > This problem still exists. What do we do? Broken code that a beginner could easily start with is not a good thing. > > Barry > >> On Feb 27, 2016, at 5:21 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: >> >> >> Jed and Emil, >> >> Yes, this is a terrible mess. The TSComputeRHSJacobianConstant() will only work for methods that do not change entries in the Jacobian matrix. And for some reason the runex3_4 and _5 which do not work are turned off in the makefile ! It looks like Rosw does work because it does not change the matrix? And yet it is still turned off? >> >> I don't know what a good fix is. If we want to keep this capability then the TSComputeRHSJacobianConstant() interface has to keep a copy of the original user provided matrix and copy it into the matrices used by the method each time. This may lead to ugly special case code. >> >> Something has to be done about (simplest thing is just tossing the entire concept) because now there is a very broken API in TS that beginning users are likely to try and get upset with. >> >> Barry >> >> >>> On Feb 26, 2016, at 11:04 AM, Dominic Meiser <dmeiser(a)txcorp.com> wrote: >>> >>> Hi all, >>> >>> I'm running into problems with TS with a time independent, linear >>> rhs. Following the user's manual, I use >>> >>> ierr=TSSetProblemType(ts, TS_LINEAR);CHKERRQ(ierr); >>> ierr=TSSetRHSFunction(ts, 0, TSComputeRHSFunctionLinear, 0);CHKERRQ(ierr); >>> ierr=TSSetRHSJacobian(ts, mat, mat, TSComputeRHSJacobianConstant, 0);CHKERRQ(ierr); >>> >>> But running the code I get: >>> >>> [dmeiser@ivyamd src]$ ./sd -ts_monitor -snes_monitor -ksp_monitor_true_residual -snes_converged_reason -ts_type beuler -ts_max_steps 20 TS dt 0.01 time 0. >>> 0 SNES Function norm 1.433485640883e+00 >>> 0 KSP preconditioned resid norm 1.433485640883e-02 true resid norm 1.433485640883e+00 ||r(i)||/||b|| 1.000000000000e+00 >>> 1 KSP preconditioned resid norm 6.528150237447e-18 true resid norm 6.122666587510e-16 ||r(i)||/||b|| 4.271173992185e-16 >>> 1 SNES Function norm 1.419007435910e+02 >>> Nonlinear solve converged due to CONVERGED_ITS iterations 1 >>> 1 TS dt 0.01 time 0.01 >>> 0 SNES Function norm 1.433342292319e+02 >>> Nonlinear solve did not converge due to DIVERGED_LINEAR_SOLVE iterations 0 >>> [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- >>> [0]PETSC ERROR: >>> [0]PETSC ERROR: TSStep has failed due to DIVERGED_NONLINEAR_SOLVE, increase -ts_max_snes_failures or make negative to attempt recovery >>> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting. >>> [0]PETSC ERROR: Petsc Development GIT revision: pre-tsfc-330-gda8116b GIT Date: 2016-02-23 15:17:17 -0600 >>> [0]PETSC ERROR: ./sd on a cplx named ivyamd.txcorp.com by dmeiser Fri Feb 26 11:46:29 2016 >>> [0]PETSC ERROR: Configure options --with-scalar-type=complex --with-x=0 --with-debugging=1 --with-ssl=0 --with-fortran-kernels=0 --with-pthread=0 --download-fftw=1 --with-mpi=1 --download-mpich=1 --with-mpiuni-fortran-binding=0 --with-fortran-interfaces >>> [0]PETSC ERROR: #1 TSStep() line 3346 in /scr_ivyamd/dmeiser/SD/petsc/src/ts/interface/ts.c >>> [0]PETSC ERROR: #2 TSSolve() line 3544 in /scr_ivyamd/dmeiser/SD/petsc/src/ts/interface/ts.c >>> [0]PETSC ERROR: #3 main() line 63 in /scr_ivyamd/dmeiser/SD/src/sd.c >>> [0]PETSC ERROR: PETSc Option Table entries: >>> [0]PETSC ERROR: -ksp_monitor_true_residual >>> [0]PETSC ERROR: -snes_converged_reason >>> [0]PETSC ERROR: -snes_monitor >>> [0]PETSC ERROR: -ts_max_steps 2 >>> [0]PETSC ERROR: -ts_monitor >>> [0]PETSC ERROR: -ts_type beuler >>> [0]PETSC ERROR: ----------------End of Error Message -------send entire error message to petsc-maint(a)mcs.anl.gov---------- >>> application called MPI_Abort(MPI_COMM_WORLD, 91) - process 0 >>> [unset]: aborting job: >>> application called MPI_Abort(MPI_COMM_WORLD, 91) - process 0 >>> >>> >>> It seems that the linear system is being solved just fine. But >>> apparently it's the wrong linear problem because the SNES >>> function norm doesn't decrease. >>> >>> I think the problem might be related to the following behaviour >>> of `scr/ts/examples/tutorials/ex3.c`: >>> >>> >>> [dmeiser@ivyamd tutorials]$ ./ex3 -nox -ts_type beuler -ts_max_steps 2 -ksp_converged_reason -snes_monitor >>> Solving a linear TS problem on 1 processor >>> Timestep 0: step size = 0.000143637, time = 0., 2-norm error = 1.01507e-15, max norm error = 3.10862e-15 >>> 0 SNES Function norm 2.018485956680e+03 >>> Linear solve converged due to CONVERGED_RTOL iterations 1 >>> 1 SNES Function norm 1.207665524289e+05 >>> Timestep 1: step size = 0.000143637, time = 0.000143637, 2-norm error = 0.000599907, max norm error = 0.000863771 >>> 0 SNES Function norm 1.195422385486e+05 >>> Linear solve did not converge due to DIVERGED_PCSETUP_FAILED iterations 0 >>> PCSETUP_FAILED due to FACTOR_NUMERIC_ZEROPIVOT >>> [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- >>> [0]PETSC ERROR: >>> [0]PETSC ERROR: TSStep has failed due to DIVERGED_NONLINEAR_SOLVE, increase -ts_max_snes_failures or make negative to attempt recovery >>> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting. >>> [0]PETSC ERROR: Petsc Development GIT revision: pre-tsfc-330-gda8116b GIT Date: 2016-02-23 15:17:17 -0600 >>> [0]PETSC ERROR: ./ex3 on a cplx named ivyamd.txcorp.com by dmeiser Fri Feb 26 11:36:28 2016 >>> [0]PETSC ERROR: Configure options --with-scalar-type=complex --with-x=0 --with-debugging=1 --with-ssl=0 --with-fortran-kernels=0 --with-pthread=0 --download-fftw=1 --with-mpi=1 --download-mpich=1 --with-mpiuni-fortran-binding=0 --with-fortran-interfaces >>> [0]PETSC ERROR: #1 TSStep() line 3346 in /scr_ivyamd/dmeiser/SD/petsc/src/ts/interface/ts.c >>> [0]PETSC ERROR: #2 TSSolve() line 3544 in /scr_ivyamd/dmeiser/SD/petsc/src/ts/interface/ts.c >>> [0]PETSC ERROR: #3 main() line 216 in /scr_ivyamd/dmeiser/SD/petsc/src/ts/examples/tutorials/ex3.c >>> [0]PETSC ERROR: PETSc Option Table entries: >>> [0]PETSC ERROR: -ksp_converged_reason >>> [0]PETSC ERROR: -nox >>> [0]PETSC ERROR: -snes_monitor >>> [0]PETSC ERROR: -ts_max_steps 2 >>> [0]PETSC ERROR: -ts_type beuler >>> [0]PETSC ERROR: ----------------End of Error Message -------send entire error message to petsc-maint(a)mcs.anl.gov---------- >>> application called MPI_Abort(MPI_COMM_WORLD, 91) - process 0 >>> [unset]: aborting job: >>> application called MPI_Abort(MPI_COMM_WORLD, 91) - process 0 >>> >>> >>> I'm aware of Jed's commit e1244c6 which I thought was supposed to >>> enable constant linear rhs without IFunction and IJacobian. Note >>> that e1244c6 also added tests that apparently produced very >>> different output (see e.g. >>> src/ts/examples/tutorials/output/ex3_4.out). >>> >>> Is constant rhs without IFunction and IJacobian no longer >>> supported or did it just breat at some point? Do I need to call >>> TSRHSJacobianSetReuse in addition to the lines of code above? Is >>> it better to treat such systems using TSComputeIFunctionLinear >>> and TSComputeIJacobianConstant? >>> >>> Thanks, >>> Dominic >>> >>> -- >>> Dominic Meiser >>> Tech-X Corporation - 5621 Arapahoe Avenue - Boulder, CO 80303 >> >

1 0

Re: [petsc-dev] MatMFFDSetType
by Barry Smith 27 Jun '16

27 Jun '16

I have added MatFDColoringSetType() in the branch barry/add-missing-apis/maint and merged the branch into next for testing. If it builds and runs everywhere I will put it into the maint branch (and master) and it will be in the next patch release. Barry > On Jun 26, 2016, at 8:49 PM, Adrian Croucher <a.croucher(a)auckland.ac.nz> wrote: > > > > On 27/06/16 13:44, Barry Smith wrote: >> Yes, you are correct, there is no MatFDSetType() and in the SetFromOptions we totally cheat http://www.mcs.anl.gov/petsc/petsc-current/src/mat/matfd/fdmatrix.c.html#Ma… and just stick the value directly into the data structure. For now you can call >> PetscOptionsSetValue(NULL,"-mat_fd_type","ds") in your code (before you start using the MatFD stuff) and we'll add the routine to the master branch. > > Thanks, that workaround does the trick for now. >> I didn't think anyone would care about the difference between ds and wp with MatFD. > It makes a noticeable difference for the problem I'm solving right now. That may be because the problem is not non-dimensionalized (as yet anyway) and so some of the solution elements are of different magnitudes than others. If I non-dimensionalize them it may possibly not matter any more about wp or ds. > > - Adrian > -- > Dr Adrian Croucher > Senior Research Fellow > Department of Engineering Science > University of Auckland, New Zealand > email: > a.croucher(a)auckland.ac.nz > > tel: +64 (0)9 923 84611 >

1 0

Re: [petsc-dev] Problem with assembly of Jacobian in ex12 & ex77
by Sander Arens 27 Jun '16

27 Jun '16

I think you forgot to set the order of the finite element spaces, which default to 0. Does it work if you try for example ex77 with the following options found here? <https://bitbucket.org/petsc/petsc/src/8bc5349a152c81f25b59a9592e1e5936c0186…> Thanks, Sander On 27 June 2016 at 15:28, <rickcha(a)googlemail.com> wrote: > Dear Petsc developer team, > > while learning the ways of the petscFE I came across a curiosity I’d like > to point out to you. The Jacobians that are created from the petscFE object > seem to have only 0-entries (apart from the obviously empty regions) in > examples 12 and 77. That might be due to a limited understanding from my > side on how exactly petscFE and the mentioned examples work, but it struck > me as unusual. > > Thanks in advance, > > M. Hartig

1 0

Problem with assembly of Jacobian in ex12 & ex77
by rickcha＠googlemail.com 27 Jun '16

27 Jun '16

Dear Petsc developer team, while learning the ways of the petscFE I came across a curiosity I’d like to point out to you. The Jacobians that are created from the petscFE object seem to have only 0-entries (apart from the obviously empty regions) in examples 12 and 77. That might be due to a limited understanding from my side on how exactly petscFE and the mentioned examples work, but it struck me as unusual. Thanks in advance, M. Hartig

1 0

Re: [petsc-dev] asm / gasm
by Mark Adams 27 Jun '16

27 Jun '16

Barry, I am having problems with your branch on ksp/ex56. First I was getting an error with the increase of the overlap that you added (I don't want that to be the default). So I set this back to 0 in gamg.c. Next, I get the same error with the vec out of range because ex56 has BCs in the matrix. I created a branch mark/fix-gamg-asm-aggs from your branch and added the dummy ASM block. But I am getting all kinds of errors. I get errors in get MatSubMatrix and a funny error when it tries to free a prefix string. I've attached the valgrind output of one process, of 8. There is a lot of noise with string operations but it looks like the first non-string thing is in MatGetSubMatrices. ex56 is elasticity. Maybe there is a block size problem?

1 0

Re: [petsc-dev] asm / gasm
by Mark Adams 27 Jun '16

27 Jun '16

Garth: we have something for you to try: branch barry/fix-gamg-asm-aggs add: -pc_gamg_use_agg_asm -mg_levels_sub_pc_type lu > > > > I added some code to add a block on each processor for any singletons, > because the MIS code strips these (so, yes, not a true MIS). I should do > this for users that put a live equation in a singleton, like a non-homo > Diri BC. I can add that you your branch. Let me know if that is OK. > > I do not understand this. Do you mean a variable that is only coupled > to itself? So the row and column for the variable have only an entry on the > diagonal? Yes, a BC that has not been removed from the matrix. > Or only the rows have an entry on the diagonal? What do you mean "strips > it out"? Do you mean it appears in NO aggregate with MIS? > > Yes, because I don't want these as coarse grid variables. They do not need a coarse grid correction. An app can have millions of these and I don't want them around (eg, they are low rank when you have more null space vectors than block size) > Rather than "adding a block on each processor for singletons won't it be > better that MIS doesn't "strip these out" but instead puts them each in > their own little (size 1) aggregate? Then they will automatically get their > own blocks? > I would then have to strip them out again for the prolongator and having a full blown ASM block for every BC vertex would be a mess. I just prefer to strip them out. > > > > > > In addition Fande is adding error checking to PCGASM so if you pass > it badly formatted subdomain information (like was passed from GAMG) it > will generate a very useful error message instead of just chugging along > with gibberish. > > > > Barry > > > > > > Mark, my confusion came from that fact that a single MPI process owns > each of the aggs; that is the list of degrees of freedom for each agg is > all on one process. > > > > NO, NO, NO > > > > This is exactly what PCASM needs but NOT what PCGASM needs. > > > > > > My aggregates span processor subdomains. > > > > The MIS aggregator is such (simple greedy) that an aggregate assigned to > a process can span only one layer of vertices into a neighbor. (The HEM > coarsener is more sophisticated and can deal with Jed's canonical thin > wire, for instance, and can span forever, sort of.) > > > > So the code now is giving you aggregates that span processors (ie, not > local indices). I am puzzled that this works. Am I misunderstanding you? > You are very clear here. Puzzled. > > I mean that ALL the indices for any single aggregate are ALL stored on > the same process; in the same aggregate list! I don't mean that the indices > for an aggregate can only be for variables that live on that process. I > think we are in understanding here, just bad communication. > > OK, good. I will add the fix for BCs and add this to ksp/ex56 and test it. Thanks, Mark > > > > > > I can change this, if ASM can not deal with this. I will just drop > off-processor indices and add non-covered indices to my new singleton > aggregate. > > > > Mark > >

1 0

Re: [petsc-dev] asm / gasm
by Barry Smith 26 Jun '16

26 Jun '16

> On Jun 26, 2016, at 9:28 PM, Mark Adams <mfadams(a)lbl.gov> wrote: > > Thanks Barry, ... we are still not communicating (see below). > > On Sun, Jun 26, 2016 at 11:56 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > I have changed the GAMG aggs support to use PCASM and not PCGASM. I don't see how it ever worked with PCGASM; very strange. Maybe it was written and tested with PCASM and then later someone changed to PCGASM and did not test it again. > > I really thought it was ASM but memory can be strange sometimes ... > > > I was wrong when I said it would be difficult to change to PCASM, it was much simpler to change to PCASM then I thought; this is because I assumed it had been done correctly for PCGASM while in fact though it used PCGASM in the code it followed the model for what PCASM needed :-(. > > Anyways I have a branch barry/fix-gamg-asm-aggs that now works for my ex10 tests I will continue to clean up the branch, fixing naming styles and test with more difficult GAMG examples. And add proper nightly tests for this functionality which it never had before. > > Great thanks. > > Which ex10 is it? and how to I run your tests? src/ksp/ksp/examples/tutorials petscmpiexec -valgrind -n 2 ./ex10 -f0 ~/Datafiles/matrices/poisson2 -ksp_type cg -ksp_monitor_short -ksp_rtol 1.e-8 -pc_type gamg -pc_gamg_type agg -pc_gamg_agg_nsmooths 1 -pc_gamg_coarse_eq_limit 100 -pc_gamg_reuse_interpolation true -pc_gamg_square_graph 1 -pc_gamg_threshold 0.0 -ksp_converged_reason -use_mat_nearnullspace true -mg_levels_ksp_max_it 2 -mg_levels_ksp_type chebyshev -mg_levels_esteig_ksp_type cg -mg_levels_esteig_ksp_max_it 10 -mg_levels_ksp_chebyshev_esteig 0,0.05,0,1.05 -mg_levels_pc_type sor -mat_block_size 1 -pc_gamg_use_agg_asm -mg_levels_pc_type asm the data files come from ftp.mcs.anl.gov in the directory pub/petsc/datafiles/matrices > > I added some code to add a block on each processor for any singletons, because the MIS code strips these (so, yes, not a true MIS). I should do this for users that put a live equation in a singleton, like a non-homo Diri BC. I can add that you your branch. Let me know if that is OK. I do not understand this. Do you mean a variable that is only coupled to itself? So the row and column for the variable have only an entry on the diagonal? Or only the rows have an entry on the diagonal? What do you mean "strips it out"? Do you mean it appears in NO aggregate with MIS? Rather than "adding a block on each processor for singletons won't it be better that MIS doesn't "strip these out" but instead puts them each in their own little (size 1) aggregate? Then they will automatically get their own blocks? > > > In addition Fande is adding error checking to PCGASM so if you pass it badly formatted subdomain information (like was passed from GAMG) it will generate a very useful error message instead of just chugging along with gibberish. > > Barry > > > Mark, my confusion came from that fact that a single MPI process owns each of the aggs; that is the list of degrees of freedom for each agg is all on one process. > > NO, NO, NO > > This is exactly what PCASM needs but NOT what PCGASM needs. > > > My aggregates span processor subdomains. > > The MIS aggregator is such (simple greedy) that an aggregate assigned to a process can span only one layer of vertices into a neighbor. (The HEM coarsener is more sophisticated and can deal with Jed's canonical thin wire, for instance, and can span forever, sort of.) > > So the code now is giving you aggregates that span processors (ie, not local indices). I am puzzled that this works. Am I misunderstanding you? You are very clear here. Puzzled. I mean that ALL the indices for any single aggregate are ALL stored on the same process; in the same aggregate list! I don't mean that the indices for an aggregate can only be for variables that live on that process. I think we are in understanding here, just bad communication. > > I can change this, if ASM can not deal with this. I will just drop off-processor indices and add non-covered indices to my new singleton aggregate. > > Mark

1 0

Re: [petsc-dev] asm / gasm
by Mark Adams 26 Jun '16

26 Jun '16

Thanks Barry, ... we are still not communicating (see below). On Sun, Jun 26, 2016 at 11:56 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > I have changed the GAMG aggs support to use PCASM and not PCGASM. I > don't see how it ever worked with PCGASM; very strange. Maybe it was > written and tested with PCASM and then later someone changed to PCGASM and > did not test it again. > I really thought it was ASM but memory can be strange sometimes ... > > I was wrong when I said it would be difficult to change to PCASM, it was > much simpler to change to PCASM then I thought; this is because I assumed > it had been done correctly for PCGASM while in fact though it used PCGASM > in the code it followed the model for what PCASM needed :-(. > > Anyways I have a branch barry/fix-gamg-asm-aggs that now works for my > ex10 tests I will continue to clean up the branch, fixing naming styles > and test with more difficult GAMG examples. And add proper nightly tests > for this functionality which it never had before. > Great thanks. Which ex10 is it? and how to I run your tests? I added some code to add a block on each processor for any singletons, because the MIS code strips these (so, yes, not a true MIS). I should do this for users that put a live equation in a singleton, like a non-homo Diri BC. I can add that you your branch. Let me know if that is OK. > > In addition Fande is adding error checking to PCGASM so if you pass it > badly formatted subdomain information (like was passed from GAMG) it will > generate a very useful error message instead of just chugging along with > gibberish. > > Barry > > > Mark, my confusion came from that fact that a single MPI process owns > each of the aggs; that is the list of degrees of freedom for each agg is > all on one process. NO, NO, NO > This is exactly what PCASM needs but NOT what PCGASM needs. > > My aggregates span processor subdomains. The MIS aggregator is such (simple greedy) that an aggregate assigned to a process can span only one layer of vertices into a neighbor. (The HEM coarsener is more sophisticated and can deal with Jed's canonical thin wire, for instance, and can span forever, sort of.) So the code now is giving you aggregates that span processors (ie, not local indices). I am puzzled that this works. Am I misunderstanding you? You are very clear here. Puzzled. I can change this, if ASM can not deal with this. I will just drop off-processor indices and add non-covered indices to my new singleton aggregate. Mark

1 0

Re: [petsc-dev] MatMFFDSetType
by Adrian Croucher 26 Jun '16

26 Jun '16

On 27/06/16 13:44, Barry Smith wrote: > Yes, you are correct, there is no MatFDSetType() and in the > SetFromOptions we totally cheat > http://www.mcs.anl.gov/petsc/petsc-current/src/mat/matfd/fdmatrix.c.html#Ma… > and just stick the value directly into the data structure. For now you > can call > PetscOptionsSetValue(NULL,"-mat_fd_type","ds") in your code (before you start using the MatFD stuff) and we'll add the routine to the master branch. Thanks, that workaround does the trick for now. > I didn't think anyone would care about the difference between ds and wp with MatFD. It makes a noticeable difference for the problem I'm solving right now. That may be because the problem is not non-dimensionalized (as yet anyway) and so some of the solution elements are of different magnitudes than others. If I non-dimensionalize them it may possibly not matter any more about wp or ds. - Adrian -- Dr Adrian Croucher Senior Research Fellow Department of Engineering Science University of Auckland, New Zealand email: a.croucher(a)auckland.ac.nz tel: +64 (0)9 923 84611

1 0