If you want to use the PETSc to compute the jacobian via finite differencing then here's the sequence of calls
ierr = SNESComputeDefaultJacobian(ts_snes,CV_Y,&J,&J,&flag);CHKERRQ(ierr); ierr = MatGetColoring(J,MATCOLORINGSL,&iscoloring);CHKERRQ(ierr); ierr = MatFDColoringCreate(J,iscoloring,&matfdcoloring);CHKERRQ(ierr); ierr = ISColoringDestroy(iscoloring);CHKERRQ(ierr); ...you missed this ierr = MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode (*)(void))f,(void*)&appctx);CHKERRQ(ierr); ierr = MatFDColoringSetFromOptions(matfdcoloring);CHKERRQ(ierr); ierr = TSSetRHSJacobian(ts,J,J,TSDefaultComputeJacobianColor,matfdcoloring); if you are using your own jacobian evaluation then use SNESSetJacobian followed by SNESComputeJacobian instead of SNESComputeDefaultJacobian. Shri ----- "Xuan YU" <xxy113 at psu.edu> wrote: > ierr = MatCreateSeqAIJ(PETSC_COMM_SELF,N,N,10,PETSC_NULL,&J);CHKERRQ(ierr); ierr = MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); ierr = MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); ierr = SNESComputeJacobian(ts_snes,CV_Y,&J,&J,&flag);CHKERRQ(ierr); ierr = MatGetColoring(J,MATCOLORINGSL,&iscoloring);CHKERRQ(ierr); ierr = MatFDColoringCreate(J,iscoloring,&matfdcoloring);CHKERRQ(ierr); ierr = MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode (*)(void))f,(void*)&appctx);CHKERRQ(ierr); ierr = MatFDColoringSetFromOptions(matfdcoloring);CHKERRQ(ierr); ierr = TSSetRHSJacobian(ts,J,J,TSDefaultComputeJacobianColor,matfdcoloring); > These are the Jacobian related codes. > > > > On Jul 7, 2010, at 1:51 PM, Satish Balay wrote: total: nonzeros=1830 > mallocs used during MatSetValues calls =1830 > > Looks like you are zero-ing out the non-zero structure - before > assembling the matrix. > > Are you calling MatZeroRows() or MatZeroEntries() or something else - > before assembling the matrix? > > Satish > > On Wed, 7 Jul 2010, Xuan YU wrote: > > I made a change: ierr = > MatCreateSeqAIJ(PETSC_COMM_SELF,N,N,5,PETSC_NULL,&J);CHKERRQ(ierr); > > Time of the code did not change much, and got the info: > Matrix Object: > type=seqaij, rows=1830, cols=1830 > total: nonzeros=1830, allocated nonzeros=36600 > total number of mallocs used during MatSetValues calls =1830 > not using I-node routines > > > > On Jul 7, 2010, at 12:51 PM, Satish Balay wrote: > > total: nonzeros=1830, allocated nonzeros=29280 > total number of mallocs used during MatSetValues calls =1830 > > There is something wrong with your preallocation or matrix > assembly. You should see zero mallocs for efficient assembly. > > http://www.mcs.anl.gov/petsc/petsc-as/documentation/faq.html#efficient-assembly > > satish > > > On Wed, 7 Jul 2010, Xuan YU wrote: > > Hi, > > I finite difference Jacobian approximation for my TS model. The size of > the > vector is 1830. I got the following info with(-ts_view): > > type: beuler > maximum steps=50 > maximum time=50 > total number of nonlinear solver iterations=647 > total number of linear solver iterations=647 > SNES Object: > type: ls > line search variant: SNESLineSearchCubic > alpha=0.0001, maxstep=1e+08, minlambda=1e-12 > maximum iterations=50, maximum function evaluations=10000 > tolerances: relative=1e-08, absolute=1e-50, solution=1e-08 > total number of linear solver iterations=50 > total number of function evaluations=51 > KSP Object: > type: gmres > GMRES: restart=30, using Classical (unmodified) Gram-Schmidt > Orthogonalization with no iterative refinement > GMRES: happy breakdown tolerance 1e-30 > maximum iterations=10000, initial guess is zero > tolerances: relative=1e-05, absolute=1e-50, divergence=10000 > left preconditioning > using PRECONDITIONED norm type for convergence test > PC Object: > type: ilu > ILU: out-of-place factorization > 0 levels of fill > tolerance for zero pivot 1e-12 > using diagonal shift to prevent zero pivot > matrix ordering: natural > factor fill ratio given 1, needed 1 > Factored matrix follows: > Matrix Object: > type=seqaij, rows=1830, cols=1830 > package used to perform factorization: petsc > total: nonzeros=1830, allocated nonzeros=1830 > total number of mallocs used during MatSetValues calls =0 > not using I-node routines > linear system matrix = precond matrix: > Matrix Object: > type=seqaij, rows=1830, cols=1830 > total: nonzeros=1830, allocated nonzeros=29280 > total number of mallocs used during MatSetValues calls =1830 > not using I-node routines > > > 50 output time step takes me 11.877s. So I guess there is something not > appropriate with my Jacobian Matrix. Could you please tell me how to speed > up > my code? > > Thanks! > > Xuan YU > xxy113 at psu.edu > > > > > > > > Xuan YU (??) > xxy113 at psu.edu > > > > > > > > Xuan YU ( ?? ) xxy113 at psu.edu > > -------------- next part -------------- An HTML attachment was scrubbed... URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20100707/31079bb5/attachment.htm>
