[petsc-users] Mat Error
Dharmendar Reddy
dharmareddy84 at gmail.com
Thu Jan 10 13:08:12 CST 2013
Hello,
I am seeing an error in my code when i switched from pets 3.2-p5
to petsc 3.3-p5. Please have a look at the error information below. I do
not pass any petsc specific command line options.
The following is the subroutine which leads to error message.
subroutine initBCSolver_solver(this,comm,SolverType,NumDof,ierr,EqnCtx)
use SolverUtils_m
implicit none
#include "finclude/petsc.h"
class(Solver_t), intent(inout) :: this
MPI_Comm :: comm
integer, intent(in) :: SolverType
integer, intent(in) :: NumDof
class(*),pointer,intent(inout),optional :: EqnCtx
PetscErrorCode :: ierr
PetscReal :: rtol
PetscReal :: abstol
PetscReal :: dtol
PetscInt :: maxits
! Local Variables
PetscInt :: psize ! problem size
PetscScalar :: pfive !
! Set solver type
!if(SolverType /= LINEAR .and. SolverType /= NonLinear) then
! print*,'Error:: Solver must be linear or Nonlinear',solverType
! stop
!end if
this%comm = comm
this%SolverType = SolverType
! Set problem size
this%numDof = numDof
psize = this%numDof ! 3
! intiate Storage vectors
call VecCreateSeq(this%comm,psize,this%sol,ierr)
call VecSetOption(this%sol, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE,ierr)
call VecDuplicate(this%sol,this%fxn,ierr)
this%flag_vec_sol = .true.
this%flag_vec_fxn = .true.
select case(this%SolverType)
case(Linear)
! create the linear operator
call MatCreate(this%comm,this%A,ierr)
call MatSetSizes(this%A,PETSC_DECIDE,PETSC_DECIDE,psize,psize,ierr)
call MatSetFromOptions(this%A,ierr)
! create storage for rhs
call VecDuplicate(this%sol,this%rhs,ierr)
! Create KSP context
call KSPCreate(this%comm,this%ksp,ierr)
this%flag_ksp_ksp = .true.
rtol = PETSC_DEFAULT_DOUBLE_PRECISION
abstol = PETSC_DEFAULT_DOUBLE_PRECISION
dtol = PETSC_DEFAULT_DOUBLE_PRECISION
maxIts = PETSC_DEFAULT_INTEGER
call KSPSetTolerances(this%ksp,rtol,abstol,dtol,maxIts,ierr);
call KSPSetFromOptions(this%ksp,ierr)
call
KSPSetOperators(this%ksp,this%A,this%A,DIFFERENT_NONZERO_PATTERN,ierr)
case(nonLinear)
call MatCreate(this%comm,this%Jac,ierr)
this%flag_mat_Jac = .true.
call
MatSetSizes(this%Jac,PETSC_DECIDE,PETSC_DECIDE,psize,psize,ierr)
call MatSetFromOptions(this%Jac,ierr)
! create snes context
call SNESCreate(this%comm,this%snes,ierr)
this%flag_snes_snes = .true.
! Set function evaluation routine and vector
! use the defualt FormFunction_snes
call SNESSetFunction(this%snes,this%fxn,FormFunction_snes,this,ierr)
! Set Jacobian matrix data structure and Jacobian evaluation
routine
call
SNESSetJacobian(this%snes,this%Jac,this%Jac,FormJacobian_snes,this,ierr)
!call
SNESSetJacobian(this%snes,this%Jac,this%Jac,PETSC_NULL_FUNCTION,PETSC_NULL_OBJECT,ierr)
call
SNESSetTolerances(this%snes,1e-15*25.6E-3,PETSC_DEFAULT_DOUBLE_PRECISION,PETSC_DEFAULT_DOUBLE_PRECISION,25000,100000,ierr);
call SNESSetFromOptions(this%snes,ierr)
end select
!call VecDuplicate(this%sol,this%bOhmic,ierr)
pfive = 0.0
call VecSet(this%sol,pfive,ierr) ! set Intial Guess
end subroutine initBCSolver_solver
end module BoundaryConditionUtils_m
[0]PETSC ERROR: --------------------- Error Message
------------------------------------
[0]PETSC ERROR: Object is in wrong state!
[0]PETSC ERROR: Must call MatXXXSetPreallocation() or MatSetUp() on
argument 1 "mat" before MatGetFactorAvailable()!
[0]PETSC ERROR:
------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.3.0, Patch 5, Sat Dec 1 15:10:41
CST 2012
[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:
/home1/00924/Reddy135/projects/utgds/test/IIIVTunnelFET2D/PoisTest on a
sandybrid named login2.stampede.tacc.utexas.edu by Reddy135 Thu Jan 10
12:50:19 2013
[0]PETSC ERROR: Libraries linked from
/opt/apps/intel13/mvapich2_1_9/petsc/3.3/sandybridge-debug/lib
[0]PETSC ERROR: Configure run at Mon Jan 7 14:42:13 2013
[0]PETSC ERROR: Configure options --with-x=0 -with-pic
--with-external-packages-dir=/opt/apps/intel13/mvapich2_1_9/petsc/3.3/externalpackages
--with-mpi-compilers=1 --with-mpi-dir=/opt/apps/intel13/mvapich2/1.9
--with-scalar-type=real --with-dynamic-loading=0 --with-shared-libraries=1
--with-spai=1 --download-spai --with-hypre=1 --download-hypre
--with-mumps=1 --download-mumps --with-scalapack=1 --download-scalapack
--with-blacs=1 --download-blacs --with-spooles=1 --download-spooles
--with-superlu=1 --download-superlu --with-superlu_dist=1
--download-superlu_dist --with-parmetis=1 --download-parmetis
--with-metis=1 --download-metis --with-debugging=yes
--with-blas-lapack-dir=/opt/apps/intel/13/composer_xe_2013.1.117/mkl/lib/intel64
--with-mpiexec=mpirun_rsh --COPTFLAGS= --CXXOPTFLAGS= --FOPTFLAGS=
[0]PETSC ERROR:
------------------------------------------------------------------------
[0]PETSC ERROR: MatGetFactorAvailable() line 3921 in
/opt/apps/intel13/mvapich2_1_9/petsc/3.3/src/mat/interface/matrix.c
[0]PETSC ERROR: PCGetDefaultType_Private() line 26 in
/opt/apps/intel13/mvapich2_1_9/petsc/3.3/src/ksp/pc/interface/precon.c
[0]PETSC ERROR: PCSetFromOptions() line 181 in
/opt/apps/intel13/mvapich2_1_9/petsc/3.3/src/ksp/pc/interface/pcset.c
[0]PETSC ERROR: KSPSetFromOptions() line 287 in
/opt/apps/intel13/mvapich2_1_9/petsc/3.3/src/ksp/ksp/interface/itcl.c
[0]PETSC ERROR: SNESSetFromOptions() line 678 in
/opt/apps/intel13/mvapich2_1_9/petsc/3.3/src/snes/interface/snes.c
--
-----------------------------------------------------
Dharmendar Reddy Palle
Graduate Student
Microelectronics Research center,
University of Texas at Austin,
10100 Burnet Road, Bldg. 160
MER 2.608F, TX 78758-4445
e-mail: dharmareddy84 at gmail.com
Phone: +1-512-350-9082
United States of America.
Homepage: https://webspace.utexas.edu/~dpr342
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130110/9561b514/attachment.html>
More information about the petsc-users
mailing list