text:
MODULE Step3DMod

  USE bellhopMod
  USE Bdry3DMod
  USE sspMod
  IMPLICIT NONE
CONTAINS

  SUBROUTINE Step3D( ray0, ray2 )

    ! Does a single step along the ray
    ! x denotes the ray coordinate, ( x, y, z )
    ! t denotes the scaled tangent to the ray (previously (xi, eta, zeta) )
    ! c * t would be the unit tangent

    USE RayNormals

    ! rays
    TYPE( ray3DPt ) :: ray0, ray1, ray2
    INTEGER         :: iSegx0, iSegy0, iSegz0
    REAL  (KIND=8 ) :: gradc0( 3 ), gradc1( 3 ), gradc2( 3 ), &
         c0, cimag0, csq0, cxx0, cyy0, czz0, cxy0, cxz0, cyz0, cnn0, cmn0, cmm0, &
         c1, cimag1, csq1, cxx1, cyy1, czz1, cxy1, cxz1, cyz1, cnn1, cmn1, cmm1, &
         c2, cimag2,       cxx2, cyy2, czz2, cxy2, cxz2, cyz2, c_mat0( 2, 2 ), c_mat1( 2, 2 ), &
         urayt0( 3 ), urayt1( 3 ), h, halfh, hw0, hw1, w0, w1, rho
    REAL  (KIND=8 ) :: gradcjump( 3 ), csjump, cn1jump, cn2jump, tBdry( 3 ), nBdry( 3 ), RM, R1, R2, Tg, Th, &
         rayt( 3 ), rayn1( 3 ), rayn2( 3 )
    REAL  (KIND=8 ) :: e1( 3 ), e2( 3 )                              ! ray normals for ray-centered coordinates
    REAL  (KIND=8 ) :: p_tilde_in(  2 ), p_hat_in(  2 ), q_tilde_in(  2 ), q_hat_in(  2 ), &
         p_tilde_out( 2 ), p_hat_out( 2 ), RotMat( 2, 2 )

    ! The numerical integrator used here is a version of the polygon (a.k.a. midpoint, leapfrog, or Box method), and similar
    ! to the Heun (second order Runge-Kutta method).
    ! However, it's modified to allow for a dynamic step change, while preserving the second-order accuracy).

    ! *** Phase 1 (an Euler step)

    CALL EvaluateSSP3D( ray0%x, c0, cimag0, gradc0, cxx0, cyy0, czz0, cxy0, cxz0, cyz0, rho, freq, 'TAB' )
    CALL RayNormal( ray0%t, ray0%phi, c0, e1, e2 ) ! Compute ray normals e1 and e2
    CALL Get_c_partials( cxx0, cxy0, cxz0, cyy0, cyz0, czz0, e1, e2, cnn0, cmn0, cmm0 ) ! Compute second partials of c along ray normals

    csq0   = c0 * c0

    iSegx0 = iSegx     ! make note of current layer
    iSegy0 = iSegy
    iSegz0 = iSegz

    h = Beam%deltas       ! initially set the step h, to the basic one, deltas
    urayt0 = c0 * ray0%t  ! unit tangent

    CALL ReduceStep3D( ray0%x, urayt0, iSegx0, iSegy0, iSegz0, h ) ! reduce h to land on boundary

    halfh = 0.5 * h       ! first step of the modified polygon method is a half step

    ray1%x = ray0%x    + halfh * urayt0
    ray1%t = ray0%t    - halfh * gradc0 / csq0

    ! ray1%f    = ray0%f    + halfh * ( c0 * ray0%DetP - cnn0 / csq0 * ray0%DetQ )
    ! ray1%g    = ray0%g    + halfh * ( c0 * ray0%DetP - cmm0 / csq0 * ray0%DetQ )
    ! ray1%h    = ray0%h    + halfh * (                - cmn0 / csq0 * ray0%DetQ )
    ! ray1%DetP = ray0%DetP + halfh / csq0 * ( -cmm0 * ray0%f - cnn0 * ray0%g + 2.0 * cmn0 * ray0%h )
    ! ray1%DetQ = ray0%DetQ + halfh * c0 * ( ray0%f + ray0%g )

    ray1%phi  = ray0%phi  + halfh * ( 1.0 / c0 ) * ray0%t( 3 ) * &
         ( ray0%t( 2 ) * gradc0( 1 ) - ray0%t( 1 ) * gradc0( 2 ) ) / &
         ( ray0%t( 1 ) ** 2 + ray0%t( 2 ) ** 2 )

    c_mat0( 1, : ) = -[ cnn0, cmn0 ] / csq0
    c_mat0( 2, : ) = -[ cmn0, cmm0 ] / csq0

    ray1%p_tilde = ray0%p_tilde + halfh * MATMUL( c_mat0, ray0%q_tilde )
    ray1%q_tilde = ray0%q_tilde + halfh *         c0    * ray0%p_tilde

    ray1%p_hat   = ray0%p_hat   + halfh * MATMUL( c_mat0, ray0%q_hat )
    ray1%q_hat   = ray0%q_hat   + halfh *         c0    * ray0%p_hat

    ! *** Phase 2

    CALL EvaluateSSP3D( ray1%x, c1, cimag1, gradc1, cxx1, cyy1, czz1, cxy1, cxz1, cyz1, rho, freq, 'TAB' )
    rayt = c1 * ray1%t           ! unit tangent to ray
    CALL RayNormal_unit( rayt, ray2%phi, e1, e2 )
    CALL Get_c_partials( cxx1, cxy1, cxz1, cyy1, cyz1, czz1, e1, e2, cnn1, cmn1, cmm1 ) ! Compute second partials of c along ray normals

    csq1   = c1 * c1
    urayt1 = c1 * ray1%t   ! unit tangent

    CALL ReduceStep3D( ray0%x, urayt1, iSegx0, iSegy0, iSegz0, h ) ! reduce h to land on boundary

    ! use blend of f' based on proportion of a full step used.
    w1  = h / ( 2.0d0 * halfh )
    w0  = 1.0d0 - w1
    hw0 = h * w0
    hw1 = h * w1

    ray2%x = ray0%x    + hw0 * urayt0        + hw1 * urayt1
    ray2%t = ray0%t    - hw0 * gradc0 / csq0 - hw1 * gradc1 / csq1

    ! ERROR: need to do hw0 and hw1 blend here as well !!!
    ! ray2%f    = ray0%f    + h * ( c1 * ray1%DetP - cnn1 / csq1 * ray1%DetQ )
    ! ray2%g    = ray0%g    + h * ( c1 * ray1%DetP - cmm1 / csq1 * ray1%DetQ )
    ! ray2%h    = ray0%h    + h * (                - cmn1 / csq1 * ray1%DetQ )
    ! ray2%DetP = ray0%DetP + h / csq1 * ( -cmm1 * ray1%f - cnn1 * ray1%g + 2.0 * cmn1 * ray1%h )
    ! ray2%DetQ = ray0%DetQ + h * c1 * ( ray1%f + ray1%g )

    ray2%phi  = ray0%phi  + h * ( 1.0 / c1 ) * ray1%t( 3 ) * &
         ( ray1%t( 2 ) * gradc1( 1 ) - ray1%t( 1 ) * gradc1( 2 ) ) / &
         ( ray1%t( 1 ) ** 2 + ray1%t( 2 ) ** 2 )
    ray2%tau  = ray0%tau + hw0 / CMPLX( c0, cimag0, KIND=8 ) + hw1 / CMPLX( c1, cimag1, KIND=8 )

    ray2%Amp       = ray0%Amp
    ray2%Phase     = ray0%Phase
    ray2%NumTopBnc = ray0%NumTopBnc
    ray2%NumBotBnc = ray0%NumBotBnc

    c_mat1( 1, : ) = -[ cnn1, cmn1 ] / csq1
    c_mat1( 2, : ) = -[ cmn1, cmm1 ] / csq1

    ray2%p_tilde = ray0%p_tilde + hw0 * MATMUL( c_mat0, ray0%q_tilde ) + hw1 * MATMUL( c_mat1, ray1%q_tilde )
    ray2%q_tilde = ray0%q_tilde + hw0 *         c0    * ray0%p_tilde   + hw1 *         c1    * ray1%p_tilde
    ray2%p_hat   = ray0%p_hat   + hw0 * MATMUL( c_mat0, ray0%q_hat   ) + hw1 * MATMUL( c_mat1, ray1%q_hat )
    ray2%q_hat   = ray0%q_hat   + hw0 *         c0    * ray0%p_hat     + hw1 *         c1    * ray1%p_hat

    ! *** If we crossed an interface, apply jump condition ***

    CALL EvaluateSSP3D( ray2%x, c2, cimag2, gradc2, cxx2, cyy2, czz2, cxy2, cxz2, cyz2, rho, freq, 'TAB' )

    ray2%c = c2

    IF ( iSegx /= iSegx0 .OR. &
         iSegy /= iSegy0 .OR. &
         iSegz /= iSegz0 ) THEN

       gradcjump =  gradc2 - gradc0  ! this is precise only for c-linear layers

