On Fri, Feb 10, 2012 at 19:51, Dmitry Karpeev <karpeev@mcs.anl.gov> wrote:

How does this explanation square with the following code from KSPSetUp?

 if (ksp->dmActive && !ksp->setupstage) {

    /* first time in so build matrix and vector data structures using DM */

    if (!ksp->vec_rhs) {ierr = DMCreateGlobalVector(ksp->dm,&ksp->vec_rhs);CHKERRQ(ierr);}

    if (!ksp->vec_sol) {ierr = DMCreateGlobalVector(ksp->dm,&ksp->vec_sol);CHKERRQ(ierr);}

    ierr = DMCreateMatrix(ksp->dm,MATAIJ,&A);CHKERRQ(ierr);

    ierr = KSPSetOperators(ksp,A,A,stflg);CHKERRQ(ierr);  

    ierr = PetscObjectDereference((PetscObject)A);CHKERRQ(ierr); 

 }

Look down the page a little, DMComputeJacobian() is what actually assembles.

    if (ksp->setupstage != KSP_SETUP_NEWRHS) {

      ierr = KSPGetOperators(ksp,&A,&A,PETSC_NULL);CHKERRQ(ierr);

      ierr = DMComputeJacobian(ksp->dm,PETSC_NULL,A,A,&stflg);CHKERRQ(ierr);

      ierr = KSPSetOperators(ksp,A,A,stflg);CHKERRQ(ierr); 

    }