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