!!! what if we cross isegx, isegy, or isegz at the same time?
       IF ( iSegz /= iSegz0 ) THEN
          nBdry = [ 0D0, 0D0, -SIGN( 1.0D0, ray2%t( 3 ) ) ]   ! inward normal to layer
       ELSE IF ( iSegx /= iSegx0 ) THEN
          nBdry = [ -SIGN( 1.0D0, ray2%t( 1 ) ), 0D0, 0D0 ]   ! inward normal to x-sgement
       ELSE
          nBdry = [ 0D0, -SIGN( 1.0D0, ray2%t( 2 ) ), 0D0 ]   ! inward normal to y-sgement
       END IF
       CALL CurvatureCorrection

    END IF

  CONTAINS
    SUBROUTINE CurvatureCorrection

      ! correct p-q due to jumps in the gradient of the sound speed

      USE cross_products

      Th    = DOT_PRODUCT( ray2%t, nBdry )   ! component of ray tangent, normal to boundary
      tBdry = ray2%t - Th * nBdry            ! tangent, along the boundary, in the reflection plane
      tBdry = tBdry / NORM2( tBdry )         ! unit boundary tangent
      Tg    = DOT_PRODUCT( ray2%t, tBdry )   ! component of ray tangent, along the boundary

      rayt = c2 * ray2%t                     ! unit tangent to ray

      rayn2 = cross_product( rayt, nBdry )   ! ray tangent x boundary normal gives refl. plane normal
      rayn2 = rayn2 / NORM2( rayn2 )         ! unit normal
      rayn1 = -cross_product( rayt, rayn2 )  ! ray tangent x refl. plane normal is first ray normal

      ! normal and tangential derivatives of the sound speed
      cn1jump = DOT_PRODUCT( gradcjump, rayn1 )
      cn2jump = DOT_PRODUCT( gradcjump, rayn2 )
      csjump  = DOT_PRODUCT( gradcjump, rayt  )

      RM = Tg / Th   ! this is tan( alpha ) where alpha is the angle of incidence
      R1 = RM * ( 2 * cn1jump - RM * csjump ) / c2 ** 2
      R2 = RM * cn2jump / c2 ** 2

      ! *** curvature correction ***

      CALL RayNormal_unit( rayt, ray2%phi, e1, e2 )   ! Compute ray normals e1 and e2

      RotMat( 1, 1 ) = DOT_PRODUCT( rayn1, e1 )
      RotMat( 1, 2 ) = DOT_PRODUCT( rayn1, e2 )
      RotMat( 2, 1 ) = -RotMat( 1, 2 )             ! DOT_PRODUCT( rayn2, e1 )
      RotMat( 2, 2 ) = DOT_PRODUCT( rayn2, e2 )

      ! rotate p-q values in e1, e2 system, onto rayn1, rayn2 system

      p_tilde_in = RotMat( 1, 1 ) * ray2%p_tilde + RotMat( 1, 2 ) * ray2%p_hat
      p_hat_in   = RotMat( 2, 1 ) * ray2%p_tilde + RotMat( 2, 2 ) * ray2%p_hat

      q_tilde_in = RotMat( 1, 1 ) * ray2%q_tilde + RotMat( 1, 2 ) * ray2%q_hat
      q_hat_in   = RotMat( 2, 1 ) * ray2%q_tilde + RotMat( 2, 2 ) * ray2%q_hat

      ! here's the actual curvature change

      p_tilde_out = p_tilde_in - q_tilde_in * R1 - q_hat_in * R2
      p_hat_out   = p_hat_in   - q_tilde_in * R2

      ! rotate p back to e1, e2 system, q does not change
      ! Note RotMat^(-1) = RotMat^T

      ray2%p_tilde = RotMat( 1, 1 ) * p_tilde_out + RotMat( 2, 1 ) * p_hat_out
      ray2%p_hat   = RotMat( 1, 2 ) * p_tilde_out + RotMat( 2, 2 ) * p_hat_out

    END SUBROUTINE CurvatureCorrection

    !**********************************************************************!

    SUBROUTINE Get_c_partials( cxx, cxy, cxz, cyy, cyz, czz, e1, e2, cnn, cmn, cmm )

      ! Computes the second partials of c along ray normals

      REAL (KIND=8), INTENT( IN  ) :: cxx, cxy, cxz, cyy, cyz, czz  ! curvature of sound speed (cartesian)
      REAL (KIND=8), INTENT( IN  ) :: e1( 3 ), e2( 3 )              ! principal normals
      REAL (KIND=8), INTENT( OUT ) :: cnn, cmn, cmm                 ! curvature of sound speed (ray-centered)

      cnn = cxx * e1( 1 )**2 + cyy * e1( 2 )**2 + czz * e1( 3 )**2 + 2.0 * cxy * e1( 1 ) * e1( 2 ) + &
           2.0 * cxz * e1( 1 ) * e1( 3 ) + 2.0 * cyz * e1( 2 ) * e1( 3 )

      cmn = cxx * e1( 1 ) * e2( 1 ) + cyy * e1( 2 ) * e2( 2 ) + czz * e1( 3 ) * e2( 3 ) + &
           cxy * ( e1( 1 ) * e2( 2 ) + e2( 1 ) * e1( 2 ) ) + cxz * ( e1( 1 ) * e2( 3 ) + e2( 1 ) * e1( 3 ) ) +  &
           cyz * ( e1( 2 ) * e2( 3 ) + e2( 2 ) * e1( 3 ) )

      cmm = cxx * e2( 1 )**2 + cyy * e2( 2 )**2 + czz * e2( 3 )**2 + 2.0 * cxy * e2( 1 ) * e2( 2 ) + &
           2.0 * cxz * e2( 1 ) * e2( 3 ) + 2.0 * cyz * e2( 2 ) * e2( 3 )

      RETURN
    END SUBROUTINE Get_c_partials
  END SUBROUTINE Step3D

  ! **********************************************************************!

  SUBROUTINE ReduceStep3D( x0, urayt, iSegx0, iSegy0, iSegz0, h )

    INTEGER,       INTENT( IN    ) :: iSegx0, iSegy0, iSegz0
    REAL (KIND=8), INTENT( IN    ) :: x0( 3 ), urayt( 3 )  ! ray coordinate and tangent
    REAL (KIND=8), INTENT( INOUT ) :: h                    ! reduced step size
    REAL (KIND=8) :: d( 3 ), d0( 3 ), tri_n( 3 )
    REAL (KIND=8) :: x( 3 )                              ! ray coordinate after full trial step
    REAL (KIND=8) :: h1, h2, h3, h4, h5, h6, h7          ! step sizes
    REAL (KIND=8) :: xSeg( 2 ), ySeg( 2 )

    ! Detect interface or boundary crossing and reduce step, if necessary, to land on that crossing.
    ! Keep in mind possibility that user put source right on an interface
    ! and that multiple events can occur (crossing interface, top, and bottom in a single step).

    x = x0 + h * urayt ! make a trial step

    ! interface crossing in depth
    ! Step reduction is not done for the top or bottom layer
    ! Instead the SSP is extrapolated
    ! This prevents problems when the boundaries are outside the domain of the SSP
    h1 = huge( h1 )
    IF ( ABS( urayt( 3 ) ) > EPSILON( h1 ) ) THEN
       IF      ( SSP%z( iSegz0     ) > x(  3 ) .AND. iSegz0     > 1  ) THEN
          h1 = ( SSP%z( iSegz0     ) - x0( 3 ) ) / urayt( 3 )
          ! write( *, * ) 'layer crossing', iSegz0, h1
       ELSE IF ( SSP%z( iSegz0 + 1 ) < x(  3 ) .AND. iSegz0 + 1 < SSP%Nz ) THEN
          h1 = ( SSP%z( iSegz0 + 1 ) - x0( 3 ) ) / urayt( 3 )
          ! write( *, * ) 'layer crossing', iSegz0, h1
       END IF
    END IF

    ! top crossing
    h2 = huge( h2 )
    d  = x - Topx              ! vector from top to ray
    IF ( DOT_PRODUCT( Topn, d )  > EPSILON( h1 ) ) THEN
       d0 = x0 - Topx   ! vector from top    node to ray origin
       h2 = -DOT_PRODUCT( d0, Topn ) / DOT_PRODUCT( urayt, Topn )
       ! write( *, * ) 'top crossing'
    END IF

    ! bottom crossing
    h3 = huge( h3 )
    d  = x - Botx              ! vector from bottom to ray
    IF ( DOT_PRODUCT( Botn, d ) > EPSILON( h1 ) ) THEN
       d0 = x0 - Botx   ! vector from bottom node to ray origin
       h3 = -DOT_PRODUCT( d0, Botn ) / DOT_PRODUCT( urayt, Botn )
       ! write( *, * ) 'bottom crossing'
    END IF

    ! top/bottom segment crossing in x
    h4 = huge( h4 )
    xSeg( 1 ) = MAX( xTopSeg( 1 ), xBotSeg( 1 ) )
    xSeg( 2 ) = MIN( xTopSeg( 2 ), xBotSeg( 2 ) )

    IF ( SSP%Type == 'H' ) THEN   ! ocean segment
       xSeg( 1 ) = MAX( xSeg( 1 ), SSP%Seg%x( iSegx0     ) )
       xSeg( 2 ) = MIN( xSeg( 2 ), SSP%Seg%x( iSegx0 + 1 ) )
    END IF

    IF ( ABS( urayt( 1 ) ) > EPSILON( h1 ) ) THEN
       IF       ( x(  1 ) < xSeg( 1 ) ) THEN
          h4 = -( x0( 1 ) - xSeg( 1 ) ) / urayt( 1 )
          ! write( *, * ) 'segment crossing in x'
       ELSE IF  ( x(  1 ) > xSeg( 2 ) ) THEN
          h4 = -( x0( 1 ) - xSeg( 2 ) ) / urayt( 1 )
          ! write( *, * ) 'segment crossing in x'
       END IF
    END IF

    ! top/bottom segment crossing in y
    h5 = huge( h5 )
    ySeg( 1 ) = MAX( yTopSeg( 1 ), yBotSeg( 1 ) )
    ySeg( 2 ) = MIN( yTopSeg( 2 ), yBotSeg( 2 ) )

    IF ( SSP%Type == 'H' ) THEN   ! ocean segment
       ySeg( 1 ) = MAX( ySeg( 1 ), SSP%Seg%y( iSegy0     ) )
       ySeg( 2 ) = MIN( ySeg( 2 ), SSP%Seg%y( iSegy0 + 1 ) )
    END IF

    IF ( ABS( urayt( 2 ) ) > 1000.0 * EPSILON( h1 ) ) THEN   !!! why 1000.0 here and not for x-crossing?
       IF       ( x(  2 ) < ySeg( 1 ) ) THEN
          h5 = -( x0( 2 ) - ySeg( 1 ) ) / urayt( 2 )
          ! write( *, * ) 'segment crossing in y'
       ELSE IF  ( x(  2 ) > ySeg( 2 ) ) THEN
          h5 = -( x0( 2 ) - ySeg( 2 ) ) / urayt( 2 )
          ! write( *, * ) 'segment crossing in y'
       END IF
    END IF

    ! triangle crossing within a top segment
    h6    = huge( h6 )
    d     = x  - Topx   ! vector from bottom node to ray end
    d0    = x0 - Topx   ! vector from bottom node to ray origin
    tri_n = [ -Top_deltay, Top_deltax, 0.0d0 ]

    IF ( ( DOT_PRODUCT( tri_n, d0 ) > 0.0d0 .AND. DOT_PRODUCT( tri_n, d ) <= 0.0d0 ) .OR. &
         ( DOT_PRODUCT( tri_n, d0 ) < 0.0d0 .AND. DOT_PRODUCT( tri_n, d ) >= 0.0d0 )  ) THEN
       h6 = -DOT_PRODUCT( d0, tri_n ) / DOT_PRODUCT( urayt, tri_n )
       ! write( *, * ) 'diagonal crossing'
    END IF

    ! triangle crossing within a bottom segment
    h7    = huge( h7 )
    d     = x  - Botx   ! vector from bottom node to ray end
    d0    = x0 - Botx   ! vector from bottom node to ray origin
    tri_n = [ -Bot_deltay, Bot_deltax, 0.0d0 ]

    IF ( ( DOT_PRODUCT( tri_n, d0 ) > 0.0d0 .AND. DOT_PRODUCT( tri_n, d ) <= 0.0d0 ) .OR. &
         ( DOT_PRODUCT( tri_n, d0 ) < 0.0d0 .AND. DOT_PRODUCT( tri_n, d ) >= 0.0d0 )  ) THEN
       h7 = -DOT_PRODUCT( d0, tri_n ) / DOT_PRODUCT( urayt, tri_n )
       ! write( *, * ) 'diagonal crossing'
    END IF

    h = MIN( h, h1, h2, h3, h4, h5, h6, h7 )  ! take limit set by shortest distance to a crossing

    IF ( h < 1.0d-4 * Beam%deltas ) THEN        ! is it taking an infinitesimal step?
       h = 1.0d-5 * Beam%deltas                 ! make sure we make some motion
       iSmallStepCtr = iSmallStepCtr + 1   ! keep a count of the number of sequential small steps
    ELSE
       iSmallStepCtr = 0                   ! didn't do a small step so reset the counter
    END IF
  END SUBROUTINE ReduceStep3D

