[petsc-users] matshell for spectral methods in fortran

[Barry Smith]

   You didn't include any code related to creating or setting the MATSHELL.  What goes wrong with your Fortran code.


> On May 23, 2017, at 2:05 PM, Luv Sharma <luvsharma11 at gmail.com> wrote:
> 
> Dear PETSc team,
> 
> I am working on a code which solves mechanical equilibrium using spectral methods.
> I want to make use of the matshell to get the action J*v.
> 
> I have been able to successfully implement it using petsc4py. But having difficulties to get it working in a fortran code.
> I am using petsc-3.7.6. 
> 
> Below is a stripped down version of the existing fortran code (module). Can you please help me in figuring out how the right way to do it a code with following structure?
> 
> !--------------------------------------------------------------------------------------------------
> module spectral_mech_basic
> 
>  implicit none
>  private
> #include <petsc/finclude/petsc.h90>
> 
> ! *PETSc data here*
>  DM ..
>  SNES ..
>  ..
> contains
> 
> !--------------------------------------------------------------------------------------------------
> subroutine basicPETSc_init
> 
>  external :: &
>    *petsc functions here*
> 
> ! initialize solver specific parts of PETSc
>  call SNESCreate(PETSC_COMM_WORLD,snes,ierr); CHKERRQ(ierr)
>  call SNESSetOptionsPrefix(snes,'mech_',ierr);CHKERRQ(ierr)  
>  call DMDACreate3d(PETSC_COMM_WORLD, &
>         DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, &                                  
>         DMDA_STENCIL_BOX, &                                                                         
>         grid(1),grid(2),grid(3), &                                                          
>         1 , 1, worldsize, &
>         9, 0, &                                                                                    
>         grid(1),grid(2),localK, &                                                             
>         da,ierr)                                                                                  
>  CHKERRQ(ierr)
>  call SNESSetDM(snes,da,ierr); CHKERRQ(ierr)
>  call DMCreateGlobalVector(da,solution_vec,ierr); CHKERRQ(ierr)                                    
>  call DMDASNESSetFunctionLocal(da,INSERT_VALUES,BasicPETSC_formResidual,PETSC_NULL_OBJECT,ierr)    
>  CHKERRQ(ierr) 
>  call SNESSetDM(snes,da,ierr); CHKERRQ(ierr)                                                       
>  call SNESGetKSP(snes,ksp,ierr); CHKERRQ(ierr)                                                    
>  call KSPGetPC(ksp,pc,ierr); CHKERRQ(ierr)                                                          
>  call PCSetType(pc,PCNONE,ierr); CHKERRQ(ierr)                                                     
>  call SNESSetFromOptions(snes,ierr); CHKERRQ(ierr)                                                 
> 
> ! init fields                 
>  call DMDAVecGetArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr)                                    
>  call DMDAVecRestoreArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr)                               
> 
> end subroutine basicPETSc_init
> !--------------------------------------------------------------------------------------------------
> 
> type(tSolutionState) function &
>   basicPETSc_solution(incInfoIn,timeinc,timeinc_old,stress_BC,rotation_BC)
>  implicit none
> ! PETSc Data
>  PetscErrorCode :: ierr   
>  SNESConvergedReason :: reason
>  external :: &
>    SNESSolve, &
> ! solve BVP 
>  call SNESSolve(snes,PETSC_NULL_OBJECT,solution_vec,ierr)
>  CHKERRQ(ierr)
> end function BasicPETSc_solution
> !--------------------------------------------------------------------------------------------------
> !> @brief forms the basic residual vector
> !--------------------------------------------------------------------------------------------------
> subroutine BasicPETSC_formResidual(in,x_scal,f_scal,dummy,ierr)
> 
>  implicit none
>  DMDALocalInfo, dimension(DMDA_LOCAL_INFO_SIZE) :: &
>    in
>  PetscScalar, dimension(3,3, &
>    XG_RANGE,YG_RANGE,ZG_RANGE), intent(in) :: &
>    x_scal
>  PetscScalar, dimension(3,3, &
>    X_RANGE,Y_RANGE,Z_RANGE), intent(out) :: &
>    f_scal
> ! constructing residual
> …..
> ….
> 
>  f_scal = tensorField_real(1:3,1:3,1:grid(1),1:grid(2),1:grid3)
> 
> end subroutine BasicPETSc_formResidual
> !--------------------------------------------------------------------------------------------------
> 
> end module spectral_mech_basic
> !--------------------------------------------------------------------------------------------------
> 
> Best regards,
> Luv
> 
>> On 3 Nov 2016, at 01:17, Barry Smith <bsmith at mcs.anl.gov> wrote:
>> 
>> 
>>  Is anyone away of cases where PETSc has been used with spectral methods?
>> 
>>   Thanks
>> 
>>    Barry
>>