[petsc-users] Questions abt ex22f

TAY wee-beng zonexo at gmail.com
Sun Apr 22 14:18:55 CDT 2012


Hi,

I have attached the ex22f.F file. The changes I added are given in bold:

...

PetscErrorCode   ierr
       DM               da
       KSP              ksp
*Vec              x,b*
       external         ComputeRHS,ComputeMatrix

....

call KSPSetUp(ksp,ierr)
       call KSPSolve(ksp,PETSC_NULL_OBJECT,PETSC_NULL_OBJECT,ierr)
*call KSPGetSolution(ksp,x,ierr)
       call VecView(x,ierr)*
       call KSPDestroy(ksp,ierr)
       call DMDestroy(da,ierr)

The error is:

[0]PETSC ERROR: 
------------------------------------------------------------------------
[0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, 
probably memory access out of range
[0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
[0]PETSC ERROR: or see 
http://www.mcs.anl.gov/petsc/petsc-as/documentation/faq.html#valgrind[0]PETSC 
ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to 
find memory corruption errors
[0]PETSC ERROR: likely location of problem given in stack below
[0]PETSC ERROR: ---------------------  Stack Frames 
------------------------------------
[0]PETSC ERROR: Note: The EXACT line numbers in the stack are not available,
[0]PETSC ERROR:       INSTEAD the line number of the start of the function
[0]PETSC ERROR:       is given.
[0]PETSC ERROR: [0] VecView line 735 
/home/wtay/Codes/petsc-3.2-p5/src/vec/vec/interface/vector.c
[0]PETSC ERROR: --------------------- Error Message 
------------------------------------
[0]PETSC ERROR: Signal received!
[0]PETSC ERROR: 
------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.2.0, Patch 5, Sat Oct 29 
13:45:54 CDT 2011
[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: ./ex22f on a arch-linu named hpc12 by wtay Sun Apr 22 
21:11:39 2012
[0]PETSC ERROR: Libraries linked from 
/home/wtay/Lib/petsc-3.2-p5_mumps_debug/lib
[0]PETSC ERROR: Configure run at Sun Nov 27 15:39:26 2011
[0]PETSC ERROR: Configure options --with-mpi-dir=/opt/openmpi-1.5.3/ 
--with-blas-lapack-dir=/opt/intel_xe_2011/mkl/lib/intel64/ 
--with-debugging=1 --download-hypre=1 
--prefix=/home/wtay/Lib/petsc-3.2-p5_mumps_debug COPTFLAGS=-O0 
FOPTFLAGS=-O0 --download-mumps=1 --download-parmetis=1 
--download-scalapack=1 --download-blacs=1
[0]PETSC ERROR: 
------------------------------------------------------------------------
[0]PETSC ERROR: User provided function() line 0 in unknown directory 
unknown file
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD
with errorcode 59.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.

Yours sincerely,

TAY wee-beng


On 22/4/2012 9:06 PM, Jed Brown wrote:
>
> Run in a debugger and/or use --with-debugging=1 so that the error 
> trace has more information. You could also show us the exact code that 
> you used.
>
> On Apr 22, 2012 2:03 PM, "TAY wee-beng" <zonexo at gmail.com 
> <mailto:zonexo at gmail.com>> wrote:
>
>     Hi,
>
>     I added "Vec x,b" after "KSP ksp"
>     and then "call KSPGetSolution(ksp, x, ierr)"
>
>     I wanted to see the output so I added "call VecView(x,ierr)" but I
>     got this error:
>
>     [0]PETSC ERROR:
>     ------------------------------------------------------------------------
>     [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation
>     Violation, probably memory access out of range
>     [0]PETSC ERROR: Try option -start_in_debugger or
>     -on_error_attach_debugger
>     [0]PETSC ERROR: or see
>     http://www.mcs.anl.gov/petsc/petsc-as/documentation/faq.html#valgrind[0]PETSC
>     ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X
>     to find memory corruption errors
>     [0]PETSC ERROR: configure using --with-debugging=yes, recompile,
>     link, and run
>     [0]PETSC ERROR: to get more information on the crash.
>     [0]PETSC ERROR: --------------------- Error Message
>     ------------------------------------
>     [0]PETSC ERROR: Signal received!
>     [0]PETSC ERROR:
>     ------------------------------------------------------------------------
>     [0]PETSC ERROR: Petsc Release Version 3.2.0, Patch 5, Sat Oct 29
>     13:45:54 CDT 2011
>     [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: ./ex22f on a arch-linu named hpc12 by wtay Sun Apr
>     22 21:02:14 2012
>     [0]PETSC ERROR: Libraries linked from
>     /home/wtay/Lib/petsc-3.2-p5_mumps_rel/lib
>     [0]PETSC ERROR: Configure run at Sun Nov 27 15:18:15 2011
>     [0]PETSC ERROR: Configure options
>     --with-mpi-dir=/opt/openmpi-1.5.3/
>     --with-blas-lapack-dir=/opt/intel_xe_2011/mkl/lib/intel64/
>     --with-debugging=0 --download-hypre=1
>     --prefix=/home/wtay/Lib/petsc-3.2-p5_mumps_rel COPTFLAGS=-O3
>     FOPTFLAGS=-O3 --download-mumps=1 --download-parmetis=1
>     --download-scalapack=1 --download-blacs=1
>     [0]PETSC ERROR:
>     ------------------------------------------------------------------------
>     [0]PETSC ERROR: User provided function() line 0 in unknown
>     directory unknown file
>     --------------------------------------------------------------------------
>     MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD
>     with errorcode 59.
>
>
>     Yours sincerely,
>
>     TAY wee-beng
>
>
>     On 22/4/2012 2:53 PM, Matthew Knepley wrote:
>>     On Sun, Apr 22, 2012 at 3:31 AM, TAY wee-beng <zonexo at gmail.com
>>     <mailto:zonexo at gmail.com>> wrote:
>>
>>         Hi,
>>
>>         I am using petsc-dev 2012-04-20.
>>
>>         Btw, I'm referring to :
>>
>>         http://www.mcs.anl.gov/petsc/petsc-dev/src/ksp/ksp/examples/tutorials/ex22f.F.html
>>
>>         Part of the code is :
>>
>>         call KSPSetFromOptions(ksp,ierr)
>>                call KSPSetUp(ksp,ierr)
>>                call KSPSolve(ksp,PETSC_NULL_OBJECT,PETSC_NULL_OBJECT,ierr)
>>                call KSPDestroy(ksp,ierr)
>>                call DMDestroy(da,ierr)
>>                call PetscFinalize(ierr)
>>
>>
>>
>>         Unlike other codes like ex29c or ex45c, there isn't a "call
>>         KSPGetSolution(ksp,x,ierr)"
>>
>>
>>     You need to declare "Vec x", and then you can call
>>     KSPGetSolution(ksp, x, ierr)
>>
>>        Matt
>>
>>         Also I want to add "call VecView(x,ierr)" to print out the
>>         results, which is usally added after the above.
>>
>>         Thank you
>>
>>         Yours sincerely,
>>
>>         TAY wee-beng
>>
>>
>>         On 22/4/2012 1:14 AM, Matthew Knepley wrote:
>>>         On Sat, Apr 21, 2012 at 6:31 PM, TAY wee-beng
>>>         <zonexo at gmail.com <mailto:zonexo at gmail.com>> wrote:
>>>
>>>             Hi,
>>>
>>>             May I know if ex22f is complete? I can't find :
>>>
>>>             call KSPGetSolution(ksp,x,ierr)
>>>
>>>             If I entered it, it says x not found.
>>>
>>>
>>>         This is correct in petsc-dev. What version are you using?
>>>
>>>           Thanks,
>>>
>>>             Matt
>>>
>>>             Thank you!
>>>
>>>             -- 
>>>             Yours sincerely,
>>>
>>>             TAY wee-beng
>>>
>>>
>>>
>>>
>>>         -- 
>>>         What most experimenters take for granted before they begin
>>>         their experiments is infinitely more interesting than any
>>>         results to which their experiments lead.
>>>         -- Norbert Wiener
>>
>>
>>
>>
>>     -- 
>>     What most experimenters take for granted before they begin their
>>     experiments is infinitely more interesting than any results to
>>     which their experiments lead.
>>     -- Norbert Wiener
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20120422/e6cf1dee/attachment-0001.htm>
-------------- next part --------------
!
!   Laplacian in 3D. Modeled by the partial differential equation
!
!   Laplacian u = 1,0 < x,y,z < 1,
!
!   with boundary conditions
!
!   u = 1 for x = 0, x = 1, y = 0, y = 1, z = 0, z = 1.
!
!   This uses multigrid to solve the linear system

      program main
      implicit none

! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!                    Include files
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!
!     petscsys.h  - base PETSc routines      petscvec.h - vectors
!     petscmat.h - matrices
!     petscksp.h    - Krylov subspace methods  petscpc.h  - preconditioners

#include <finclude/petscsys.h>
#include <finclude/petscvec.h>
#include <finclude/petscmat.h>
#include <finclude/petscpc.h>
#include <finclude/petscksp.h>
#include <finclude/petscdmda.h>
      PetscErrorCode   ierr
      DM               da
      KSP              ksp
      Vec              x,b
      external         ComputeRHS,ComputeMatrix
      PetscInt i1,i3

      call  PetscInitialize(PETSC_NULL_CHARACTER,ierr)

      i3 = -3
      i1 = 1
      call KSPCreate(MPI_COMM_WORLD,ksp,ierr)
      call DMDACreate3d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,            &
     &    DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,                               &
     &    DMDA_STENCIL_STAR,i3,i3,i3,PETSC_DECIDE,PETSC_DECIDE,                        &
     &    PETSC_DECIDE,i1,i1,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                        &
     &                PETSC_NULL_INTEGER,da,ierr)  
      call DMSetFunction(da,ComputeRHS,ierr)
      call DMSetJacobian(da,ComputeMatrix,ierr)   
      call KSPSetDM(ksp,da,ierr)

      call KSPSetFromOptions(ksp,ierr)
      call KSPSetUp(ksp,ierr)
      call KSPSolve(ksp,PETSC_NULL_OBJECT,PETSC_NULL_OBJECT,ierr)
      call KSPGetSolution(ksp,x,ierr)
      call VecView(x,ierr)
      call KSPDestroy(ksp,ierr)
      call DMDestroy(da,ierr)
      call PetscFinalize(ierr)

      end

      subroutine ComputeRHS(da,x,b,ierr)
      implicit none

#include <finclude/petscsys.h>
#include <finclude/petscvec.h>
#include <finclude/petscdmda.h>

      PetscErrorCode  ierr
      PetscInt mx,my,mz
      PetscScalar  h
      Vec          x,b
      DM           da

      call DMDAGetInfo(da,PETSC_NULL_INTEGER,mx,my,mz,                        &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                 &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                 &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                 &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                 &
     &               PETSC_NULL_INTEGER,ierr)
      h    = 1.d0/((mx-1)*(my-1)*(mz-1))

      call VecSet(b,h,ierr)
      return
      end

    
      subroutine ComputeMatrix(da,x,JJ,jac,str,ierr)
      implicit none

#include <finclude/petscsys.h>
#include <finclude/petscvec.h>
#include <finclude/petscmat.h>
#include <finclude/petscdmda.h>

      Mat          jac,JJ
      PetscErrorCode    ierr

      DM           da
      PetscInt    i,j,k,mx,my,mz,xm
      PetscInt    ym,zm,xs,ys,zs,i1,i7
      PetscScalar  v(7),Hx,Hy,Hz
      PetscScalar  HxHydHz,HyHzdHx
      PetscScalar  HxHzdHy
      MatStencil   row(4),col(4,7)
      Vec          x
      MatStructure str
      i1 = 1
      i7 = 7
      call DMDAGetInfo(da,PETSC_NULL_INTEGER,mx,my,mz,                       &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                &
     &               PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,                &
     &               PETSC_NULL_INTEGER,ierr)

      Hx = 1.d0 / (mx-1)
      Hy = 1.d0 / (my-1)
      Hz = 1.d0 / (mz-1)
      HxHydHz = Hx*Hy/Hz
      HxHzdHy = Hx*Hz/Hy
      HyHzdHx = Hy*Hz/Hx
      call DMDAGetCorners(da,xs,ys,zs,xm,ym,zm,ierr)        
  
      do 10,k=zs,zs+zm-1    
        do 20,j=ys,ys+ym-1
          do 30,i=xs,xs+xm-1
          row(MatStencil_i) = i  
          row(MatStencil_j) = j
          row(MatStencil_k) = k 
          if (i.eq.0 .or. j.eq.0 .or. k.eq.0 .or. i.eq.mx-1 .or.         &
     &         j.eq.my-1 .or. k.eq.mz-1) then
            v(1) = 2.d0*(HxHydHz + HxHzdHy + HyHzdHx)
            call MatSetValuesStencil(jac,i1,row,i1,row,v,INSERT_VALUES,    &
     &                               ierr)
          else 
            v(1) = -HxHydHz
             col(MatStencil_i,1) = i
             col(MatStencil_j,1) = j
             col(MatStencil_k,1) = k-1
            v(2) = -HxHzdHy
             col(MatStencil_i,2) = i
             col(MatStencil_j,2) = j-1
             col(MatStencil_k,2) = k
            v(3) = -HyHzdHx
             col(MatStencil_i,3) = i-1 
             col(MatStencil_j,3) = j
             col(MatStencil_k,3) = k
            v(4) = 2.d0*(HxHydHz + HxHzdHy + HyHzdHx)
             col(MatStencil_i,4) = i 
             col(MatStencil_j,4) = j
             col(MatStencil_k,4) = k
            v(5) = -HyHzdHx
             col(MatStencil_i,5) = i+1 
             col(MatStencil_j,5) = j 
             col(MatStencil_k,5) = k
            v(6) = -HxHzdHy
             col(MatStencil_i,6) = i 
             col(MatStencil_j,6) = j+1 
             col(MatStencil_k,6) = k
            v(7) = -HxHydHz
             col(MatStencil_i,7) = i 
             col(MatStencil_j,7) = j 
             col(MatStencil_k,7) = k+1
      call MatSetValuesStencil(jac,i1,row,i7,col,v,INSERT_VALUES,               &
     &                               ierr)
          endif
 30       continue
 20     continue
 10   continue

      call MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY,ierr)
      call MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY,ierr)
      return 
      end


More information about the petsc-users mailing list