END MODULE Step3DMod

max_stars_repo_name: IvanaEscobar/sandbox

max_issues_repo_name: IvanaEscobar/sandbox

max_forks_repo_name: IvanaEscobar/sandbox


text:
cBHEADER*********************************************************************
c Copyright (c) 2008,  Lawrence Livermore National Security, LLC.
c Produced at the Lawrence Livermore National Laboratory.
c This file is part of HYPRE.  See file COPYRIGHT for details.
c
c HYPRE is free software; you can redistribute it and/or modify it under the
c terms of the GNU Lesser General Public License (as published by the Free
c Software Foundation) version 2.1 dated February 1999.
c
c $Revision: 1.7 $
cEHEADER*********************************************************************

c****************************************************************************
c HYPRE_SStruct_ls fortran interface routines
c****************************************************************************

c****************************************************************************
c                HYPRE_SStructBiCGSTAB routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABCreate
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabcreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabcreate error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABDestroy
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabdestroy(fsolver)

      integer ierr
      integer*8 fsolver

       call HYPRE_SStructBiCGSTABDestroy(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabdestroy error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetup
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetup(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

       call HYPRE_SStructBiCGSTABSetup(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetup error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSolve
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsolve(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructBiCGSTABSolve(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsolve error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetTol
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsettol(fsolver, ftol)

      integer ierr
      integer*8 fsolver
      double precision ftol

      call HYPRE_SStructBiCGSTABSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsettol error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetMinIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetminite(fsolver, fmin_iter)

      integer ierr
      integer fmin_iter
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABSetMinIter(fsolver, fmin_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetminiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetMaxIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetmaxite(fsolver, fmax_iter)

      integer ierr
      integer fmax_iter
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABSetMaxIter(fsolver, fmax_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetmaxiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetStopCrit
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetstopcr(fsolver, fstop_crit)

      integer ierr
      integer fstop_crit
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABSetStopCri(fsolver, fstop_crit, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetstopcrit error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetPrecond
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetprecon(fsolver, fprecond_id,
     1                                            fprecond)

      integer ierr
      integer fprecond_id
      integer*8 fsolver
      integer*8 fprecond

      call HYPRE_SStructBiCGSTABSetPrecond(fsolver, fprecond_id,
     1                                     fprecond, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetprecond error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetLogging
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetloggin(fsolver, flog)

      integer ierr
      integer flog
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABSetLogging(fsolver, flog, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetlogging error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABSetPrintLevel
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabsetprintl(fsolver, fprint)

      integer ierr
      integer fprint
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABSetPrintLe(fsolver, fprint, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabsetprintlevel error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABGetNumIterations
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabgetnumite(fsolver, fnumiter)

      integer ierr
      integer fnumiter
      integer*8 fsolver

      call HYPRE_SStructBiCGSTABGetNumIter(fsolver, fnumiter,
     1                                     ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabgetnumiterations error = ',
     1                                          ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABGetFinalRelativeResidualNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabgetfinalr(fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructBiCGSTABGetFinalRe(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabgetfinalrelative error = ',
     1                                                             ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructBiCGSTABGetResidual
c--------------------------------------------------------------------------
      subroutine fhypre_sstructbicgstabgetresidu(fsolver, fresidual)

      integer ierr
      integer*8 fsolver
      integer*8 fresidual

      call HYPRE_SStructBiCGSTABGetResidua(fsolver, fresidual, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructbicgstabgetresidual error = ', ierr
      endif

      return
      end

c****************************************************************************
c                HYPRE_SStructGMRES routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESCreate
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmrescreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructGMRESCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmrescreate error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESDestroy
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmresdestroy(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructGMRESDestroy(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmresdestroy error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetup
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetup(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructGMRESSetup(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetup error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSolve
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressolve(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructGMRESSolve(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressolve error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetKDim
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetkdim(fsolver, fkdim)

      integer ierr
      integer fkdim
      integer*8 fsolver

      call HYPRE_SStructGMRESSetKDim(fsolver, fkdim, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetkdim error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetTol
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressettol(fsolver, ftol)

      integer ierr
      integer*8 fsolver
      double precision ftol

      call HYPRE_SStructGMRESSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressettol error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetMinIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetminiter(fsolver, fmin_iter)

      integer ierr
      integer fmin_iter
      integer*8 fsolver

      call HYPRE_SStructGMRESSetMinIter(fsolver, fmin_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetminiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetMaxIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetmaxiter(fsolver, fmax_iter)

      integer ierr
      integer fmax_iter
      integer*8 fsolver

      call HYPRE_SStructGMRESSetMaxIter(fsolver, fmax_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetmaxiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetStopCrit
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetstopcrit(fsolver, fstop_crit)

      integer ierr
      integer fstop_crit
      integer*8 fsolver

      call HYPRE_SStructGMRESSetStopCrit(fsolver, fstop_crit, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetstopcrit error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetPrecond
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetprecond(fsolver, fprecond_id,
     1                                         fprecond)

      integer ierr
      integer fprecond_id
      integer*8 fsolver
      integer*8 fprecond

      call HYPRE_SStructGMRESSetPrecond(fsolver, fprecond_id, fprecond,
     1                                  ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetprecond error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetLogging
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetlogging(fsolver, flog)

      integer ierr
      integer flog
      integer*8 fsolver

      call HYPRE_SStructGMRESSetLogging(fsolver, flog, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetlogging error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESSetPrintLevel
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmressetprintleve(fsolver, flevel)

      integer ierr
      integer flevel
      integer*8 fsolver

      call HYPRE_SStructGMRESSetPrintLevel(fsolver, flevel, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmressetprintlevel error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESGetNumIterations
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmresgetnumiterat(fsolver, fnumiter)

      integer ierr
      integer fnumiter
      integer*8 fsolver

      call HYPRE_SStructGMRESGetNumIterati(fsolver, fnumiter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmresgetnumiteration error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESGetFinalRelativeResidualNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmresgetfinalrela(fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructGMRESGetFinalRelat(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmresgetfinalrelative error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructGMRESGetResidual
c--------------------------------------------------------------------------
      subroutine fhypre_sstructgmresgetresidual(fsolver, fresidual)

      integer ierr
      integer*8 fsolver
      integer*8 fresidual

      call HYPRE_SStructGMRESGetResidual(fsolver, fresidual, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructgmresgetresidual error = ', ierr
      endif

      return
      end

c****************************************************************************
c                HYPRE_SStructInterpreter routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructPVectorSetRandomValues
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpvectorsetrandomv(fsolver, fseed)

      integer ierr
      integer*8 fsolver
      integer*8 fseed

      call HYPRE_SStructPVectorSetRandomVa(fsolver, fseed, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpvectorsetrandomvalues error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructVectorSetRandomValues
c--------------------------------------------------------------------------
      subroutine fhypre_sstructvectorsetrandomva(fsolver, fseed)

      integer ierr
      integer*8 fsolver
      integer*8 fseed

      call HYPRE_SStructVectorSetRandomVal(fsolver, fseed, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructvectorsetrandomvalues error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSetRandomValues
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsetrandomvalues(fsolver, fseed)

      integer ierr
      integer*8 fsolver
      integer*8 fseed

      call HYPRE_SStructSetRandomValues(fsolver, fseed, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsetrandomvalues error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSetupInterpreter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsetupinterpreter(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSetupInterpreter(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsetupinterpreter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSetupMatvec
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsetupmatvec(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSetupMatvec(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsetupmatvec error = ', ierr
      endif

      return
      end

c****************************************************************************
c                HYPRE_SStructFAC routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructFACCreate
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfaccreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructFACCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfaccreate error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACDestroy2
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacdestroy2(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructFACDestroy2(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacdestroy2 error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetup2
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetup2(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructFACSetup2(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetup2 error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSolve3
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsolve3(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructFACSolve3(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsolve3 error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetTol
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsettol(fsolver, ftol)

      integer ierr
      integer*8 fsolver
      integer*8 ftol

      call HYPRE_SStructFACSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsettol error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetPLevels
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetplevels(fsolver, fnparts, fplevels)

      integer ierr
      integer*8 fsolver
      integer*8 fnparts
      integer*8 fplevels

      call HYPRE_SStructFACSetPLevels(fsolver, fnparts, fplevels, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetplevels error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetPRefinements
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetprefinement(fsolver, fnparts,
     1                                            frfactors)

      integer ierr
      integer*8 fsolver
      integer*8 fnparts
      integer*8 frfactors(3)

      call HYPRE_SStructFACSetPRefinements(fsolver, fnparts, frfactors,
     1                                     ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetprefinements error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetMaxLevels
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetmaxlevels(fsolver, fmaxlevels)

      integer ierr
      integer*8 fsolver
      integer*8 fmaxlevels

      call HYPRE_SStructFACSetMaxLevels(fsolver, fmaxlevels, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetmaxlevels error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetMaxIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetmaxiter(fsolver, fmaxiter)

      integer ierr
      integer*8 fsolver
      integer*8 fmaxiter

      call HYPRE_SStructFACSetMaxIter(fsolver, fmaxiter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetmaxiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetRelChange
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetrelchange(fsolver, frelchange)

      integer ierr
      integer*8 fsolver
      integer*8 frelchange

      call HYPRE_SStructFACSetRelChange(fsolver, frelchange, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetrelchange error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetZeroGuess
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetzeroguess(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructFACSetZeroGuess(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetzeroguess error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetNonZeroGuess
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetnonzerogues(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructFACSetNonZeroGuess(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetnonzeroguess error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetRelaxType
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetrelaxtype(fsolver, frelaxtype)

      integer ierr
      integer*8 fsolver
      integer*8 frelaxtype

      call HYPRE_SStructFACSetRelaxType(fsolver, frelaxtype, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetrelaxtype error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetNumPreRelax
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetnumprerelax(fsolver, fnumprerelax)

      integer ierr
      integer*8 fsolver
      integer*8 fnumprerelax

      call HYPRE_SStructFACSetNumPreRelax(fsolver, fnumprerelax, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetnumprerelax error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetNumPostRelax
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetnumpostrela(fsolver,
     1                                            fnumpostrelax)

      integer ierr
      integer*8 fsolver
      integer*8 fnumpostrelax

      call HYPRE_SStructFACSetNumPostRelax(fsolver, fnumpostrelax, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetnumpostrelax error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetCoarseSolverType
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetcoarsesolve(fsolver,
     1                                            fcsolvertype)

      integer ierr
      integer*8 fsolver
      integer*8 fcsolvertype

      call HYPRE_SStructFACSetCoarseSolver(fsolver, fcsolvertype,
     1                                      ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetcoarsesolvertype error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACSetLogging
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacsetlogging(fsolver, flog)

      integer ierr
      integer flog
      integer*8 fsolver

      call HYPRE_SStructFACSetLogging(fsolver, flog, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacsetlogging error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACGetNumIterations
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacgetnumiteratio(fsolver, fnumiter)

      integer ierr
      integer fnumiter
      integer*8 fsolver

      call HYPRE_SStructFACGetNumIteration(fsolver, fnumiter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacgetnumiteration error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructFACGetFinalRelativeResidualNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructfacgetfinalrelati(fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructFACGetFinalRelativ(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructfacgetfinalrelative error = ', ierr
      endif

      return
      end

c****************************************************************************
c                HYPRE_SStructPCG routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructPCGCreate
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgcreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructPCGCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgcreate error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGDestroy
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgdestroy(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructPCGDestroy(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgdestroy error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetup
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsetup(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructPCGSetup(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsetup error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSolve
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsolve(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructPCGSolve(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsolve error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetTol
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsettol(fsolver, ftol)

      integer ierr
      integer*8 fsolver
      double precision ftol

      call HYPRE_SStructPCGSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsettol error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetMaxIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsetmaxiter(fsolver, fmax_iter)

      integer ierr
      integer fmax_iter
      integer*8 fsolver

      call HYPRE_SStructPCGSetMaxIter(fsolver, fmax_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsetmaxiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetTwoNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsettwonorm(fsolver, ftwo_norm)

      integer ierr
      integer ftwo_norm
      integer*8 fsolver

      call HYPRE_SStructPCGSetTwoNorm(fsolver, ftwo_norm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsettwonorm error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetRelChange
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsetrelchange(fsolver, frel_change)

      integer ierr
      integer frel_change
      integer*8 fsolver

      call HYPRE_SStructPCGSetRelChange(fsolver, frel_change, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsetrelchange error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetPrecond
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsetprecond(fsolver, fprecond_id,
     1                                       fprecond)

      integer ierr
      integer fprecond_id
      integer*8 fsolver
      integer*8 fprecond

      call HYPRE_SStructPCGSetPrecond(fsolver, fprecond_id, fprecond,
     1                                ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsetprecond error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetLogging
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsetlogging(fsolver, flog)

      integer ierr
      integer flog
      integer*8 fsolver

      call HYPRE_SStructPCGSetLogging(fsolver, flog, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsetlogging error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGSetPrintLevel
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcgsetprintlevel(fsolver, flevel)

      integer ierr
      integer flevel
      integer*8 fsolver

      call HYPRE_SStructPCGSetPrintLevel(fsolver, flevel, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcgsetprintlevel error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGGetNumIterations
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcggetnumiteratio(fsolver, fnumiter)

      integer ierr
      integer fnumiter
      integer*8 fsolver

      call HYPRE_SStructPCGGetNumIteration(fsolver, fnumiter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcggetnumiteration error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGGetFinalRelativeResidualNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcggetfinalrelati(fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructPCGGetFinalRelativ(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcggetfinalrelative error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructPCGGetResidual
c--------------------------------------------------------------------------
      subroutine fhypre_sstructpcggetresidual(fsolver, fresidual)

      integer ierr
      integer*8 fsolver
      integer*8 fresidual

      call HYPRE_SStructPCGGetResidual(fsolver, fresidual, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructpcggetresidual error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructDiagScaleSetup
c--------------------------------------------------------------------------
      subroutine fhypre_sstructdiagscalesetup(fsolver, fA, fy, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fy
      integer*8 fx

      call HYPRE_SStructDiagScaleSetup(fsolver, fA, fy, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructdiagscalesetup error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructDiagScale
c--------------------------------------------------------------------------
      subroutine fhypre_sstructdiagscale(fsolver, fA, fy, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fy
      integer*8 fx

      call HYPRE_SStructDiagScale(fsolver, fA, fy, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructdiagscale error = ', ierr
      endif

      return
      end

c****************************************************************************
c                HYPRE_SStructSplit routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructSplitCreate
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitcreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructSplitCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitcreate error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitDestroy
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitdestroy(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSplitDestroy(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitdestroy error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSetup
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsetup(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructSplitSetup(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsetup error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSolve
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsolve(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructSplitSolve(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsolve error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSetTol
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsettol(fsolver, ftol)

      integer ierr
      integer*8 fsolver
      double precision ftol

      call HYPRE_SStructSplitSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsettol error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSetMaxIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsetmaxiter(fsolver, fmax_iter)

      integer ierr
      integer fmax_iter
      integer*8 fsolver

      call HYPRE_SStructSplitSetMaxIter(fsolver, fmax_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsetmaxiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSetZeroGuess
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsetzeroguess(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSplitSetZeroGuess(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsetzeroguess error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSetNonZeroGuess
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsetnonzerogu(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSplitSetNonZeroGue(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsetnonzeroguess error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitSetStructSolver
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitsetstructsol(fsolver, fssolver)

      integer ierr
      integer fssolver
      integer*8 fsolver

      call HYPRE_SStructSplitSetStructSolv(fsolver, fssolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitsetstructsolve error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitGetNumIterations
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitgetnumiterat(fsolver, fnumiter)

      integer ierr
      integer fnumiter
      integer*8 fsolver

      call HYPRE_SStructSplitGetNumIterati(fsolver, fnumiter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitgetnumiteration error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSplitGetFinalRelativeResidualNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsplitgetfinalrela(fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructSplitGetFinalRelat(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsplitgetfinalrelative error = ', ierr
      endif

      return
      end

c****************************************************************************
c                HYPRE_SStructSYSPFMG routines
c****************************************************************************

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGCreate
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgcreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructSysPFMGCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgcreate error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGDestroy
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgdestroy(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSysPFMGDestroy(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgdestroy error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetup
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetup(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructSysPFMGSetup(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetup error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSolve
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsolve(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructSysPFMGSolve(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsolve error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetTol
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsettol(fsolver, ftol)

      integer ierr
      integer*8 fsolver
      double precision ftol

      call HYPRE_SStructSysPFMGSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsettol error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetMaxIter
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetmaxiter(fsolver, fmax_iter)

      integer ierr
      integer fmax_iter
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetMaxIter(fsolver, fmax_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetmaxiter error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetRelChange
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetrelchang(fsolver, frel_change)

      integer ierr
      integer frel_change
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetRelChang(fsolver, frel_change, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetrelchange error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetZeroGuess
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetzerogue(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetZeroGues(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetzeroguess error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetNonZeroGuess
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetnonzero(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetNonZeroG(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetnonzeroguess error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetRelaxType
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetrelaxty(fsolver, frelax_type)

      integer ierr
      integer frelax_type
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetRelaxTyp(fsolver, frelax_type, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetrelaxtype error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetNumPreRelax
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetnumprer(fsolver,
     1                                            fnum_pre_relax)

      integer ierr
      integer fnum_pre_relax
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetNumPreRe(fsolver, fnum_pre_relax,
     1                                        ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetnumprerelax error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetNumPostRelax
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetnumpost(fsolver,
     1                                            fnum_post_relax)

      integer ierr
      integer fnum_post_relax
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetNumPostR(fsolver, fnum_post_relax,
     1                                        ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetnumpostrelax error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetSkipRelax
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetskiprel(fsolver, fskip_relax)

      integer ierr
      integer fskip_relax
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetSkipRela(fsolver, fskip_relax, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetskiprelax error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetDxyz
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetdxyz(fsolver, fdxyz)

      integer ierr
      integer*8 fsolver
      double precision fdxyz

      call HYPRE_SStructSysPFMGSetDxyz(fsolver, fdxyz, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetdxyz error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetLogging
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetlogging(fsolver, flog)

      integer ierr
      integer flog
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetLogging(fsolver, flog, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetlogging error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGSetPrintLevel
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmgsetprintle(fsolver, fprint_level)

      integer ierr
      integer fprint_level
      integer*8 fsolver

      call HYPRE_SStructSysPFMGSetPrintLev(fsolver, fprint_level,
     1                                       ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmgsetprintlevel error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGGetNumIterations
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmggetnumiter(fsolver, fnumiter)

      integer ierr
      integer fnumiter
      integer*8 fsolver

      call HYPRE_SStructSysPFMGGetNumItera(fsolver, fnumiter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmggetnumiteration error = ', ierr
      endif

      return
      end

c--------------------------------------------------------------------------
c HYPRE_SStructSysPFMGGetFinalRelativeResidualNorm
c--------------------------------------------------------------------------
      subroutine fhypre_sstructsyspfmggetfinalre(fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructSysPFMGGetFinalRel(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructsyspfmggetfinalrelative error = ', ierr
      endif

      return
      end

c*****************************************************************************
c*
c* HYPRE_SStructMaxwell interface
c*
c*****************************************************************************

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellCreate
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellcreate(fcomm, fsolver)

      integer ierr
      integer*8 fcomm
      integer*8 fsolver

      call HYPRE_SStructMaxwellCreate(fcomm, fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellcreate = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellDestroy
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwelldestroy(fsolver)

      integer ierr
      integer*8 fsolver

      call HYPRE_SStructMaxwellDestroy(fsolver, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwelldestroy = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetup
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetup (fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructMaxwellSetup(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetup = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSolve
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsolve (fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructMaxwellSolve(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsolve = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSolve2
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsolve2(fsolver, fA, fb, fx)

      integer ierr
      integer*8 fsolver
      integer*8 fA
      integer*8 fb
      integer*8 fx

      call HYPRE_SStructMaxwellSolve2(fsolver, fA, fb, fx, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsolve2 = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_MaxwellGrad
c*--------------------------------------------------------------------------
      subroutine fhypre_maxwellgrad (fgrid, fT)

      integer ierr
      integer*8 fgrid
      integer*8 fT

      call HYPRE_MaxwellGrad(fgrid, fT, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellgrad = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetGrad
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetgrad (fsolver, fT)

      integer ierr
      integer*8 fsolver
      integer*8 fT

      call HYPRE_SStructMaxwellSetGrad(fsolver, fT, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetgrad = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetRfactors
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetrfactor (fsolver,frfactors)

      integer ierr
      integer*8 fsolver
      integer*8 frfactors(3)

      call HYPRE_SStructMaxwellSetRfactors(fsolver, frfactors, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetrfactors = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetTol
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsettol (fsolver, ftol)

      integer ierr
      integer*8 fsolver
      double precision ftol

      call HYPRE_SStructMaxwellSetTol(fsolver, ftol, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsettol = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetConstantCoef
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetconstan (fsolver,
     1                                            fconstant_coef)

      integer ierr
      integer*8 fsolver
      integer fconstant_coef

      call HYPRE_SStructMaxwellSetConstant(fsolver, fconstant_coef,
     1                                     ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetconstantcoef = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetMaxIter
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetmaxiter (fsolver, fmax_iter)

      integer ierr
      integer*8 fsolver
      integer fmax_iter

      call HYPRE_SStructMaxwellSetMaxIter(fsolver, fmax_iter, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetmaxiter = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetRelChange
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetrelchan (fsolver, frel_change)

      integer ierr
      integer*8 fsolver
      integer frel_change

      call HYPRE_SStructMaxwellSetRelChang(fsolver, frel_change, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetrelchange = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetNumPreRelax
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetnumprer (fsolver,
     1                                            fnum_pre_relax)

      integer ierr
      integer*8 fsolver
      integer fnum_pre_relax

      call HYPRE_SStructMaxwellSetNumPreRe(fsolver, fnum_pre_relax,
     1                                     ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetnumprerelax = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetNumPostRelax
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetnumpost (fsolver,
     1                                            fnum_post_relax)

      integer ierr
      integer*8 fsolver
      integer fnum_post_relax

      call HYPRE_SStructMaxwellSetNumPostR(fsolver, fnum_post_relax,
     1                                     ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetnumpostrelax = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetLogging
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetlogging (fsolver, flogging)

      integer ierr
      integer*8 fsolver
      integer flogging

      call HYPRE_SStructMaxwellSetLogging(fsolver, flogging, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetlogging = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellSetPrintLevel
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellsetprintle (fsolver, fprint_level)

      integer ierr
      integer*8 fsolver
      integer flogging

      call HYPRE_SStructMaxwellSetPrintLev(fsolver, fprint_level,
     1                                     ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellsetprintlevel = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellPrintLogging
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellprintloggi (fsolver, fmyid)

      integer ierr
      integer*8 fsolver
      integer flogging

      call HYPRE_SStructMaxwellPrintLoggin(fsolver, fmyid, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellprintlogging = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellGetNumIterations
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellgetnumiter (fsolver,
     1                                            fnum_iterations)

      integer ierr
      integer*8 fsolver
      integer fnum_iterations

      call HYPRE_SStructMaxwellGetNumItera(fsolver,
     1                                     fnum_iterations, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellgetnumiterations = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellGetFinalRelativeResidualNorm
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellgetfinalre (fsolver, fnorm)

      integer ierr
      integer*8 fsolver
      double precision fnorm

      call HYPRE_SStructMaxwellGetFinalRel(fsolver, fnorm, ierr)

      if(ierr .ne. 0) then
         print *,
     1      'fhypre_sstructmaxwellgetfinalrelativeresidualnorm = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellPhysBdy
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellphysbdy (fgrid_l, fnum_levels,
     1                                         frfactors,
     2                                         fBdryRanks_ptr,
     3                                         fBdryRanksCnt_ptr)

      integer ierr
      integer*8 fgrid_l
      integer*8 frfactors
      integer*8 fBdryRanks_ptr
      integer*8 fBdryRanksCnt_ptr
      integer fnum_levels

      call HYPRE_SStructMaxwellPhysBdy(fgrid_l, fnum_levels, frfactors,
     1                                 fBdryRanks_ptr,
     2                                 fBdryRanksCnt_ptr, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellphysbdy = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellEliminateRowsCols
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwelleliminater (fparA, fnrows, frows)

      integer ierr
      integer*8 fparA
      integer*8 frows
      integer*8 fnrows

      call HYPRE_SStructMaxwellEliminateRo(fparA, fnrows, frows,
     1                                         ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwelleliminaterows = ', ierr
      endif

      return
      end

c*--------------------------------------------------------------------------
c* HYPRE_SStructMaxwellZeroVector
c*--------------------------------------------------------------------------
      subroutine fhypre_sstructmaxwellzerovector (fb, frows, fnrows)

      integer ierr
      integer*8 fb
      integer*8 frows
      integer*8 fnrows

      call HYPRE_SStructMaxwellZeroVector(fb, frows, fnrows, ierr)

      if(ierr .ne. 0) then
         print *, 'fhypre_sstructmaxwellzerovector = ', ierr
      endif

      return
      end

ext: f

max_issues_repo_name: leyel/BoSSS

max_forks_repo_name: leyel/BoSSS


text:
MODULE IO
    USE chemistry
    USE physicscore

    CHARACTER (LEN=100) :: abundSaveFile, abundLoadFile, outputFile, columnFile, outFile
    LOGICAL :: columnOutput=.False.,fullOutput=.False.,readAbunds=.False.,writeAbunds=.False.
    CHARACTER (LEN=15),ALLOCATABLE :: outSpecies(:)
    INTEGER :: nout
    INTEGER, ALLOCATABLE :: outIndx(:)

    INTEGER, PARAMETER :: outputId=10,columnId=11,abundLoadID=71,abundSaveID=72,outID=74,debugId=79,inputId=21
CONTAINS
    !Reads input reaction and species files as well as the final step of previous run if this is phase 2
    SUBROUTINE fileSetup
        IMPLICIT NONE
        INQUIRE(UNIT=columnId, OPENED=columnOutput)
        IF (columnOutput) WRITE(columnId,333) specName(outIndx)
        333 FORMAT("Time,Density,gasTemp,av,zeta,",(999(A,:,',')))

        INQUIRE(UNIT=outputId, OPENED=fullOutput)
        IF (fullOutput) THEN
            WRITE(outputId,334) fc,fo,fn,fs
            WRITE(outputId,*) "Radfield ", radfield
            WRITE(outputId,335) specName
        END IF
        335 FORMAT("Time,Density,gasTemp,av,zeta,point,",(999(A,:,',')))
        334 FORMAT("Elemental abundances, C:",1pe15.5e3," O:",1pe15.5e3," N:",1pe15.5e3," S:",1pe15.5e3)

        INQUIRE(UNIT=abundLoadID, OPENED=readAbunds)
        INQUIRE(UNIT=abundSaveID, OPENED=writeAbunds)
    END SUBROUTINE fileSetup

    SUBROUTINE readInputAbunds
        !read start file if choosing to use abundances from previous run
        IF (readAbunds) THEN
            DO l=1,points
                READ(abundLoadID,*) fhe,fc,fo,fn,fs,fmg
                READ(abundLoadID,*) abund(:nspec,l)
                REWIND(abundLoadID)
            END DO
        END IF
    END SUBROUTINE readInputAbunds

    !Writes physical variables and fractional abundances to output file, called every time step.
    SUBROUTINE output

        IF (fullOutput) THEN
            WRITE(outputId,8020) timeInYears,density(dstep),gasTemp(dstep),av(dstep),zeta,dstep,abund(:neq-1,dstep)
            8020 FORMAT(1pe11.3,',',1pe11.4,',',0pf8.2,',',1pe11.4,',',1pe11.4,',',I4,',',(999(1pe15.5,:,',')))
        END IF

        !Every 'writestep' timesteps, write the chosen species out to separate file
        !choose species you're interested in by looking at parameters.f90
        IF (writeCounter==writeStep .and. columnOutput) THEN
            writeCounter=1
            WRITE(columnId,8030) timeInYears,density(dstep),gasTemp(dstep),av(dstep),zeta,abund(outIndx,dstep)
            8030  FORMAT(1pe11.3,',',1pe11.4,',',0pf8.2,',',1pe11.4,',',1pe11.4,',',(999(1pe15.5,:,',')))
        ELSE
            writeCounter=writeCounter+1
        END IF
    END SUBROUTINE output

    SUBROUTINE finalOutput
        IF (writeAbunds) THEN
            DO dstep=1,points
                WRITE(abundSaveID,*) fhe,fc,fo,fn,fs,fmg
                WRITE(abundSaveID,8010) abund(:neq-1,dstep)
            8010  FORMAT((999(1pe15.5,:,',')))
            END DO
        END IF
    END SUBROUTINE finalOutput

    SUBROUTINE closeFiles
        CLOSE(outputId)
        CLOSE(columnId)
        CLOSE(abundSaveID)
        CLOSE(abundLoadID)
    END SUBROUTINE closeFiles

    SUBROUTINE debugout
        OPEN(debugId,file='output/debuglog',status='unknown')       !debug file.
        WRITE(debugId,*) "Integrator failed, printing relevant debugging information"
        WRITE(debugId,*) "dens",density(dstep)
        WRITE(debugId,*) "density in integration array",abund(nspec+1,dstep)
        WRITE(debugId,*) "Av", av(dstep)
        WRITE(debugId,*) "Temp", gasTemp(dstep)
        DO i=1,nreac
            if (rate(i) .ge. huge(i)) write(debugId,*) "Rate(",i,") is potentially infinite"
        END DO
    END SUBROUTINE debugout

    SUBROUTINE simpleDebug(message)
    !A very simply subroutine for debugging, just write a bunch of variables to screen
    !so we can check they're all as expected.
    !Argument message is a string to print before the variables
        CHARACTER(LEN=*) :: message
        WRITE(*,*) message
        WRITE(*,*)"endAtFinalDensity",endAtFinalDensity
        WRITE(*,*)"freefall",freefall
        WRITE(*,*)"initialDens",initialDens
        WRITE(*,*)"finalDens",finalDens
        WRITE(*,*)"initialTemp",initialTemp
        WRITE(*,*)"finalTime",finalTime
        WRITE(*,*)"rout",rout
        WRITE(*,*)"baseAv",baseAv
        WRITE(*,*) "freezeFactor",freezeFactor
        WRITE(*,*) "abstol_factor",abstol_factor
        WRITE(*,*) "neq",neq
        WRITE(*,*) "density abund",abund(neq,1)
        WRITE(*,*) "density arr",density
        WRITE(*,*) "gasTemp",gasTemp
        WRITE(*,*) "coldens",coldens
        WRITE(*,*) "av",av
    END SUBROUTINE simpleDebug

END MODULE IO

max_stars_repo_name: jonholdship/UCLCHEM

max_issues_repo_name: jonholdship/UCLCHEM

max_forks_repo_name: jonholdship/UCLCHEM


text:
	SUBROUTINE TI_MTCH ( match, dattim, times, ntime, minute,
     +			     ipos, iret )
C************************************************************************
C* TI_MTCH								*
C*									*
C* This subroutine searches for a particular date/time in a list of	*
C* date/times.  The position of the nearest date/time, determined by	*
C* time match scheme used, is returned in IPOS.  If there is no match	*
C* possible, IPOS is set to 0.						*
C*									*
C* TI_MTCH  ( MATCH, DATTIM, TIMES, NTIME, MINUTE, IPOS, IRET )		*
C*									*
C* Input parameters:							*
C*	MATCH		INTEGER		Flag for time match scheme	*
C*					   1 = exact only		*
C*					   2 = closest before		*
C*					   3 = closest after		*
C*					   4 = closest before or after	*
C*	DATTIM		CHAR*		Date/time			*
C*	TIMES (NTIME)	CHAR*		Sorted list of date/times	*
C*	NTIME		INTEGER		Number of date/times in list	*
C*	MINUTE		INTEGER		Minutes difference for match	*
C*					   0 = ignore			*
C*									*
C* Output parameters:							*
C*	IPOS		INTEGER		Position of nearest date/time	*
C*				 	  0 = not found			*
C*	IRET		INTEGER		Return code			*
C**									*
C* Log:									*
C* S. Jacobs/NCEP	 6/99	Created					*
C* R. Curtis/EAI	10/00	Changed LLMXTM to MXNMFL in tsort array *
C************************************************************************
	INCLUDE		'GEMPRM.PRM'
C*
	CHARACTER*(*)	dattim, times (*)
C*
	CHARACTER	tsort(MXNMFL+1)*20, ttemp*20, ttemp2*20
C------------------------------------------------------------------------
	iret = 0
C
C*	Initialize ipos to not found.
C
	ipos = 0
C
C*	Check the number of minutes for the difference range. If the
C*	value is 0, set the number of minutes to a large value.
C
	IF  ( minute .eq. 0 )  THEN
C
C*	    Set the range in minutes to more than 1 year.
C
	    imin = 600000
	  ELSE
	    imin = minute
	END IF
C
C*	Try for an exact match.
C
	CALL ST_FIND ( dattim, times, ntime, ipos, ier )
C
C*	If the time is not found, and the match type is not 1,
C*	continue checking.
C
	IF  ( ( ipos .eq. 0 ) .and. ( match .ne. 1 ) )  THEN
C
C*	    Add all the times, including the search time,
C*	    to a new array.
C
	    DO  i = 1, ntime
		tsort (i) = times (i)
	    END DO
	    tsort(ntime+1) = dattim
	    nsort = ntime + 1
C
C*	    Sort the new array of times.
C
	    CALL TI_SORT ( nsort, tsort, tsort, ier )
C
C*	    Find the location of the search time in the new array.
C
	    CALL ST_FIND ( dattim, tsort, nsort, jpos, ier )
C
C*	    Find the closest time before the search time.
C
	    IF  ( match .eq. 2 )  THEN
C
C*		If the search time is the first in the array,
C*		there is no valid match to be found.
C
		IF  ( jpos .ne. 1 )  THEN
		    ttemp = tsort (jpos-1)
C
C*		    Check the time difference range.
C
		    CALL TI_DIFF ( dattim, ttemp, nmin, ier )
		    IF  ( ABS(nmin) .le. imin )  THEN
C
C*			Find the closest time in the original array.
C
			CALL ST_FIND ( ttemp, times, ntime, ipos, ier )
		    END IF
		END IF
C
C*	      Find the closest time after the search time.
C
	      ELSE IF  ( match .eq. 3 )  THEN
C
C*		If the search time is the last in the array,
C*		there is no valid match to be found.
C
		IF  ( jpos .ne. nsort )  THEN
		    ttemp = tsort (jpos+1)
C
C*		    Check the time difference range.
C
		    CALL TI_DIFF ( dattim, ttemp, nmin, ier )
		    IF  ( ABS(nmin) .le. imin )  THEN
C
C*			Find the closest time in the original array.
C
			CALL ST_FIND ( ttemp, times, ntime, ipos, ier )
		    END IF
		END IF
C
C*	      Find the closest time before or after the search time.
C
	      ELSE IF  ( match .eq. 4 )  THEN
C
C*		If the search time is the first in the array,
C*		only check the after case.
C
		IF  ( jpos .ne. 1 )  THEN
		    ttemp = tsort (jpos-1)
		    CALL TI_DIFF ( dattim, ttemp, mmin, ier )
		  ELSE
		    mmin = 610000
		END IF
C
C*		If the search time is the last in the array,
C*		only check the before case.
C
		IF  ( jpos .ne. nsort )  THEN
		    ttemp2 = tsort (jpos+1)
		    CALL TI_DIFF ( dattim, ttemp2, nmin, ier )
		  ELSE
		    nmin = 610000
		END IF
C
C*		Check the before and after time differences to
C*		find the closest.
C
		IF  ( ABS(mmin) .lt. ABS(nmin) )  THEN
		    IF  ( ABS(mmin) .le. imin )  THEN
C
C*			Find the closest time in the original array.
C
			CALL ST_FIND (ttemp, times, ntime, ipos, ier)
		    END IF
		  ELSE
		    IF  ( ABS(nmin) .le. imin )  THEN
C
C*			Find the closest time in the original array.
C
			CALL ST_FIND (ttemp2, times, ntime, ipos, ier)
		    END IF
		END IF
	    END IF
C
	END IF
C*
	RETURN
	END

max_stars_repo_name: oxelson/gempak

max_issues_repo_name: oxelson/gempak

max_forks_repo_name: oxelson/gempak


text:
module constants
   use mpi ! Necessary because of use of MPI_DOUBLE_PRECISION

   implicit none

   private

   integer,  parameter :: xp = selected_real_kind(8)
   integer,  parameter :: MPI_REAL_XP = MPI_DOUBLE_PRECISION

   real(xp), parameter :: pi = 4.0_xp*atan(1.0_xp)
   real(xp), parameter :: G  = 1._xp

   integer,  parameter :: npoints = 2**10
   real(xp), parameter :: epsilon2 = npoints**(-2._xp/3._xp) / 400._xp ! We take 1/20th of the initial caracteristic distance

   real(xp), parameter :: dt = 1.e-3_xp             ! Timestep
   real(xp), parameter :: maxtime = 8.e0_xp         ! Maximum time (ad hoc)
   integer,  parameter :: maxiter = int(maxtime/dt) ! Number of iteration

   integer,  parameter :: dump_freq = 10            ! Frequency at which the system is sampled

   logical,  parameter :: flag_diag = .true.        ! Choose force/energy computation subroutine complexity: false=n², true=n²/2
   logical,  parameter :: flag_memory = .true.      ! Choose mpi logic: false=full copy of pos. arrays, true=low mem. impact
   ! The following parameter describes how much we want to spare memory. If each process uses n bytes of memory for its own data, a
   ! memory_factor = 1 allow the code to use 2*n bytes for the interactions computations. If memory_factor = k, it uses n*(1+1/k)
   ! bytes. This option should be a divider of the number N of particles per process, since it will split this one into
   ! N/memory_factor blocks of the same size.
   integer,  parameter :: memory_factor = 1         ! Only useful if flag_memory = .true.
   integer,  parameter :: MASTER = 0                ! Choose which rank should be chosen as master node

   integer :: un, una ! Unit numbers for file operation

   public :: xp, MPI_REAL_XP, &
             pi, G, &
             npoints, epsilon2, &
             dt, maxtime, maxiter, &
             dump_freq, &
             flag_diag, flag_memory, &
             memory_factor, MASTER, &
             un, una

contains

end module constants

max_stars_repo_name: cphyc/n-body

max_issues_repo_name: cphyc/n-body

max_forks_repo_name: cphyc/n-body


text:
! { dg-do run }
! PR 26499  Test write with rewind sequences to make sure buffering and
! end-of-file conditions are handled correctly.  Derived from test case by Dale
! Ranta.  Submitted by Jerry DeLisle <jvdelisle@gcc.gnu.org>.
      program test
      dimension idata(1011)
      idata = -42
      open(unit=11,form='unformatted')
        idata(1)   = -705
        idata(  1011)   = -706
       write(11)idata
        idata(1)   = -706
        idata(  1011)   = -707
       write(11)idata
        idata(1)   = -707
        idata(  1011)   = -708
       write(11)idata
       read(11,end=        1000 )idata
       call abort()
 1000  continue
       rewind 11
       read(11,end=        1001 )idata
        if(idata(1).ne. -705.or.idata(  1011).ne. -706)call abort()
 1001  continue
       close(11,status='keep')
      open(unit=11,form='unformatted')
      rewind 11
      read(11)idata
      if(idata(1).ne.-705)then
      call abort()
      endif
      read(11)idata
      if(idata(1).ne.-706)then
      call abort()
      endif
      read(11)idata
      if(idata(1).ne.-707)then
      call abort()
      endif
      close(11,status='delete')
      stop
      end

max_stars_repo_name: maurizioabba/rose

max_issues_repo_name: sujankh/rose-matlab

max_forks_repo_name: sujankh/rose-matlab


module internal

  use, intrinsic :: ISO_C_Binding

  use internal_c

  implicit none

  public

  integer, parameter :: USER_NAME_MAX = 32

  contains

  pure function c_f_string(c_string) result(f_string)

    character(kind=C_CHAR), dimension(:), intent(in) :: c_string
    character(len=:),       allocatable              :: f_string
    integer :: i

    allocate(character(len=size(c_string))::f_string)
    arrayToString: do i=1,len(f_string)
      if (c_string(i) /= C_NULL_CHAR) then
        f_string(i:i)=c_string(i)
      else
        f_string = f_string(:i-1)
        exit
      endif
    enddo arrayToString

  end function c_f_string

  pure function f_c_string(f_string) result(c_string)

    character(len=*), intent(in)                       :: f_string
    character(kind=C_CHAR), dimension(len(f_string)+1) :: c_string
    integer :: i

    do i=1,len(f_string)
      c_string(i)=f_string(i:i)
    enddo
    c_string(i) = C_NULL_CHAR

  end function f_c_string

  function getuser()

    character(len=:), allocatable :: getuser

    character(kind=C_CHAR), dimension(USER_NAME_MAX+1) :: user
    integer(C_INT)                                     :: stat

    call getuser_c(user,stat)
    if(stat /= 0) &
      error stop 'getuser: cannot determine user name'

    getuser = c_f_string(user)

  end function getuser

  function gethome()

    character(len=:), allocatable :: gethome

    character(kind=C_CHAR), dimension(:), allocatable :: home
    integer(C_INT)                                    :: stat

    allocate(home(PATH_MAX()))
    call gethome_c(home,stat)
    if(stat /= 0) &
      error stop 'gethome: cannot determine home directory'

    gethome = c_f_string(home)

  end function gethome

  ! ToDo: Check whether the header can be used in Fortran (true parameter!)
  function PATH_MAX()

    integer(C_INT) :: PATH_MAX

    PATH_MAX = PATH_MAX_c()

  end function PATH_MAX

end module internal