Skip to content

Coloring#160

Merged
bendudson merged 12 commits intomasterfrom
coloring
Feb 22, 2016
Merged

Coloring#160
bendudson merged 12 commits intomasterfrom
coloring

Conversation

@bendudson
Copy link
Contributor

Adds support for algebraic constraints and finite difference Jacobian approximation using coloring to the IMEX-BDF2 solver

New example IMEX/drift-wave-constraint

Some documentation in README.md files and the user manual

Would allow PETSc to use coloring to speed up finite
difference calculation of the system Jacobian.
Seems to work well for the IMEX/diffusion-nl problem

Doesn't converge for the IMEX/drift-wave problem.
If model was not split, then previously the RHS function
was treated as a convective (explicit) part.

Now switched so that unsplit models are treated as
entirely diffusive (implicit)
When using coloring to calculate the Jacobian,
re-use the Jacobian 4 times by default.
o IMEX-BDF2 solver can now handle constraints as part
  of the implicit solve. This can be used to solve
  for potential, for example.

o Added examples/IMEX/drift-wave-constraint
  which solves for the potential phi as a constraint
  rather than inside the RHS.

o Moved generic and duplicated variable loading code
  from IDA solver to Solver base class.

o Changed error handling in Solver::save_vars().
  Now throws exception with useful error message, rather than
  returning error code.
The linear flag in the diffusive() RHS function call
is now set when the Jacobian is being calculated
using coloring.
Conflicts:
	src/solver/impls/imex-bdf2/imex-bdf2.cxx
Added a small section to the user manual
on the IMEX-BDF2 solver, and added some example
code for solving with constraints.
@d7919
Copy link
Member

d7919 commented Jan 29, 2016

I think line 590 of imex-bdf2.cxx may need modifying to support PETSc<3.4.

In 3.4 SNESDefaultComputeJacobianColor was renamed
to SNESComputeJacobianDefaultColor
@d7919
Copy link
Member

d7919 commented Feb 3, 2016

I find that using the matrix_free=true approach leads to errors unless I specify -pc_type on the run line to set the preconditioner to none. I think it's possible to pass in a different matrix to SNESSetJacobian for use with preconditioning but I think probably all we want to do is set the PETSc default preconditioner to none if matrix_free=true and use_precon=false || !hasUserPrecon(). This would be around line 660 I think.

@bendudson
Copy link
Contributor Author

Hi David,

I found it worked with preconditioning. Is this for the IMEX/diffusion-nl
and IMEX/drift-wave-constraint examples, because those were working well
for me. Is it a PETSc version issue, as I was testing with 3.4.?

Ben

On 3 February 2016 at 09:46, David Dickinson notifications@github.com
wrote:

I find that using the matrix_free=true approach leads to errors unless I
specify -pc_type on the run line to set the preconditioner to none. I
think it's possible to pass in a different matrix to SNESSetJacobian for
use with preconditioning but I think probably all we want to do is set the
PETSc default preconditioner to none if matrix_free=true and use_precon=false
|| !hasUserPrecon(). This would be around line 660 I think.


Reply to this email directly or view it on GitHub
#160 (comment).

@d7919
Copy link
Member

d7919 commented Feb 3, 2016

This is actually for the gyro-fluid model that I'm working on. I'm using v3.4.5 and I'm running with solver:matrix_free=true, solver:use_precon=false and -info and I get the following output

