[petsc-users] my code runs too slow
Xuan YU
xxy113 at psu.edu
Wed Jul 7 13:54:15 CDT 2010
Reordering does work. But it went back to beginning. Info is shown
below. I want to follow the example(http://www.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-2.3.3/src/ts/examples/tutorials/ex7.c.html
) ierr =
TSSetRHSJacobian
(ts,J,J,TSDefaultComputeJacobianColor,matfdcoloring);CHKERRQ(ierr); I
thought this code can compute Jacobian. Actually, I have not write the
Jacobian function yet. Does it necessary to speed up the code?
TS Object:
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=45750
total number of mallocs used during MatSetValues calls =1830
not using I-node routines
On Jul 7, 2010, at 2:36 PM, Satish Balay wrote:
> What function are you using to compute jacobian? Perhaps you are
> missing MatAssembly() calls in that function?
>
>
> Or reordering the code as follows should work - but the correct
> location for MatAssembly() calls are in the routine you are using for
> assembling the jacobian.
>
> Satish
>
>>>>>>>>>>
>
> ierr = MatCreateSeqAIJ(PETSC_COMM_SELF,N,N,
> 10,PETSC_NULL,&J);CHKERRQ(ierr);
> ierr = SNESComputeJacobian(ts_snes,CV_Y,&J,&J,&flag);CHKERRQ(ierr);
>
> ierr = MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
> ierr = MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);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);
>
> <<<<<<
>
>
>
>
> On Wed, 7 Jul 2010, Xuan YU wrote:
>
>> Sorry, this is the complete one:
>>
>> [0]PETSC ERROR: --------------------- Error Message
>> ------------------------------------
>> [0]PETSC ERROR: Object is in wrong state!
>> [0]PETSC ERROR: Not for unassembled matrix!
>> [0]PETSC ERROR:
>> ------------------------------------------------------------------------
>> [0]PETSC ERROR: Petsc Development HG revision:
>> fdc7be12de37b8a400b6e561de9850bcc7e79f4f HG Date: Fri Jul 02
>> 21:17:36 2010
>> -0500
>> [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: ./pihm on a arch-linu named lionxi.rcc.psu.edu by
>> xxy113 Wed
>> Jul 7 14:32:11 2010
>> [0]PETSC ERROR: Libraries linked from
>> /gpfs/home/xxy113/soft/petsc-dev/arch-linux-gnu-c-debug/lib
>> [0]PETSC ERROR: Configure run at Sun Jul 4 14:09:16 2010
>> [0]PETSC ERROR: Configure options --download-f-blas-lapack=1
>> --download-mpich=1
>> [0]PETSC ERROR:
>> ------------------------------------------------------------------------
>> [0]PETSC ERROR: MatGetColoring() line 491 in src/mat/color/color.c
>> [0]PETSC ERROR: VecRestoreArrayPrivate3() line 201 in /work/petsc/
>> pihm.c
>>
>> On Jul 7, 2010, at 2:28 PM, Satish Balay wrote:
>>
>>>>> ierr =
>>>>> SNESComputeJacobian(ts_snes,CV_Y,&J,&J,&flag);CHKERRQ(ierr);
>>>
>>> Perhaps there is a function thats set to compute jacobinan thats not
>>> assembling the matrix properly?
>>>
>>> Also - please dont' truncate error messages - when you send
>>> them. Incomplete info is not useful.
>>>
>>> Satish
>>>
>>> On Wed, 7 Jul 2010, Matthew Knepley wrote:
>>>
>>>> You have to assemble
>>>>
>>>> a) after setting values
>>>>
>>>> b) before using the matrix
>>>>
>>>> Please consult the user examples where this is done correctly and
>>>> the
>>>> manual
>>>> section which explains the assembly process.
>>>>
>>>> Matt
>>>>
>>>> On Wed, Jul 7, 2010 at 8:13 PM, Xuan YU <xxy113 at psu.edu> wrote:
>>>>
>>>>>
>>>>> On Jul 7, 2010, at 2:06 PM, Satish Balay wrote:
>>>>>
>>>>> On Wed, 7 Jul 2010, Xuan YU 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);
>>>>>
>>>>>
>>>>> This assembly removes the unused space here. Since no values are
>>>>> inserted - it squezes out all of the allocated space. Perhaps
>>>>> you just
>>>>> need to remove these 2 calls as the actual matrix is assembled
>>>>> further
>>>>> down the code.
>>>>>
>>>>>
>>>>> I removed these 2
>>>>> But got Error Message
>>>>>
>>>>> [0]PETSC ERROR: Object is in wrong state!
>>>>> [0]PETSC ERROR: Not for unassembled matrix!
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> Satish
>>>>>
>>>>> 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
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> 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/e44e44a4/attachment-0001.htm>
More information about the petsc-users
mailing list