[0] PCSetUp(): Setting up new PC
[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: No support for this operation for this object type!
[0]PETSC ERROR: Matrix type mffd does not support getting diagonal block!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.4.5, Jun, 29, 2014 
[0]PETSC ERROR: See docs/changes/index.html for recent updates.
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
[0]PETSC ERROR: See docs/index.html for manual pages.
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Libraries linked from src/petsc-3.4.5/arch-linux2-cxx-opt_REAL/lib
[0]PETSC ERROR: Configure run at Mon Jul 13 12:39:41 2015
[0]PETSC ERROR: Configure options  --with-clanguage=cxx --download-ptscotch --download-pastix --download-superlu_dist --download-parmetis --download-metis --download-mumps --download-hypre --download-scalapack --download-blacs --with-fftw --with-fftw-dir=/usr --with-netcdf --with-netcdf-dir=/usr --with-mpi=yes --with-precision=double --with-scalar-type=real --with-shared-libraries=0 --with-debugging=no --with-make-np=64 --COPTFLAGS="-O3 -flto=64 -fwhole-program " --CXXOPTFLAGS="-O3 -flto=64 -fwhole-program " --FOPTFLAGS=-O3
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: MatGetDiagonalBlock() line 161 in src/mat/interface/matrix.c
[0]PETSC ERROR: PCSetUp_BJacobi() line 126 in src/ksp/pc/impls/bjacobi/bjacobi.c
[0]PETSC ERROR: PCSetUp() line 890 in src/ksp/pc/interface/precon.c
[0]PETSC ERROR: KSPSetUp() line 278 in src/ksp/ksp/interface/itfunc.c
[0]PETSC ERROR: KSPSolve() line 399 in src/ksp/ksp/interface/itfunc.c
[0]PETSC ERROR: SNESSolve_NEWTONLS() line 220 in src/snes/impls/ls/ls.c
[0]PETSC ERROR: SNESSolve() line 3636 in src/petsc-3.4.5/src/snes/interface/snes.c

@d7919
Copy link
Member

d7919 commented Feb 3, 2016

Ah sorry, it looks like this was from an unclean build.

Also introduces higher order schemes (up to fourth order).

This commit involves a substantial restructuring of the original
imex-bdf code in order to make it easier to change the order of the
scheme.

Introduced imex-bdf3 and imex-bdf4 schemes.

Adaptivity possible (set adaptive = true). This requires use of
significantly more complicated coefficients (see code for reference) but
general structure of scheme essentially unchanged. Error estimated by
comparing solution with requested order scheme and one order
lower. Rough argument then used to decide how to change the
timestep. This argument works ok for requested order >=3, a little bit
trickier when order = 2.

Some things to note are:
1. Changing the timestep frequently can significantly increase the
number of implicit solves required. The nadapt option aims to help
control this.
2. Recommend order >=3 when adaptive.

Whilst this commit introduces a fully tested and working feature there
are some outstanding elements that are intended to be addressed in the
future:
1. Add documentation of all input parameters (comments in code should be
enough for now).
2. Catch failures of the SNESSolve call to allow a repeat with a smaller
timestep (if desired) rather than aborting as done currently.
3. Develop automatic decision making for how often to check the error
when adaptive.
4. Investigate alternative schemes (e.g. TVB rather than BDF).

Some outstanding issues/concerns also exist:
1. MMS time test (non-adaptive) gives 2nd order for 3rd and 4th order
schemes in some situations (depends on relative size of explicit and
implicit ddt values).
2. Use of two SNES solve object could possibly lead to a drift in
results used for the error checking over time (most likely leading to a
slightly smaller timestep than needed). This is because if the KSP
tolerances are not tight enough the SNES solve can be slightly dependent
on the previous result, which means the high and low order solvers, which
have slightly different solutions from the previous stage, could give
different results.
3. Adaptivity may struggle if the explicit scheme is limiting the
timestep. This isn't tested but could be a slight issue. May want to
introduce a cfl like parameter which allows us to add some safety margin
to our timestep estimate.
bendudson added a commit that referenced this pull request Feb 22, 2016
Major changes/improvements to the IMEX-BDF2 solver

- Adaptive timestepping and higher order schemes (BDF3)
- Coloring for faster finite difference Jacobian calculation
@bendudson bendudson merged commit 1a6346f into master Feb 22, 2016
@d7919 d7919 deleted the coloring branch October 19, 2016 12:15
@d7919 d7919 restored the coloring branch October 19, 2016 12:15
@d7919 d7919 deleted the coloring branch October 19, 2016 12:15
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants